Errores estándar para predicción de lazo usando R

Question

Estoy tratando de utilizar un modelo de lazo para la predicción, y necesito para estimar errores estándar. Seguramente alguien ya ha escrito un paquete para hacer esto. Pero en cuanto puedo ver, ninguno de los paquetes en CRAN que hacer predicciones mediante un lazo devolverá errores estándar para las predicciones.

Mi pregunta es: ¿hay un paquete o algún código de R para calcular errores estándar para las predicciones de lazo?

Answer 1

Kyung et al. (2010), "Penalizado de regresión, los errores estándar, y Bayesiano lazo", Análisis Bayesiano , 5, 2, sugieren que no hay un consenso sobre estadísticamente la validez del método de cálculo de los errores estándar para el lazo de predicciones. Tibshirani parece estar de acuerdo (lámina 43) que los errores estándar son todavía un problema no resuelto.

Answer 2

En una nota relacionada, que puede ser de ayuda, Tibshirani y sus colegas han propuesto una significación de la prueba para el lazo. El (todavía no publicada) manuscrito está disponible, y que se titulaba "la importancia de la prueba para el lazo".^†

† Ahora publicado: Lockhart(2014), "la significación de la prueba para el lazo", Ann. Estatismo., 42, 2

Answer 3

Bayesiano LAZO es la única alternativa para el problema del cálculo de los errores estándar. Los errores estándar se calculan automáticamente en Bayesiano LAZO...puede implementar Bayesiano LAZO muy fácilmente mediante Muestreo de Gibbs esquema de...

Bayesiano LAZO necesidades antes de que las distribuciones de ser asignados a los parámetros del modelo. En LAZO modelo, tenemos el objetivo de la función $||\mathbf{y}-\mathbf{X}\boldsymbol{\beta}||_2^2 + \lambda||\boldsymbol{\beta}||_1$ con $\lambda$ como el parámetro de regularización. Aquí tenemos $\ell_1$-norma para $\boldsymbol{\beta}$ entonces, un tipo especial de estado de la distribución que se necesita para esto, LAPLACE de distribución de una escala mezcla de una distribución normal con distribución exponencial como la mezcla de densidad. Basado en el completo condicional posteriores de cada uno de los parámetros que se puedan deducir.

Entonces uno puede utilizar Muestreo de Gibbs para la simulación de la cadena. Ver Park & Cassella (2008), "El Bayesiano Lazo", JASA, 103, 482.

Hay tres desventajas inherentes de frecuentista LAZO:

1. Uno tiene que elegir $\lambda$ por validación cruzada o por otros medios.

2. Los errores estándar son difíciles de calcular como la de LARS y otros algoritmos de producir estimaciones puntuales para $\boldsymbol{\beta}$.

3. La estructura jerárquica del problema en cuestión no puede ser codificada con frecuentista modelo, que es bastante fácil en el marco Bayesiano.

Answer 4

Sandipan Karmakar respuesta te dice qué hacer, esto debe ayudarle en el "cómo":

> library(monomvn)
>
> ## following the lars diabetes example
> data(diabetes)
> str(diabetes)
'data.frame':   442 obs. of  3 variables:
 $ x : AsIs [1:442, 1:10] 0.038075.... -0.00188.... 0.085298.... -0.08906.... 0.005383.... ...
      ..- attr(*, "dimnames")=List of 2
      .. ..$ : NULL
  .. ..$ : chr  "age" "sex" "bmi" "map" ...

 $ y : num  151 75 141 206 135 97 138 63 110 310 ...

[...]

> ## Bayesian Lasso regression
> reg_blas <- with(diabetes, blasso(x, y))
t=100, m=8
t=200, m=5
t=300, m=8
t=400, m=8
t=500, m=7
t=600, m=8
t=700, m=8
t=800, m=8
t=900, m=5
> 
> ## posterior mean beta (setting those with >50% mass at zero to exactly zero)
> (beta <- colMeans(reg_blas$beta) * (colMeans(reg_blas$beta != 0)  > 0.5))
      b.1       b.2       b.3       b.4       b.5       b.6       b.7       b.8 
   0.0000 -195.9795  532.7136  309.1673 -101.1288    0.0000 -196.4315    0.0000 
      b.9      b.10 
 505.4726    0.0000 
> 
> ## n x nsims matrix of realizations from the posterior predictive:
> post_pred_y <- with(reg_blas, X %*% t(beta))
> 
> ## predictions:
> y_pred <- rowMeans(post_pred_y)
> head(y_pred)
[1]  52.772443 -78.690610  24.234753   9.717777 -23.360369 -45.477199
> 
> ## sd of y:
> sd_y <- apply(post_pred_y, 1, sd)
> head(sd_y)
[1] 6.331673 6.756569 6.031290 5.236101 5.657265 6.150473
> 
> ## 90% credible intervals
> ci_y <- t(apply(post_pred_y, 1, quantile, probs=c(0.05, 0.95)))
> head(ci_y)
             5%       95%
[1,]  42.842535  62.56743
[2,] -88.877760 -68.47159
[3,]  14.933617  33.85679
[4,]   1.297094  18.01523
[5,] -32.709132 -14.13260
[6,] -55.533807 -35.77809

Answer 5

Para añadir a las respuestas anteriores, el problema parece ser que incluso un bootstrap es probablemente insuficiente como la estimación de la penalizado modelo es sesgada y arranque de sólo hablar de la varianza - ignorando el sesgo de la estimación. Esto está muy bien resumida en la viñeta para el sancionado el paquete en la Página 18.

Si se utiliza para la predicción, sin embargo, ¿por qué es un error estándar del modelo requiere? No se puede cruzar validar o bootstrap adecuadamente y producir un error estándar alrededor de una métrica relativa a la predicción como MSE?