Ridge penalizou GLMs usando aumento de linha?

Eu li que a regressão de crista pode ser alcançada simplesmente adicionando linhas de dados à matriz de dados original, onde cada linha é construída usando 0 para as variáveis ​​dependentes e a raiz quadrada de $k$ ou zero para as variáveis ​​independentes. Uma linha extra é então adicionada para cada variável independente.

Eu queria saber se é possível obter uma prova para todos os casos, incluindo regressão logística ou outros GLMs.

A regressão de Ridge minimiza .$\sum_{i=1}^n (y_i-x_i^T\beta)^2+\lambda\sum_{j=1}^p\beta_j^2$

(Geralmente, é necessária uma constante, mas não reduzida. Nesse caso, ela está incluída no e nos preditores - mas se você não deseja reduzi-la, não possui uma linha correspondente para a pseudo observação. você quer reduzi-lo, você faz tem uma linha para ele. Eu vou escrevê-lo como se ele não é contado no p , e não encolhido, como é o caso mais complicado. o outro caso é uma mudança trivial a partir deste. )$\beta$$p$

Podemos escrever o segundo termo como pseudo-observações se pudermos escrever cada "y" e cada um dos vetores ( p + 1 ) correspondentes "x", de modo que$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$

Então

$(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$

Isso funciona para regressão linear. Não funciona para a regressão logística, porque a regressão logística comum não minimiza uma soma de resíduos quadrados.

[A regressão de Ridge não é a única coisa que pode ser feita por meio de truques de pseudo-observação - eles aparecem em vários outros contextos]

A generalização desta receita para os GLMs não é realmente difícil, pois os GLMs geralmente são adequados usando mínimos quadrados com ponderação iterativa . Portanto, dentro de cada iteração, pode-se substituir a etapa regular de mínimos quadrados ponderada com uma etapa de mínimos quadrados ponderada penalizada de crista para obter um GLM penalizado de crista. De fato, em combinação com penalidades adaptativas da crista, esta receita é usada para ajustar GLMs penalizados por L0 (também conhecido como melhor subconjunto, ou seja, GLMs em que o número total de coeficientes diferentes de zero é penalizado). Isso foi implementado, por exemplo, no pacote l0ara , consulte este documento e este para obter detalhes.

Também é importante notar que a maneira mais rápida e fechada de resolver uma regressão regular de cume é usar

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

para o caso em que n>=p, ou usando

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

quando p>ne para um modelo sem interceptação.

Isso é mais rápido do que usar a receita de aumento de linha , ou seja, fazer

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))))
}

Se você precisar de restrições de não-negatividade em seus coeficientes ajustados , basta fazer

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
}

que então dá um resultado um pouco mais preciso entre

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 
}

(e estritamente falando, apenas a solução nnls(A=crossprod(X)+diag(lambdas), b=crossprod(X,Y))$x é a correta).

Ainda não descobri como o caso restrito da não-negatividade poderia ser otimizado ainda mais para o p > ncaso - deixe-me saber se alguém saberia como fazer isso ... [ lmridge_nnls_largep = function (X, Y, lambda) t(X) %*% nnls(A=tcrossprod(X)+lambda*diag(nrow(X)), b=Y)$xnão funciona]