Uso del error estándar de la distribución bootstrap

Question

(ignora el código R si es necesario, ya que mi pregunta principal es independiente del lenguaje)

Si quiero ver la variabilidad de una estadística simple (por ejemplo, la media), sé que puedo hacerlo a través de la teoría como:

x = rnorm(50)

# Estimate standard error from theory
summary(lm(x~1))
# same as...
sd(x) / sqrt(length(x))

o con el bootstrap como:

library(boot)

# Estimate standard error from bootstrap
(x.bs = boot(x, function(x, inds) mean(x[inds]), 1000))
# which is simply the standard *deviation* of the bootstrap distribution...
sd(x.bs$t)

---

Sin embargo, lo que me pregunto es, ¿puede ser útil/válido(?) recurrir a la norma error de una distribución bootstrap en determinadas situaciones? La situación con la que estoy tratando es una función no lineal relativamente ruidosa, como:

# Simulate dataset
set.seed(12345)
n   = 100
x   = runif(n, 0, 20)
y   = SSasymp(x, 5, 1, -1) + rnorm(n, sd=2)
dat = data.frame(x, y)

Aquí el modelo ni siquiera converge utilizando el conjunto de datos original,

> (fit = nls(y ~ SSasymp(x, Asym, R0, lrc), dat))
Error in numericDeriv(form[[3L]], names(ind), env) : 
  Missing value or an infinity produced when evaluating the model

por lo que las estadísticas que me interesan en cambio son más estabilizado estimaciones de estos parámetros nls - quizás sus medias a través de un número de réplicas bootstrap.

# Obtain mean bootstrap nls parameter estimates
fit.bs = boot(dat, function(dat, inds)
              tryCatch(coef(nls(y ~ SSasymp(x, Asym, R0, lrc), dat[inds, ])),
                       error=function(e) c(NA, NA, NA)), 100)
pars = colMeans(fit.bs$t, na.rm=T)

En este caso, los resultados se aproximan a los que utilicé para simular los datos originales:

> pars
[1]  5.606190  1.859591 -1.390816

Una versión ploteada tiene este aspecto:

# Plot
with(dat, plot(x, y))

newx = seq(min(x), max(x), len=100)
lines(newx, SSasymp(newx, pars[1], pars[2], pars[3]))

lines(newx, SSasymp(newx, 5, 1, -1), col='red')
legend('bottomright', c('Actual', 'Predicted'), bty='n', lty=1, col=2:1)

Ahora, si quiero la variabilidad de estos estabilizado estimaciones de los parámetros, creo que puedo, asumiendo la normalidad de esta distribución bootstrap, simplemente calcular sus errores estándar:

> apply(fit.bs$t, 2, function(x) sd(x, na.rm=T) / sqrt(length(na.omit(x))))
[1] 0.08369921 0.17230957 0.08386824

¿Es éste un planteamiento sensato? ¿Existe un método general mejor para inferir los parámetros de modelos no lineales inestables como éste? (Supongo que aquí podría hacer una segunda capa de remuestreo, en lugar de basarme en la teoría para la última parte, pero eso podría llevar mucho tiempo dependiendo del modelo. Aun así, no estoy seguro de si estos errores estándar serían útiles para algo, ya que se acercarían a 0 si simplemente aumento el número de réplicas bootstrap).

Muchas gracias y, por cierto, soy ingeniero, así que disculpen que sea relativamente novato por aquí.

Answer 1

Esta pregunta plantea varios problemas. En primer lugar, está la cuestión de si los promedios bootstrap serán estimadores razonables incluso cuando algunos de los estimadores bootstrap individuales no sean computables (falta de convergencia, inexistencia de soluciones). En segundo lugar, dado que los estimadores bootstrapped son sensibles, se plantea la cuestión de cómo obtener intervalos de confianza o quizás sólo errores estándar para estas estimaciones.

La idea de promediar estimaciones bootstrap está estrechamente relacionada, si no es que es la misma, con la agregación bootstrap, o bagging, utilizada en el aprendizaje automático para mejorar el rendimiento de predicción de predictores débiles. Véase ESL Sección 8.7. En determinados casos $-$ también para estimar parámetros $-$ el promedio de las estimaciones bootstrap puede reducir la varianza del estimador resultante en comparación con el simple uso del estimador en el conjunto de datos original.

El objetivo de la pregunta es, sin embargo, producir estimaciones incluso en casos en los que el algoritmo para calcular las estimaciones pueda fallar ocasionalmente o en los que el estimador sea ocasionalmente indefinido. Como planteamiento general existe un problema:

• Si se promedian las estimaciones bootstrapped y se descartan ciegamente las muestras bootstrapped para las que las estimaciones no son computables, en general se obtendrán resultados sesgados.

La gravedad del problema general depende de varios factores. Por ejemplo, la frecuencia con la que la estimación no es computable y si la distribución condicional de la muestra dado que la estimación es no computable difiere de la distribución condicional de la muestra dado que la estimación es computable. Yo no recomendaría utilizar el método.

Para la segunda parte de la pregunta necesitamos un poco de notación. Si $X$ denota nuestro conjunto de datos original, $\hat{\theta}$ nuestro estimador (supongamos por simplicidad que tiene valor real y se le permite tomar el valor NA) tal que $\hat{\theta}(X)$ es la estimación para el conjunto de datos original, y $Y$ denota una única muestra bootstrap, entonces el promedio bootstrap está calculando efectivamente el estimador $$\tilde{\theta}(X) = E(\hat{\theta}(Y) \mid X, A(X))$$ donde $A(X)$ denota el acontecimiento, en función de $X$ sobre la cual $\hat{\theta}(Y) \neq \text{NA}$ . Es decir, calculamos la expectativa condicional del estimador en una muestra bootstrapada $-$ condicionamiento de la muestra original $X$ y el evento, $A(X)$ que el estimador es computable para la muestra bootstrap. El cálculo bootstrap real es una aproximación basada en el muestreo de $\tilde{\theta}(X)$ .

La sugerencia de la pregunta es calcular la desviación típica empírica de los estimadores bootstrap, que es una estimación de la desviación típica de $\hat{\theta}(Y)$ condicionado a $X$ y $A(X)$ . La desviación típica deseada, el error típico, es la desviación típica de $\tilde{\theta}(X)$ . No se puede obtener lo segundo de lo primero. No veo otra forma obvia y general que utilizar una segunda capa de bootstrapping para obtener una estimación fiable del error estándar.

La discusión sobre la estimación del error típico es independiente de la forma en que el condicionamiento sobre $A(X)$ afecta al sesgo del estimador $\tilde{\theta}(X)$ . Si el efecto es grave, incluso con estimaciones correctas del error estándar, un intervalo de confianza será engañoso.

Editar :

El papel muy bonito Estimación y precisión tras la selección del modelo de Efron ofrece un método general para estimar el error estándar de un estimador embolsado sin utilizar una segunda capa de bootstrapping. El artículo no trata explícitamente de los estimadores que en ocasiones no son computables.