Ridge ukarał GLM za pomocą powiększania rzędów?

Czytałem, że regresję grzbietu można osiągnąć, po prostu dodając wiersze danych do oryginalnej macierzy danych, gdzie każdy wiersz jest konstruowany przy użyciu 0 dla zmiennych zależnych i pierwiastka kwadratowego $k$ lub zero dla zmiennych niezależnych. Następnie dodaje się jeden dodatkowy wiersz dla każdej niezależnej zmiennej.

Zastanawiałem się, czy można uzyskać dowód dla wszystkich przypadków, w tym dla regresji logistycznej lub innych GLM.

$\sum_{i=1}^n (y_i-x_i^T\beta)^2+\lambda\sum_{j=1}^p\beta_j^2$

$\beta$$p$

$p$$(p+1)$

$(y_{n+j}-x_{n+j}^T\beta)^2=\lambda\beta_j^2\,,\quad j=1,\ldots,p$

$y_{n+j}=0$$x_{n+j,j}=\sqrt{\lambda}$$x_{n+j,k}=0$$x_{n+j,0}=0$

Następnie

$(y_{n+j}-[x_{n+j,0}\beta_0+x_{n+j,1}\beta_1+x_{n+j,2}\beta_2+...+x_{n+j,p}\beta_p])^2=\lambda\beta_j^2$

Działa to w przypadku regresji liniowej. Nie działa w przypadku regresji logistycznej, ponieważ zwykła regresja logistyczna nie minimalizuje sumy kwadratów reszt.

[Regresja grzbietu nie jest jedyną rzeczą, którą można zrobić za pomocą takich pseudoobserwacyjnych sztuczek - pojawiają się one w wielu innych kontekstach]

Uogólnienie tego przepisu na GLM rzeczywiście nie jest trudne, ponieważ GLM są zwykle dopasowane przy użyciu iteracyjnie przeważonych najmniejszych kwadratów . Dlatego w ramach każdej iteracji można zastąpić regularny ważony krok najmniejszych kwadratów krokiem karnym ważonym krokiem najmniejszych kwadratów, aby uzyskać punkt GLM karany kalenicą. W rzeczywistości, w połączeniu z adaptacyjnymi karami kalenicowymi, przepis ten stosuje się, aby dopasować GLM z karą L0 (czyli najlepszy podzbiór, tj. GLM, w których karana jest całkowita liczba niezerowych współczynników). Zostało to zaimplementowane na przykład w pakiecie L0ara , zobacz ten dokument i ten, aby uzyskać szczegółowe informacje.

Warto również zauważyć, że stosuje się najszybszy zamknięty sposób rozwiązywania regularnej regresji grzbietu

lmridge_solve = function (X, y, lambda, intercept = TRUE) {
  if (intercept) {
    lambdas = c(0, rep(lambda, ncol(X)))
    X = cbind(1, X)
  } else { lambdas = rep(lambda, ncol(X)) }
  solve(crossprod(X) + diag(lambdas), crossprod(X, y))[, 1]
}

w przypadku n>=p, gdy lub przy użyciu

lmridge_solve_largep = function (X, Y, lambda) (t(X) %*% solve(tcrossprod(X)+lambda*diag(nrow(X)), Y))[,1]

kiedy p>ni dla modelu bez przechwytywania.

Jest to szybsze niż stosowanie przepisu powiększania wierszy , tj. Robienie

lmridge_rbind = function (X, y, lambda, intercept = TRUE) {
  if (intercept) {
    lambdas = c(0, rep(lambda, ncol(X)))
    X = cbind(1, X)
  } else { lambdas = rep(lambda, ncol(X)) }
  qr.solve(rbind(X, diag(sqrt(lambdas))), c(y, rep(0, ncol(X))))
}

Jeśli będziesz potrzebować ograniczeń nieujemności względem dopasowanych współczynników , możesz to zrobić

library(nnls)

nnlmridge_solve = function (X, y, lambda, intercept = TRUE) {
  if (intercept) {
    lambdas = c(0, rep(lambda, ncol(X)))
    X = cbind(1, X)
  } else { lambdas = rep(lambda, ncol(X)) }
  nnls(A=crossprod(X)+diag(lambdas), b=crossprod(X,Y))$x
}

co daje nieco dokładniejszy wynik btw niż

nnlmridge_rbind = function (X, y, lambda, intercept = TRUE) {
  if (intercept) {
    lambdas = c(0, rep(lambda, ncol(X)))
    X = cbind(1, X)
  } else { lambdas = rep(lambda, ncol(X)) }
  nnls(A=rbind(X,diag(sqrt(lambdas))), b=c(Y,rep(0,ncol(X))))$x 
}

(i ściśle mówiąc, tylko rozwiązanie nnls(A=crossprod(X)+diag(lambdas), b=crossprod(X,Y))$x jest wtedy poprawne).

Jeszcze nie zorientowałem się, w jaki sposób można zoptymalizować przypadek ograniczony nieegatywnością dla tej p > nsprawy - daj mi znać, jeśli ktoś będzie wiedział, jak to zrobić ... [ lmridge_nnls_largep = function (X, Y, lambda) t(X) %*% nnls(A=tcrossprod(X)+lambda*diag(nrow(X)), b=Y)$xnie działa]