8.4 Selección de variables explicativas
Cuando se dispone de un conjunto grande de posibles variables explicativas suele ser especialmente importante determinar cuales de estas deberían ser incluidas en el modelo de regresión. Si alguna de las variables no contiene información relevante sobre la respuesta no se debería incluir (se simplificaría la interpretación del modelo, aumentaría la precisión de la estimación y se evitarían problemas como la multicolinealidad). Se trataría entonces de conseguir un buen ajuste con el menor número de variables explicativas posible.
Para actualizar un modelo (p.e. eliminando o añadiendo variables) se puede emplear la función update:
modelo.completo <- lm(fidelida ~ . , data = datos)
summary(modelo.completo)##
## Call:
## lm(formula = fidelida ~ ., data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.3351 -2.0733 0.5224 2.9218 6.7106
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -9.5935 4.8213 -1.990 0.0496 *
## velocida -0.6023 1.9590 -0.307 0.7592
## precio -1.0771 2.0283 -0.531 0.5967
## flexprec 3.4616 0.3997 8.660 1.62e-13 ***
## imgfabri -0.1735 0.6472 -0.268 0.7892
## servconj 9.0919 3.8023 2.391 0.0189 *
## imgfvent 1.5596 0.9221 1.691 0.0942 .
## calidadp 0.4874 0.3451 1.412 0.1613
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.281 on 91 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.7869, Adjusted R-squared: 0.7705
## F-statistic: 48 on 7 and 91 DF, p-value: < 2.2e-16modelo.reducido <- update(modelo.completo, . ~ . - imgfabri)
summary(modelo.reducido)##
## Call:
## lm(formula = fidelida ~ velocida + precio + flexprec + servconj +
## imgfvent + calidadp, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -13.2195 -2.0022 0.4724 2.9514 6.8328
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -9.9900 4.5656 -2.188 0.0312 *
## velocida -0.5207 1.9254 -0.270 0.7874
## precio -1.0017 1.9986 -0.501 0.6174
## flexprec 3.4709 0.3962 8.761 9.23e-14 ***
## servconj 8.9111 3.7230 2.394 0.0187 *
## imgfvent 1.3699 0.5883 2.329 0.0221 *
## calidadp 0.4844 0.3432 1.411 0.1615
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.26 on 92 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.7867, Adjusted R-squared: 0.7728
## F-statistic: 56.56 on 6 and 92 DF, p-value: < 2.2e-16Para obtener el modelo “óptimo” lo ideal sería evaluar todos los modelos posibles.
8.4.1 Búsqueda exhaustiva
La función regsubsets del paquete leaps permite seleccionar los mejores modelos
fijando el número de variables explicativas.
Por defecto, evalúa todos los modelos posibles con un determinado número de
parámetros (variando desde 1 hasta un máximo de nvmax=8)
y selecciona el mejor (nbest=1).
library(leaps)
res <- regsubsets(fidelida ~ . , data = datos)
summary(res)## Subset selection object
## Call: regsubsets.formula(fidelida ~ ., data = datos)
## 7 Variables (and intercept)
## Forced in Forced out
## velocida FALSE FALSE
## precio FALSE FALSE
## flexprec FALSE FALSE
## imgfabri FALSE FALSE
## servconj FALSE FALSE
## imgfvent FALSE FALSE
## calidadp FALSE FALSE
## 1 subsets of each size up to 7
## Selection Algorithm: exhaustive
## velocida precio flexprec imgfabri servconj imgfvent calidadp
## 1 ( 1 ) " " " " " " " " "*" " " " "
## 2 ( 1 ) " " " " "*" " " "*" " " " "
## 3 ( 1 ) " " " " "*" " " "*" "*" " "
## 4 ( 1 ) " " " " "*" " " "*" "*" "*"
## 5 ( 1 ) " " "*" "*" " " "*" "*" "*"
## 6 ( 1 ) "*" "*" "*" " " "*" "*" "*"
## 7 ( 1 ) "*" "*" "*" "*" "*" "*" "*"# names(summary(res))Al representar el resultado se obtiene un gráfico con los mejores modelos ordenados
según el criterio determinado por el argumento scale = c("bic", "Cp", "adjr2", "r2").
Por ejemplo, en este caso, empleando el coeficiente de determinación ajustado, obtendríamos:
plot(res, scale = "adjr2")
En este caso (considerando que una mejora del 2% no es significativa), el modelo resultante sería:
lm(fidelida ~ servconj + flexprec, data = hatco)##
## Call:
## lm(formula = fidelida ~ servconj + flexprec, data = hatco)
##
## Coefficients:
## (Intercept) servconj flexprec
## -3.462 7.829 3.402Notas:
Si se emplea alguno de los criterios habituales, el mejor modelo con un determinado número de variables no depende del criterio empleado. Pero estos criterios pueden diferir al comparar modelos con distinto número de variables explicativas.
Si el número de variables explicativas es grande, en lugar de emplear una búsqueda exhaustiva se puede emplear un criterio por pasos, mediante el argumento
method = c("backward", "forward", "seqrep"), pero puede ser recomendable emplear el paqueteMASSpara obtener directamente el modelo final.
8.4.2 Selección por pasos
Si el número de variables es grande (no sería práctico evaluar todas las posibilidades) se suele utilizar alguno (o varios) de los siguientes métodos:
Selección progresiva (forward): Se parte de una situación en la que no hay ninguna variable y en cada paso se incluye una aplicando un criterio de entrada (hasta que ninguna de las restantes lo verifican).
Eliminación progresiva (backward): Se parte del modelo con todas las variables y en cada paso se elimina una aplicando un criterio de salida (hasta que ninguna de las incluidas lo verifican).
Regresión paso a paso (stepwise): El más utilizado, se combina un criterio de entrada y uno de salida. Normalmente se parte sin ninguna variable y en cada paso puede haber una inclusión y una exclusión (forward/backward).
La función stepAIC del paquete MASS permite seleccionar el modelo por pasos,
hacia delante o hacia atrás según criterio AIC o BIC (también esta disponible una función step del paquete base stats con menos opciones).
La función stepwise del paquete RcmdrMisc es una interfaz de stepAIC
que facilita su uso:
library(MASS)
library(RcmdrMisc)
modelo <- stepwise(modelo.completo, direction = "forward/backward", criterion = "BIC")##
## Direction: forward/backward
## Criterion: BIC
##
## Start: AIC=437.24
## fidelida ~ 1
##
## Df Sum of Sq RSS AIC
## + servconj 1 3813.6 4013.2 375.71
## + velocida 1 3558.5 4268.2 381.81
## + flexprec 1 2615.5 5211.3 401.57
## + imgfvent 1 556.9 7269.9 434.53
## + imgfabri 1 394.2 7432.5 436.72
## <none> 7826.8 437.24
## + calidadp 1 325.8 7501.0 437.63
## + precio 1 46.2 7780.6 441.25
##
## Step: AIC=375.71
## fidelida ~ servconj
##
## Df Sum of Sq RSS AIC
## + flexprec 1 2175.6 1837.6 302.97
## + precio 1 831.5 3181.7 357.32
## + velocida 1 772.3 3240.9 359.15
## + calidadp 1 203.8 3809.4 375.15
## <none> 4013.2 375.71
## + imgfvent 1 74.8 3938.4 378.44
## + imgfabri 1 2.3 4010.9 380.25
## - servconj 1 3813.6 7826.8 437.24
##
## Step: AIC=302.97
## fidelida ~ servconj + flexprec
##
## Df Sum of Sq RSS AIC
## + imgfvent 1 129.8 1707.7 300.31
## <none> 1837.6 302.97
## + imgfabri 1 69.3 1768.3 303.76
## + calidadp 1 50.7 1786.9 304.80
## + precio 1 0.2 1837.4 307.56
## + velocida 1 0.0 1837.5 307.57
## - flexprec 1 2175.6 4013.2 375.71
## - servconj 1 3373.7 5211.3 401.57
##
## Step: AIC=300.31
## fidelida ~ servconj + flexprec + imgfvent
##
## Df Sum of Sq RSS AIC
## <none> 1707.7 300.31
## - imgfvent 1 129.82 1837.6 302.97
## + calidadp 1 24.70 1683.0 303.47
## + precio 1 0.96 1706.8 304.85
## + imgfabri 1 0.66 1707.1 304.87
## + velocida 1 0.41 1707.3 304.88
## - flexprec 1 2230.67 3938.4 378.44
## - servconj 1 2850.14 4557.9 392.91summary(modelo)##
## Call:
## lm(formula = fidelida ~ servconj + flexprec + imgfvent, data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -12.9301 -2.1395 0.0695 2.9632 7.4286
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -6.7761 3.1343 -2.162 0.0331 *
## servconj 7.4320 0.5902 12.592 <2e-16 ***
## flexprec 3.4503 0.3097 11.140 <2e-16 ***
## imgfvent 1.5369 0.5719 2.687 0.0085 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 4.24 on 95 degrees of freedom
## (1 observation deleted due to missingness)
## Multiple R-squared: 0.7818, Adjusted R-squared: 0.7749
## F-statistic: 113.5 on 3 and 95 DF, p-value: < 2.2e-16Los métodos disponibles son "backward/forward", "forward/backward", "backward" y "forward".
Cuando el número de variables explicativas es muy grande (o si el tamaño de la muestra es pequeño en comparación) pueden aparecer problemas al emplear los métodos anteriores (incluso pueden no ser aplicables). Una alternativa son los métodos de regularización (Ridge regression, Lasso) disponibles en el paquete glmnet.