Los métodos de envoltura se basan en algoritmos de búsqueda codiciosos, ya que evalúan todas las combinaciones posibles de las funciones y seleccionan la combinación que produce el mejor resultado para un algoritmo de ML específico. Una desventaja de este enfoque es que probar todas las combinaciones posibles de las características puede ser computacionalmente muy costoso, particularmente si el conjunto de características es muy grande.

Existen varios tipos de métodos. En este caso, vamos a ver los siguientes:

• Selección hacia delante (Forward Selection)
  
  
• Eliminación hacia atrás (Backward Elimination)
  
  
• Selección bidireccional (Step-wise Selection)
  
  

Los pasos del algoritmo son:

• Se selecciona un nivel de significancia (normalmente 0.05)
  
  
• Se entrenan tantos modelos como características hay (cada modelo con una característica diferente)
  
  
• Se selecciona la característica con el p-value más bajo (la más significativa)
  
  
• Se vuelven a entrenar modelos añadiendo una característica cada vez (de las que no se han seleccionado en el paso anterior) al modelo anterior
  
  
• Volvemos a seleccionar la característica con el p-value más bajo y la incorporamos a nuestro modelo original y volvemos al paso anterior
  
  
• El algoritmo para cuando se alcanza el límite de las características, o no encontramos ninguna con un p-value inferior al nivel de significancia

Vamos a ver un ejemplo con el dataset de precios de las viviendas de Boston de sklearn.

Cargamos las librerías, los datos y separamos en test y entrenamiento.

import pandas as pd
import numpy as np

from sklearn.datasets import load_boston
from sklearn.model_selection import train_test_split

import statsmodels.api as sm

# Configuración warnings
# ==============================================================================
import warnings
warnings.filterwarnings('ignore')

boston = load_boston()
bos = pd.DataFrame(boston.data, columns = boston.feature_names)
bos['Price'] = boston.target
X = bos.drop("Price", 1)       # feature matrix
y = bos['Price']               # target feature
bos.shape

(506, 14)

Tenemos 506 registros con 14 características.

X_train, X_test, y_train, y_test = train_test_split(X, y, train_size   = 0.7)

Vamos a crear una función para recrear el algoritmo.

def forward_selection(data, target, significance_level=0.05):
    initial_features = data.columns.tolist()
    best_features = []
    while (len(initial_features)>0):
        remaining_features = list(set(initial_features)-set(best_features))
        new_pval = pd.Series(index=remaining_features)
        for new_column in remaining_features:
            model = sm.OLS(target, sm.add_constant(data[best_features+[new_column]])).fit()
            new_pval[new_column] = model.pvalues[new_column]
        min_p_value = new_pval.min()
        if(min_p_value<significance_level):
            best_features.append(new_pval.idxmin())
        else:
            break
    return best_features

Seleccionamos un nivel de significancia, y ejecutamos la función.

significance = 0.05
selected_features_FS = forward_selection(X_train, y_train, significance)

print(selected_features_FS)

['LSTAT', 'RM', 'PTRATIO', 'DIS', 'NOX', 'B', 'RAD', 'TAX', 'ZN', 'CRIM', 'CHAS']

Ahora ya podemos seleccionar sólo aquellas características que nos ha devuelto el algoritmo.

X_train_fs = X_train[selected_features_FS]
X_test_fs = X_test[selected_features_FS]

print(X_train.shape, X_test.shape)
print(X_train_fs.shape, X_test_fs.shape)

(354, 13) (152, 13)
(354, 11) (152, 11)

Existen algunas librerías que podemos usar para automatizar el proceso anterior, con alguna que otra diferencia. En este caso, vamos a usar la función SequentialFeatureSelector de la librería mlxtend.

Además de las librerías anteriores, necesitamos un par más para poder usar la función.

from mlxtend.feature_selection import SequentialFeatureSelector as SFS
from sklearn.linear_model import LinearRegression

from mlxtend.plotting import plot_sequential_feature_selection as plot_sfs
import matplotlib.pyplot as plt

Primero tenemos que crear nuestro transformador.

sfs = SFS(LinearRegression(),
          k_features=11,
          forward=True,
          floating=False,
          scoring = 'r2',
          cv = 0)

Al método le pasamos los siguientes parámetros:

• Modelo para evaluar las características (en nuestro caso LinearRegression())
  
  
• Número de características con las que nos queremos quedar (11)
  
  
• Si el algoritmo es de paso adelante o atrás (forward = True)
  
  
• El parámetro floating lo veremos más adelante, cuando hablemos del algoritmo selección bidireccional
  
  
• La métrica a utilizar para evaluar el modelo ($R²$)
  
  
• El último parámetro (cv = 0) tiene que ver con la validación cruzada, que veremos más adelante.

Ahora sólo nos queda entrenar al transformador (con la función fit()), y aplicarlo sobre nuestros datos (con transform()).

sfs.fit(X_train, y_train)
sfs.k_feature_names_

('CRIM',
 'ZN',
 'CHAS',
 'NOX',
 'RM',
 'DIS',
 'RAD',
 'TAX',
 'PTRATIO',
 'B',
 'LSTAT')

X_train_fs = sfs.transform(X_train)
X_test_fs = sfs.transform(X_test)

print(X_train.shape, X_test.shape)
print(X_train_fs.shape, X_test_fs.shape)

(354, 13) (152, 13)
(354, 11) (152, 11)

Como vemos, este método cambia con respecto al que hemos creado nosotros en varios aspectos, entre otros:

• La selección de características en cada paso ya no se hace en función del p-value. Ahora se entrena el modelo y se selcciona en cada paso la que mejor rendimiento ofrezca al modelo (según la métrica utilizada, en nuestro caso $R²$)
  
  
• A diferencia de nuestro método, ahora tenemos que indicar el número de características que queremos tener al final.

Este último punto puede ser problemático, ya que no podemos saber a priori qué número de características es el idóneo.

Para intentar solucionarlo, podemos indicarle al método que pruebe con un rango de ellas (lo ideal sería con todas, aunque dependiendeo del número puede ser inviable) y comprobar qué número nos da el mejor rendimiento.

sfs1 = SFS(LinearRegression(),
         k_features=(1,13),
         forward=True,
         floating=False,
         cv=0)
sfs1.fit(X_train, y_train)

Para comprobar el rendimiento de cada combinación, vamos a usar el método get_metric_dict().

sfs1.get_metric_dict()

{1: {'feature_idx': (12,),
  'cv_scores': array([0.56240129]),
  'avg_score': 0.5624012933913177,
  'feature_names': ('LSTAT',),
  'ci_bound': nan,
  'std_dev': 0.0,
  'std_err': nan},
 2: {'feature_idx': (5, 12),
  'cv_scores': array([0.65548207]),
  'avg_score': 0.6554820707483102,
  'feature_names': ('RM', 'LSTAT'),
  'ci_bound': nan,
  'std_dev': 0.0,
  'std_err': nan},
 3: {'feature_idx': (5, 10, 12),
  'cv_scores': array([0.68684585]),
  'avg_score': 0.6868458545728382,
  'feature_names': ('RM', 'PTRATIO', 'LSTAT'),
  'ci_bound': nan,
  'std_dev': 0.0,
  'std_err': nan},
 4: {'feature_idx': (3, 5, 10, 12),
  'cv_scores': array([0.70348744]),
  'avg_score': 0.7034874374066211,
  'feature_names': ('CHAS', 'RM', 'PTRATIO', 'LSTAT'),
   'PTRATIO',

También podemos visualizar los datos en una gráfica gracias a la clase plotting de la librería mlxtend.

fig1 = plot_sfs(sfs1.get_metric_dict(), kind='std_dev')
plt.title('Sequential Forward Selection (w. StdErr)')
plt.grid()
plt.show()

El eje y del gráfico anterior, indica el $R²$ del modelo.

En este caso, el proceso es el inverso del anterior. Comenzamos entrenando el modelo con todas las características, y vamos eliminando aquellas con el p-value más alto. El proceso termina cuando no existen características con p-value superior al umbral dado.

def backward_elimination(data, target, significance_level=0.05):
    features = data.columns.tolist()
    while(len(features)>0):
        features_with_constant = sm.add_constant(data[features])
        p_values = sm.OLS(target, features_with_constant).fit().pvalues[1:]
        max_p_value = p_values.max()
        if(max_p_value >= significance_level):
            excluded_feature = p_values.idxmax()
            features.remove(excluded_feature)
        else:
            break 
    return features

significance = 0.05
selected_features_BE = backward_elimination(X_train, y_train, significance)

print(selected_features_BE)

['CRIM', 'ZN', 'CHAS', 'NOX', 'RM', 'DIS', 'RAD', 'TAX', 'PTRATIO', 'B', 'LSTAT']

X_train_be = X_train[selected_features_BE]
X_test_be = X_test[selected_features_BE]

print(X_train.shape, X_test.shape)
print(X_train_be.shape, X_test_be.shape)

(354, 13) (152, 13)
(354, 11) (152, 11)

Usamos la misma librería anterior, pero en este caso, ponemos el parámetro forward = False.

sbe = SFS(LinearRegression(),
          k_features=11,
          forward=False,
          floating=False,
          scoring = 'r2',
          cv = 0)

sbe.fit(X_train, y_train)
sbe.k_feature_names_

('CRIM',
 'ZN',
 'CHAS',
 'NOX',
 'RM',
 'DIS',
 'RAD',
 'TAX',
 'PTRATIO',
 'B',
 'LSTAT')

En este caso, los pasos son los mismo que en la selección hacia delante, pero cada vez que incorporamos una característica al modelo, usamos eliminación hacia atrás para eliminar posibles características irrelevantes. Es una mezcla de los dos algoritmos vistos anteriormente.

def stepwise_selection(data, target,SL_in=0.05,SL_out = 0.05):
    initial_features = data.columns.tolist()
    best_features = []
    while (len(initial_features)>0):
        remaining_features = list(set(initial_features)-set(best_features))
        new_pval = pd.Series(index=remaining_features)
        for new_column in remaining_features:
            model = sm.OLS(target, sm.add_constant(data[best_features+[new_column]])).fit()
            new_pval[new_column] = model.pvalues[new_column]
        min_p_value = new_pval.min()
        if(min_p_value<SL_in):
            best_features.append(new_pval.idxmin())
            while(len(best_features)>0):
                best_features_with_constant = sm.add_constant(data[best_features])
                p_values = sm.OLS(target, best_features_with_constant).fit().pvalues[1:]
                max_p_value = p_values.max()
                if(max_p_value >= SL_out):
                    excluded_feature = p_values.idxmax()
                    best_features.remove(excluded_feature)
                else:
                    break 
        else:
            break
    return best_features

Necesitamos dos umbrales, el de selección hacia delante y el de eliminación hacia atrás.

s_in = 0.05
s_out = 0.05
selected_features_SS = stepwise_selection(X_train, y_train, s_in, s_out)

print(selected_features_SS)

['LSTAT', 'RM', 'PTRATIO', 'DIS', 'NOX', 'B', 'RAD', 'TAX', 'ZN', 'CRIM', 'CHAS']

Seguimos usando la misma librería, pero dejando el parámetro forward = True y floating = True.

# Sequential Forward Floating Selection(sffs)
sffs = SFS(LinearRegression(),
         k_features=11,
         forward=True,
         floating=True,
         cv=0)

sffs.fit(X_train, y_train)
sffs.k_feature_names_

('CRIM',
 'ZN',
 'CHAS',
 'NOX',
 'RM',
 'DIS',
 'RAD',
 'TAX',
 'PTRATIO',
 'B',
 'LSTAT')