Uso de LASSO desde el paquete lars (o glmnet) en R para la selección de variables

Lo siento si esta pregunta es un poco básica.

Estoy buscando utilizar la selección de variables LASSO para un modelo de regresión lineal múltiple en R. Tengo 15 predictores, uno de los cuales es categórico (¿eso causará un problema?). Después de configurar mi e utilizo los siguientes comandos:$x$$y$

model = lars(x, y)
coef(model)

Mi problema es cuando lo uso coef(model). Esto devuelve una matriz con 15 filas, con un predictor adicional agregado cada vez. Sin embargo, no hay sugerencias sobre qué modelo elegir. ¿Me he perdido algo? ¿Hay alguna manera de que pueda obtener el paquete lars para devolver solo un " mejor " modelo?

Hay otras publicaciones que sugieren usar en su glmnetlugar, pero esto parece más complicado. Un intento es el siguiente, utilizando los mismos e . ¿Me he perdido algo aquí ?: $x$$y$

cv = cv.glmnet(x, y)
model = glmnet(x, y, type.gaussian="covariance", lambda=cv$lambda.min)
predict(model, type="coefficients")

El comando final devuelve una lista de mis variables, la mayoría con un coeficiente, aunque algunas son = 0. ¿Es esta la elección correcta del " mejor " modelo seleccionado por LASSO? Si luego ajusto un modelo lineal con todas mis variables que tienen coeficientes not=0, obtengo estimaciones de coeficientes muy similares, pero ligeramente diferentes. ¿Hay alguna razón para esta diferencia? ¿Sería aceptable reajustar el modelo lineal con estas variables elegidas por LASSO y tomarlo como mi modelo final? De lo contrario, no puedo ver ningún valor p de importancia. ¿Me he perdido algo?

Hace

type.gaussian="covariance" 

¿Asegúrate de que glmnetusa regresión lineal múltiple?

¿La normalización automática de las variables afecta a los coeficientes? ¿Hay alguna manera de incluir términos de interacción en un procedimiento LASSO?

Estoy buscando usar este procedimiento más como una demostración de cómo se puede usar LASSO que para cualquier modelo que realmente se usará para cualquier inferencia / predicción importante si eso cambia algo.

Gracias por tomarse el tiempo de leer esto. Cualquier comentario general sobre LASSO / lars / glmnet también sería muy apreciado.

Usarlo glmnetes realmente fácil una vez que lo entiendes gracias a su excelente viñeta en http://web.stanford.edu/~hastie/glmnet/glmnet_alpha.html (también puedes consultar la página del paquete CRAN). En cuanto a la mejor lambda para glmnet, la regla general es usar

cvfit <- glmnet::cv.glmnet(x, y)
coef(cvfit, s = "lambda.1se")

en lugar de lambda.min.

Para hacer lo mismo, larsdebes hacerlo a mano. Aqui esta mi solucion

cv <- lars::cv.lars(x, y, plot.it = FALSE, mode = "step")
idx <- which.max(cv$cv - cv$cv.error <= min(cv$cv))
coef(lars::lars(x, y))[idx,]

Tenga en cuenta que esto no es exactamente lo mismo, porque se detiene en un nudo de lazo (cuando entra una variable) en lugar de en cualquier punto.

Tenga en cuenta que ahora glmnetes el paquete preferido, se mantiene activamente, más que nunca lars, y que ha habido preguntas sobre glmnetvs larsrespondidas anteriormente (los algoritmos utilizados difieren).

En cuanto a su pregunta de usar el lazo para elegir variables y luego ajustar OLS, es un debate en curso. Google for OLS publica Lasso y hay algunos documentos que discuten el tema. Incluso los autores de Elementos de aprendizaje estadístico admiten que es posible.

Editar : Aquí está el código para reproducir con mayor precisión lo que glmnethace enlars

  cv <- lars::cv.lars(x, y, plot.it = FALSE)
  ideal_l1_ratio <- cv$index[which.max(cv$cv - cv$cv.error <= min(cv$cv))]
  obj <- lars::lars(x, y)
  scaled_coefs <- scale(obj$beta, FALSE, 1 / obj$normx)
  l1 <- apply(X = scaled_coefs, MARGIN = 1, FUN = function(x) sum(abs(x)))
  coef(obj)[which.max(l1 / tail(l1, 1) > ideal_l1_ratio),]

Estoy volviendo a esta pregunta hace un tiempo, ya que creo que he resuelto la solución correcta.

Aquí hay una réplica usando el conjunto de datos mtcars:

library(glmnet)
`%ni%`<-Negate(`%in%')
data(mtcars)

x<-model.matrix(mpg~.,data=mtcars)
x=x[,-1]

glmnet1<-cv.glmnet(x=x,y=mtcars$mpg,type.measure='mse',nfolds=5,alpha=.5)

c<-coef(glmnet1,s='lambda.min',exact=TRUE)
inds<-which(c!=0)
variables<-row.names(c)[inds]
variables<-variables[variables %ni% '(Intercept)']

'variables' le da la lista de las variables que resuelven la mejor solución.

Quizás la comparación con la regresión por pasos de selección directa ayude (consulte el siguiente enlace a un sitio de uno de los autores http://www-stat.stanford.edu/~tibs/lasso/simple.html) Este es el enfoque utilizado en el Capítulo 3.4.4 de Los elementos del aprendizaje estadístico (disponible en línea de forma gratuita). Pensé que el Capítulo 3.6 en ese libro ayudó a comprender la relación entre mínimos cuadrados, mejor subconjunto y lazo (más un par de otros procedimientos). También me parece útil tomar la transposición del coeficiente, t (coef (modelo)) y write.csv, para poder abrirlo en Excel junto con una copia del gráfico (modelo) en el lateral. Es posible que desee ordenar por la última columna, que contiene la estimación de mínimos cuadrados. Luego puede ver claramente cómo se agrega cada variable en cada paso por partes y cómo cambian los coeficientes como resultado. Por supuesto, esta no es toda la historia, pero espero que sea un comienzo.

larsy glmnetoperar en matrices crudas. Para incluir los términos de interacción, deberá construir las matrices usted mismo. Eso significa una columna por interacción (que es por nivel por factor si tiene factores). Mire lm()para ver cómo lo hace (advertencia: habrá dragones).

Para hacerlo ahora, haga algo como: Para hacer un término de interacción manualmente, podría (pero quizás no debería , porque es lento) hacer:

int = D["x1"]*D["x2"]
names(int) = c("x1*x2")
D = cbind(D, int)

Luego, para usar esto en lars (suponiendo que tenga una ypatada):

lars(as.matrix(D), as.matrix(y))

Desearía poder ayudarte más con las otras preguntas. Encontré este porque lars me da pena y la documentación en él y en la web es muy escasa.

LARS resuelve TODO el camino de la solución. La ruta de la solución es lineal por partes: hay un número finito de puntos de "muesca" (es decir, valores del parámetro de regularización) en los que cambia la solución.

Entonces, la matriz de soluciones que está obteniendo son todas las soluciones posibles. En la lista que devuelve, también debería proporcionarle los valores del parámetro de regularización.