La terminología es, probablemente, no se utilizan de forma sistemática, por lo que la siguiente es sólo la forma en la que entiendo la pregunta original. Desde mi entender, el normal CIs calcula no son lo que se le pidió. Cada conjunto de bootstrap replica le da un intervalo de confianza, no muchos. La manera de calcular las diferentes CI-tipos a partir de los resultados de un conjunto de bootstrap replica es la siguiente:
B <- 999 # number of replicates
muH0 <- 100 # for generating data: true mean
sdH0 <- 40 # for generating data: true sd
N <- 200 # sample size
DV <- rnorm(N, muH0, sdH0) # simulated data: original sample
Ya quiero comparar los cálculos de los resultados de paquete boot, la primera vez que definir una función que será llamada para cada repetición. Sus argumentos son la muestra original, y un índice del vector de la especificación de los casos para una sola repetición. Devuelve $M^{\star}$, el plug-in de estimación de $\mu$, así como el $S_{M}^{2\star}$, el plug-in de estimación de la varianza de la media de $\sigma_{M}^{2}$. El último será necesaria únicamente para el bootstrap $t$-CI.
> getM <- function(orgDV, idx) {
+ bsM <- mean(orgDV[idx]) # M*
+ bsS2M <- (((N-1) / N) * var(orgDV[idx])) / N # S^2*(M)
+ c(bsM, bsS2M)
+ }
> library(boot) # for boot(), boot.ci()
> bOut <- boot(DV, statistic=getM, R=B)
> boot.ci(bOut, conf=0.95, type=c("basic", "perc", "norm", "stud"))
BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 999 bootstrap replicates
CALL :
boot.ci(boot.out = bOut, conf = 0.95, type = c("basic", "perc", "norm", "stud"))
Intervals :
Level Normal Basic Studentized Percentile
95% ( 95.6, 106.0 ) ( 95.7, 106.2 ) ( 95.4, 106.2 ) ( 95.4, 106.0 )
Calculations and Intervals on Original Scale
Sin usar el paquete boot usted puede simplemente utilizar replicate() para obtener un conjunto de bootstrap repeticiones.
boots <- t(replicate(B, getM(DV, sample(seq(along=DV), replace=TRUE))))
Pero sigamos con los resultados de boot.ci() para tener una referencia.
boots <- bOut$t # estimates from all replicates
M <- mean(DV) # M from original sample
S2M <- (((N-1)/N) * var(DV)) / N # S^2(M) from original sample
Mstar <- boots[ , 1] # M* for each replicate
S2Mstar <- boots[ , 2] # S^2*(M) for each replicate
biasM <- mean(Mstar) - M # bias of estimator M
La básica, percentil, y $t$-CI dependen de la distribución empírica de bootstrap estimaciones. Para obtener el $\alpha/2$ $1 - \alpha/2$ cuantiles, nos encontramos con los índices correspondientes a los ordenados vector de estimaciones bootstrap (tenga en cuenta que boot.ci() va a hacer más complicado de interpolación para encontrar los cuantiles empíricos cuando los índices no son números naturales).
(idx <- trunc((B + 1) * c(0.05/2, 1 - 0.05/2)) # indices for sorted vector of estimates
[1] 25 975
> (ciBasic <- 2*M - sort(Mstar)[idx]) # basic CI
[1] 106.21826 95.65911
> (ciPerc <- sort(Mstar)[idx]) # percentile CI
[1] 95.42188 105.98103
Para el $t$-CI, necesitamos el bootstrap $t^{\star}$ estimaciones para calcular la crítica $t$-valores. Para el normal estándar de CI, el valor crítico será sólo el $z$-valor de la distribución normal estándar.
# standard normal CI with bias correction
> zCrit <- qnorm(c(0.025, 0.975)) # z-quantiles from std-normal distribution
> (ciNorm <- M - biasM + zCrit * sqrt(var(Mstar)))
[1] 95.5566 106.0043
> tStar <- (Mstar-M) / sqrt(S2Mstar) # t*
> tCrit <- sort(tStar)[idx] # t-quantiles from empirical t* distribution
> (ciT <- M - tCrit * sqrt(S2M)) # studentized t-CI
[1] 106.20690 95.44878
Para la estimación de la cobertura de las probabilidades de estas CI-tipos, usted tendrá que ejecutar esta simulación muchas veces. Simplemente envuelva el código en una función, devuelve una lista con los CI-resultados y ejecutarlo con replicate() como se demostró en este gist.
