Używanie LASSO z pakietu lars (lub glmnet) w R do wyboru zmiennych

Przepraszam, jeśli to pytanie jest trochę podstawowe.

Chciałbym użyć selekcji zmiennych LASSO dla modelu wielokrotnej regresji liniowej w R. Mam 15 predyktorów, z których jeden jest kategoryczny (czy to spowoduje problem?). Po ustawieniu mojego i Y używam następujące polecenia:$x$$y$

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

Mój problem polega na tym, kiedy używam coef(model). Zwraca macierz z 15 wierszami, z każdym dodatkowym predyktorem dodawanym za każdym razem. Jednak nie ma sugestii, który model wybrać. Czy coś przeoczyłem? Czy istnieje sposób, aby uzyskać pakiet Larsa, aby zwrócił tylko jeden „ najlepszy ” model?

Istnieją inne posty sugerujące użycie glmnetzamiast tego, ale wydaje się to bardziej skomplikowane. Próba jest następujący, przy użyciu tego samego i y . Czy coś tu przegapiłem ?: $x$$y$

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

Ostatnie polecenie zwraca listę moich zmiennych, większość ze współczynnikiem, chociaż niektóre mają wartość = 0. Czy to właściwy wybór „ najlepszego ” modelu wybranego przez LASSO? Jeśli następnie dopasuję model liniowy ze wszystkimi moimi zmiennymi, które miały współczynniki not=0, otrzymam bardzo podobne, ale nieco inne, oszacowania współczynników. Czy jest powód tej różnicy? Czy akceptowalny byłby remont modelu liniowego z tymi zmiennymi wybranymi przez LASSO i uznanie go za mój ostateczny model? W przeciwnym razie nie widzę żadnych wartości p dla istotności. Czy coś przeoczyłem?

Robi

type.gaussian="covariance" 

upewnić się, że glmnetużywa wielu regresji liniowej?

Czy automatyczna normalizacja zmiennych wpływa w ogóle na współczynniki? Czy jest jakiś sposób na uwzględnienie warunków interakcji w procedurze LASSO?

Chciałbym użyć tej procedury bardziej jako demonstracji tego, jak można użyć LASSO, niż w jakimkolwiek modelu, który faktycznie zostanie użyty do każdego ważnego wnioskowania / prognozy, jeśli to cokolwiek zmieni.

Dziękujemy za poświęcenie czasu na przeczytanie tego. Wszelkie ogólne uwagi na temat LASSO / lars / glmnet będą również mile widziane.

Korzystanie glmnetjest bardzo łatwe, gdy się go zrozumie, dzięki doskonałej winiecie w http://web.stanford.edu/~hastie/glmnet/glmnet_alpha.html (możesz również sprawdzić stronę pakietu CRAN). Jeśli chodzi o najlepszą lambda glmnet, zasadą jest użycie

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

zamiast lambda.min.

Aby zrobić to samo lars, musisz to zrobić ręcznie. Oto moje rozwiązanie

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,]

Pamiętaj, że nie jest to dokładnie to samo, ponieważ zatrzymuje się na węźle lasso (gdy wchodzi zmienna) zamiast w dowolnym punkcie.

Należy pamiętać, że glmnetjest to teraz preferowany pakiet, jest on aktywnie utrzymywany, bardziej niż lars, i że wcześniej pojawiły się pytania dotyczące odpowiedzi glmnetvs lars(stosowane algorytmy różnią się).

Jeśli chodzi o twoje pytanie o użycie lasso do wyboru zmiennych, a następnie dopasowania OLS, jest to ciągła debata. Google for OLS po Lasso i jest kilka artykułów na ten temat. Nawet autorzy Elements of Statistics Learning przyznają, że jest to możliwe.

Edycja : Oto kod do dokładniejszego odtwarzania tego, co glmnetdziałalars

  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),]

Powracam do tego pytania jakiś czas temu, ponieważ myślę, że rozwiązałem właściwe rozwiązanie.

Oto replika korzystająca z zestawu danych 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)']

„zmienne” daje listę zmiennych, które rozwiązują najlepsze rozwiązanie.

Być może porównanie z regresją stopniową wyboru do przodu pomoże (patrz poniższy link do strony jednego z autorów http://www-stat.stanford.edu/~tibs/lasso/simple.html). Takie podejście zastosowano w rozdziale 3.4.4 elementów uczenia statystycznego (dostępne online za darmo). Pomyślałem, że rozdział 3.6 tej książki pomógł zrozumieć związek między najmniejszymi kwadratami, najlepszym podzbiorem i lasso (plus kilka innych procedur). Przydaje mi się również transpozycja współczynnika t (coef (model)) i write.csv, dzięki czemu mogę otworzyć go w programie Excel wraz z kopią wykresu (model) z boku. Możesz posortować według ostatniej kolumny, która zawiera oszacowanie najmniejszych kwadratów. Następnie możesz wyraźnie zobaczyć, w jaki sposób każda zmienna jest dodawana na każdym kroku i jak w rezultacie zmieniają się współczynniki. Oczywiście to nie jest cała historia, ale mam nadzieję, że to będzie początek.

larsi glmnetdziałają na surowych matrycach. Aby uwzględnić warunki interakcji, musisz sam zbudować macierze. Oznacza to jedną kolumnę na interakcję (która jest na poziom na czynnik, jeśli masz czynniki). Sprawdź, lm()jak to robi (ostrzeżenie: są smoki).

Aby zrobić to teraz, zrób coś takiego: Aby ręcznie dokonać określenie interakcji, to może (ale może nie powinienem , bo to spowalnia) zrobić:

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

Następnie użyj tego w Larach (zakładając, że masz do czynienia z ykopaniem):

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

Chciałbym pomóc Ci w przypadku innych pytań. Znalazłem ten, ponieważ lars sprawia mi ból, a dokumentacja w nim i w Internecie jest bardzo cienka.

LARS rozwiązuje CAŁĄ ścieżkę rozwiązania. Ścieżka rozwiązania jest fragmentarycznie liniowa - istnieje skończona liczba punktów „wycięcia” (tj. Wartości parametru regularyzacji), w których zmienia się rozwiązanie.

Zatem macierz otrzymywanych rozwiązań to wszystkie możliwe rozwiązania. Na liście, którą zwraca, powinien również podać wartości parametru regularyzacji.