Valor p de computación con sistema de arranque de R

Question

Yo uso paquete "arranque" para calcular un aproximado de 2 lados hizo valor p pero el resultado es demasiado lejos del valor de p de usar t.test. No puedo averiguar de lo que hice mal en mi código de R. Alguien por favor me puede dar un Consejo para esto

time = c(14,18,11,13,18,17,21,9,16,17,14,15,
         12,12,14,13,6,18,14,16,10,7,15,10)
group=c(rep(1:2, each=12))
sleep = data.frame(time, group)

require(boot)
diff = function(d1,i){
    d = d1[i,]
    Mean= tapply(X=d$time, INDEX=d$group, mean)
    Diff = Mean[1]-Mean[2]
    Diff
}

set.seed(1234)
b3 = boot(data = sleep, statistic = diff, R = 5000, strata=sleep$group)

pvalue = mean(abs(b3$t) > abs(b3$t0))
pvalue 

2-sided hizo p-valor (pvalue) = 0.4804 pero los 2 lados p-valor de t.test es 0.04342. Ambos valores de p son alrededor de 11 veces diferencia. ¿Cómo puede esto suceder? Muchas gracias de antemano!

Answer 1

La respuesta de Elvis se basa en permutaciones, pero en mi opinión no tiene claro lo que está mal con el original enfoque bootstrap. Permítanme hablar de una solución basada exclusivamente en bootstrap.

El problema crucial de su original de la simulación es que bootstrap siempre proporciona la VERDADERA distribución de la estadística de prueba. Sin embargo, cuando se calcula el p-valor que tiene que comparar el valor obtenido de la estadística de prueba para su distribución BAJO H0, es decir, no con la verdadera distribución!

[Vamos a dejar claro. Por ejemplo, se sabe que la prueba estadística de T de la clásica prueba de t tiene la clásica "central" distribución t de bajo H0 y un noncentral de distribución en general. Sin embargo, todo el mundo está familiarizado con el hecho de que el valor observado de T en comparación con el clásico "central" distribución t, es decir, no se intenta obtener el verdadero [noncenral] t-distribución de hacer la comparación con T.]

Su p-valor 0.4804 es tan grande, debido a que el valor observado "t0" de la prueba estadística de Media[1] Media[2] se encuentra muy cerca del centro de la bootstrapped muestra de la "t". Es natural y normalmente siempre es así [es decir, independientemente de la validez de H0], porque el bootstrap de la muestra "t" emula el de la distribución REAL de la Media[1] Media[2]. Pero, como se señaló anteriormente [y también por Elvis], lo que realmente se necesita es la distribución de la Media[1] Media[2] BAJO H0. Es obvio que

1) bajo H0 la distribución de la Media[1] Media[2] estará centrado alrededor de 0,

2) su forma no depende de la validez de H0.

Estos dos puntos implica que la distribución de la Media[1] Media[2] bajo H0 puede ser emulado por el bootstrap de la muestra "t" DESPLAZADO de modo que se centra alrededor de 0. En R:

b3.under.H0 <- b3$t - mean(b3$t)

y el correspondiente p-valor será:

mean(abs(b3.under.H0) > abs(b3$t0))

lo que le da una "muy buena" valor de 0.0232. :-)

Permítanme señalar que el punto "2)" que se menciona anteriormente se llama "traducción equivariance" de la estadística de prueba y NO tiene que aguantar en general! I. e. para algunos de la estadística de prueba, el cambio de la bootstrap "t" no te proporciona una estimación válida de la distribución de la estadística de prueba de bajo HO! Echa un vistazo a este debate y, especialmente, en la respuesta de P. Dalgaard: http://tolstoy.newcastle.edu.au/R/e6/help/09/04/11096.html

Las pruebas problema no rinde perfectamente simétrica la distribución de la estadística de prueba, pero hay que tener en cuenta que hay algunos problemas con la obtención de DOS CARAS, los valores de p en el caso de sesgada bootstrap distribución de la estadística de prueba. De nuevo, lea el enlace de arriba.

[Y por último, me gustaría utilizar el "puro" prueba de permutación en su situación; es decir, la segunda mitad de Elvis respuesta. :-)]

HTH :-)

jan.

Answer 2

Usted está usando bootstrap para generar los datos en virtud de la distribución empírica de los datos observados. Esto puede ser útil para dar un intervalo de confianza de la diferencia entre los dos medios:

> quantile(b3$t,c(0.025,0.975))
     2.5%     97.5% 
0.4166667 5.5833333 

Para obtener un $p$-valor, es necesario generar permutaciones bajo la hipótesis nula. Esto se puede hacer por ejemplo como este:

diff2 = function(d1,i){
    d = d1; 
    d$group <- d$group[i];  # randomly re-assign groups
    Mean= tapply(X=d$time, INDEX=d$group, mean)
    Diff = Mean[1]-Mean[2]
    Diff
}

> set.seed(1234)
> b4 = boot(data = sleep, statistic = diff2, R = 5000)
> mean(abs(b4$t) > abs(b4$t0))
[1] 0.046

En esta solución, el tamaño de los grupos no es fijo, de forma aleatoria a asignar un grupo a cada individuo por bootstraping a partir de la inicial de los grupos. Parece de fiar para mí, sin embargo, una más clásica solución es fijar el número de individuos de cada grupo, por lo que sólo permutar los grupos en lugar de bootstraping (esto es generalmente motivados por el diseño del experimento, donde los tamaños de los grupos son fijadas de antemano):

> R <- 10000; d <- sleep
> b5 <- numeric(R); for(i in 1:R) { 
+    d$group <- sample(d$group, length(d$group)); 
+    b5[i] <- mean(d$time[d$group==1])-mean(d$time[d$group==2]); 
+ }
> mean(abs(b5) > 3)
[1] 0.0372