¿Ridge penalizó los GLM usando el aumento de filas?

He leído que la regresión de crestas podría lograrse simplemente agregando filas de datos a la matriz de datos original, donde cada fila se construye usando 0 para las variables dependientes y la raíz cuadrada de o cero para las variables independientes. Luego se agrega una fila adicional para cada variable independiente.$k$

Me preguntaba si es posible obtener una prueba para todos los casos, incluida la regresión logística u otros GLM.

La regresión de cresta minimiza .$\sum_{i=1}^n (y_i-x_i^T\beta)^2+\lambda\sum_{j=1}^p\beta_j^2$

(A menudo se requiere una constante, pero no se reduce. En ese caso, se incluye en el$\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$

Luego

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

Esto funciona para la regresión lineal. No funciona para la regresión logística, porque la regresión logística ordinaria no minimiza una suma de residuos al cuadrado.

[La regresión de cresta no es lo único que se puede hacer a través de estos trucos de pseudoobservación: surgen en varios otros contextos]

De hecho, generalizar esta receta a GLM no es difícil, ya que los GLM generalmente se ajustan utilizando mínimos cuadrados repesados ​​de forma iterativa . Por lo tanto, dentro de cada iteración, se puede sustituir el paso regular de mínimos cuadrados ponderados con un paso de mínimos cuadrados ponderado de cresta para obtener un GLM penalizado de cresta. De hecho, en combinación con penalizaciones de cresta adaptativa, esta receta se utiliza para ajustar GLM penalizados L0 (también conocido como el mejor subconjunto, es decir, GLM donde se penaliza el número total de coeficientes distintos de cero). Esto se ha implementado, por ejemplo, en el paquete l0ara , consulte este documento y este para obtener más detalles.

También vale la pena señalar que la forma más rápida de resolver una regresión de cresta regular es 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 el caso donde n>=p, o usando

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

cuando p>ny para un modelo sin intercepción.

Esto es más rápido que usar la receta de aumento de filas , es decir, hacer

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

Si resulta que necesita restricciones de no negatividad en sus coeficientes ajustados, entonces simplemente puede hacer

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 luego da un resultado un poco más preciso por cierto que

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 
}

(y estrictamente hablando, solo la solución nnls(A=crossprod(X)+diag(lambdas), b=crossprod(X,Y))$x es la correcta).

Todavía no he descubierto cómo el caso restringido de no negatividad podría optimizarse aún más para el p > ncaso; avíseme si alguien sabe cómo hacer esto ... [ lmridge_nnls_largep = function (X, Y, lambda) t(X) %*% nnls(A=tcrossprod(X)+lambda*diag(nrow(X)), b=Y)$xno funciona]