Procesos de Remuestreo

Descargar como pdf o txt
Descargar como pdf o txt
Está en la página 1de 73

Métodos de Remuestreo

Tema 3. Estimación de errores estándar mediante


remuestreo

basado en
B. Efron, R. Tibshirani (1993). An Introduction to the bootstrap.
O. Kirchkamp (2017). Resampling methods.

Curso 2018/19

1/73
Errores estándar

I Los estadı́sticos muestrales se usan frecuentemente en Estadı́stica,


por ello es necesario conocer su precisión.

I El bootstrap permite encontrar el error estándar de los estadı́sticos


basándose en el principio de plug-in.

I Supongamos una v.a. X con media µF y varianza σF2

I Sea x = (x1 , x2 , . . . , xn ) una m.a.s procedente de la distribución F .

I La media muestral
σF2
 
x̄ ∼ µF ,
n

2/73
Cálculo de errores estándar mediante el TCL

I De este modo la esperanza de x̄ es la misma que la de la v.a. X


original, pero la varianza es igual a 1/n veces la varianza de X .

σ
I Ası́ el error estándar de x̄ es simplemente √F .
n

I En una distribución normal se espera que X sea menor que una vez
la desviación estándar de µF aproximadamente el 68 % de las
ocasiones y menor que dos desviaciones estándar alrededor del 95 %
de las veces, aplicando el TCL (teorema central del lı́mite).

3/73
4/73
Limitaciones del TCL

I La aproximación del TCL funciona bien cuando el tamaño muestral


n es grande pero la aproximación tiene limitaciones.

I Supongamos que X sigue una distribución de Bernoulli:

PF {X = 1} = p
PF {X = 0} = 1 − p

I El parámetro p es la probabilidad de éxito que está entre 0 y 1

5/73
Limitaciones del TCL

I Una m.a.s. es una sucesión de unos y ceros de modo que la suma

n
X
s= xi ∼ Bin(n, p)
i=1

s
I La media x = es igual a p
b que es el estimador plug-in de p, de
n
modo que
 
p(1 − p)
b ∼ N p,
p
n

6/73
Limitaciones del TCL

I Se toma el ejemplo de una distribución binomial con n = 25 en los


casos de p = 0,25 y p = 0,9.

I Para el caso de p = 0,9 la aproximación a la normal por el TCL no


es muy buena.

n = 25

p09 = rbinom (20000 , n , 0 .9 ) / n


p025 = rbinom (20000 , n , 0 .25 ) / n

7/73
Limitaciones del TCL

par ( mfrow = c (1 ,2))

hist ( p09 , prob =T , xlim = c (0 .6 ,1 .1 ) , col = " lightblue " ,


xlab = " p " , main = " p = 0 .9 " )
xs1 = seq (0 .6 ,1 .1 ,0 .0001 )
ys1 = dnorm ( xs1 , 0 .9 , sqrt ((0 .9 * 0 .1 ) / n ))
lines ( xs1 , ys1 , lwd =3 , col = " red " )

hist ( p025 , prob =T , xlim = c (0 ,0 .6 ) , col = " pink " ,


xlab = " p " , main = " p = 0 .25 " )
xs2 = seq (0 ,0 .6 ,0 .0001 )
ys2 = dnorm ( xs2 , 0 .25 , sqrt ((0 .25 * 0 .75 ) / n ))
lines ( xs2 , ys2 , lwd =3 , col = " red " )

8/73
9/73
Bootstrap y errores estándar

I El bootstrap permite calcular errores estándar sin que tenga


importancia lo complicado que sea el estimador que se considere.
I Los métodos bootstrap dependen del concepto de muestra
bootstrap.
I Partimos de la función de distribución empı́rica Fb que asigna
probabilidad 1/n a cada uno de los elementos de la muestra
observada.
I Una muestra bootstrap se define como una muestra aleatoria de
tamaño n extraı́da de Fb

Fb → x∗ = (x1∗ , x2∗ , . . . , xn∗ )

I La notación estrella * indica que x∗ no es el conjunto de datos


original sino una versión remuestreada de la muestra original x.

10/73
Bootstrap y errores estándar
I Alternativamente, se puede decir que una muestra bootstrap
x1∗ , x2∗ , . . . xn∗ es una muestra aleatoria de tamaño n tomada con
reemplazamiento de la muestra original (que hace el papel de
población).

I El algoritmo se denomina bootstrap no paramétrico porque


depende solo de la función de distribución empı́rica.

I Por ejemplo podrı́amos tener una muestra bootstrap como

x1∗ = x7
x2∗ = x3
x3∗ = x3
··· ··· ···
xn∗ = x2

11/73
Algoritmo Bootstrap y errores estándar

1. Seleccionar B muestras bootstrap x∗1 , x∗2 , . . . , x∗B cada una


obtenida a partir de la muestra original x con reemplazamiento.
2. Evaluar la réplica bootstrap en el estimador correspondiente

θb∗ (b) = s(x∗b )

para b = 1, 2, . . . , B.
3. Estimar el error estándar seF (θ)
b mediante
v
u
u 1 X B  2
se
bB = t θb∗ (b) − θb∗ (·)
B−1
b=1

B
1 X b∗
donde θb∗ (·) = θ (b)
B
b=1

12/73
13/73
Ejemplo de los institutos de máster en leyes
I La correlación entre GPA y LSAT es

library ( bootstrap )
( lawCor = with ( law , cor ( GPA , LSAT )))

[1] 0 .7763745

I Pero ¿cómo es de preciso el estimador del coeficiente de correlación


lineal?

I Si la distribución conjunta de ambas variables F es normal


bivariante, entonces ρ̂ (siguiendo a Efron&Tibshirani) tiene un error
estándar igual a
1 − ρb2
σ
bbρ
= √ ≈ 0,115
n−3

14/73
Ejemplo de los institutos de máster en leyes
I Ası́,
( se = (1 - lawCor ˆ2) / sqrt ( dim ( law )[1] -3))

[1] 0 .1146741

c (( lawCor - 1 .96 * se ) , min (1 , ( lawCor + 1 .96 * se )))

[1] 0 .5516133 1 .00000

I Alternativamente se puede usar la librerı́a psychometric:

psychometric :: CIr ( lawCor , dim ( law )[1])

[1] 0 .4385108 0 .9219648

15/73
Ejemplo de los institutos de máster en leyes

I Usando bootstrap se puede evitar asumir que F se distribuye como


una normal bivariante.

ssamplesize = dim ( law )[1]


ind = 1: samplesize

law.boot =
replicate (1000 , { indB = sample ( ind , replace = TRUE );
with ( law [ indB ,] , cor ( GPA , LSAT ))})

sd ( law.boot )

[1] 0 .1336493

16/73
Error estándar del coeficiente de correlación

I Desde el punto de vista clásico, el error estándar del estimador de la


correlación, se puede estimar ası́:

# Asumes SE _ r = sqrt ((1 - r ˆ2) / (n -2))

cor.test.plus = function ( x ) {
list (x ,
Stand ard.Erro r =
unname ( sqrt ((1 - x $ estimate ˆ2) / x $ parameter )))
}

17/73
Error estándar del coeficiente de correlación

library ( bootstrap )
cor.test.plus ( cor.test ( law $ GPA , law $ LSAT ))

Pearson ’s product - moment correlation

data : law $ GPA and law $ LSAT


t = 4 .4413 , df = 13 , p - value = 0 .0006651
alternative hypothesis : true correlation is not equal to 0

95 percent confidence interval :


0 .4385108 0 .9219648

sample estimates :
cor
0 .7763745

Stand ard.Error
[1] 0 .174806

18/73
Ejemplo de los institutos de máster en leyes

I ¿Cómo converge de rápido el estimador bootstrap?

samplesize = dim ( law )[1]


ind = 1: samplesize
lawBS = function ( B ) sd ( replicate (B ,
{ indB = sample ( ind , replace = TRUE );
with ( law [ indB ,] , cor ( GPA , LSAT ))}))

BStamannos = seq (200 ,5000 ,200)


BSestimas = sapply ( BStamannos , lawBS )

library ( ggplot2 )
qplot ( BStamannos , BSestimas , geom = " path " )

19/73
Ejemplo de los institutos de máster en leyes

I Se compara la distribución empı́rica muestral de la población Fb (θb∗ )


con la distribución poblacional F (θb∗ )

ind1 = 1: dim ( law )[1]


law.boot = replicate (5000 , { indB = sample ( ind1 ,
size = length ( ind1 ) , replace = TRUE );
with ( law [ indB ,] , cor ( GPA , LSAT ))})

ind2 = 1: dim ( law82 )[1]


law82.boot = replicate (5000 ,{ indB = sample ( ind2 ,
size = length ( ind2 ) , replace = TRUE );
with ( law82 [ indB ,] , cor ( GPA , LSAT ))})

21/73
Ejemplo de los institutos de máster en leyes

library ( latticeExtra )

densityplot (∼law82.boot + law.boot , plot.points = FALSE ,


auto.key = list ( columns =2 , size =3 , between =1 ,
col = c ( " red " ," blue " )) ,
par.settings = list ( superpo se.line =
list ( col = c ( " red " ," blue " ))) , col = c ( " red " ," blue " ))

22/73
23/73
Ejemplo de los institutos de máster en leyes

I Se compara el estimador bootstrap del error estándar a partir de la


muestra sebF
(b
ρ) con respecto al error estándar a partir de la
población seF (b
ρ)

sd ( law.boot )

[1] 0 .1327579

sd ( law82.boot )

[1] 0 .1309014

24/73
Bootstrap Paramétrico

I En muchas ocasiones se tienen fórmulas analı́ticas para calcular los


errores. En este caso se puede aplicar el bootstrap aprovechando que
éstas se conocen.

I Se denomina a este tipo de remuestreo como bootstrap


paramétrico y se define el estimador bootstrap del error estándar
como
seb
F
(θb∗ )
par

I donde Fbpar es un estimador de F que se obtiene a partir de un


modelo paramétrico aplicado a los datos.

25/73
Bootstrap Paramétrico

I En el ejemplo law82, en lugar de estimar la función de distribución


F mediante la función de distribución empı́rica, se puede asumir que
la población se distribuye como una normal bivariante.

I Para la media y la matriz de covarianzas de esta distribución, los


estimadores razonables serı́an respectivamente

(x̄ , ȳ )

!
(xi − x )2
P P
1 i (xi − x )(yi − y )
P i P 2
14 i (xi − x )(yi − y ) i (yi − y )

I Se denota a la población normal bivariante que se obtiene con esta


media y matriz de covarianzas como Fbpar .

26/73
Bootstrap Paramétrico

I Se denomina al estimador bootstrap paramétrico del error estándar


del parámetro como seF̂par (θb∗ ).

I En lugar de muestrear con reemplazamiento a partir de los datos


originales, se sacan B muestras de tamaño n del estimador
paramétrico de la población Fbpar

Fbpar → (x1∗ , x2∗ , . . . , xn∗ )

I Posteriormente se siguen los mismos pasos 2 y 3 del algoritmo


general del bootstrap no paramétrico: se calcula el correspondiente
estadı́stico en cada muestra bootstrap y luego se calcula la
desviación estándar de las B réplicas.

27/73
Ejemplo de los centros de estudios de máster

I En el ejemplo de los datos de los centros de estudios de máster en


leyes se tiene que si (x , y ) se distribuyen como una normal bivariante
entonces se pueden generar observaciones de este vector, definiendo

x = µx + σx z1
z1 + c · z2
y = µy + σy √
1 + c2
donde z1 , z2 ∼ N(0, 1)
s
σx2 σy2
c= 2
−1
σxy

28/73
Ejemplo de los institutos de máster en leyes

library ( bootstrap )

paraBoot = function ( datos ) {


ndatos = dim ( datos )[1]
sigma = cov ( datos )
mu = rapply ( datos , mean )

c = sqrt ( prod ( diag ( sigma )) / sigma [1 ,2]ˆ2 -1)


z1 = rnorm ( ndatos )
z2 = rnorm ( ndatos )
x = mu [1] + sqrt ( sigma [1 ,1]) * z1
y = mu [2] + sqrt ( sigma [2 ,2]) * ( z1 + c * z2 ) / sqrt (1+ c ˆ2)
cbind (x , y )
}

29/73
Ejemplo de los institutos de máster en leyes

Alternativamente, se puede programar con una librerı́a especı́fica que


trata la normal multivariante: mvrnorm.

paraBoot = function ( datos ) {


ndatos = dim ( datos )[1]
sigma = cov ( datos )
mu = rapply ( datos , mean )

x = MASS :: mvrnorm ( ndatos , mu = mu , Sigma = sigma )


return ( x )
}

30/73
Ejemplo de los institutos de máster en leyes

# Bootstrap parametrico
pBoot = replicate (5000 , cor ( paraBoot ( law ))[2 ,1])

# Aproximacion bootstrap
sd ( pBoot )

[1] 0 .1214904

samplesize = dim ( law )[1]


# Aproximacion asintotica
(1 - cor ( law )[2 ,1]ˆ2) / sqrt ( samplesize -3)

[1] 0 .1146741

31/73
Ejemplo de los institutos de máster en leyes
library ( latticeExtra )

ind1 = 1: dim ( law )[1]


law.boot = replicate (5000 , { indB = sample ( ind1 ,
size = dim ( law )[1] , replace = TRUE );
with ( law [ indB ,] , cor ( GPA , LSAT ))})

ind2 = 1: dim ( law82 )[1]


law82.boot = replicate (5000 ,{ indB = sample ( ind2 ,
size = dim ( law82 )[1] , replace = TRUE );
with ( law82 [ indB ,] , cor ( GPA , LSAT ))})

densityplot (∼pBoot + law82.boot + law.boot , plot.points = FALSE ,


auto.key = list ( columns =3 , size =2 , between =1 ,
col = c ( " orange " ," darkgreen " ," blue " )) ,
par.settings = list ( superpo se.line =
list ( col = c ( " orange " ," darkgreen " ," blue " ))) ,
col = c ( " orange " ," darkgreen " ," blue " ))

32/73
33/73
Bootstrap Paramétrico

I La mayor parte de los errores estándar son aproximaciones basadas


en la distribución normal.

I Estas aproximaciones se parecen a los resultados que se obtienen con


el bootstrap paramétrico cuando se hace remuestreo de la
distribución normal.

I Cuando se usa bootstrap paramétrico se obtienen resultados más


precisos que en las aproximaciones asintóticas cuando éstas existen.

34/73
Aplicación a datos multivariantes

I Ejemplo: Se tienen unos datos sobre calificaciones en 5 asignaturas


de 88 alumnos (ver el libro sobre Análisis Multivariante, de Mardia,
Kent and Bibby, (1979)):
I mec: mechanics
I vec: vectors
I alg: algebra
I ana: analysis
I sta: statistics

library ( bootstrap )
data ( scor )
plot ( scor )

35/73
36/73
Aplicación a datos multivariantes
I El vector de medias y la correspondiente matriz de covarianzas son:

colMeans ( scor )

mec vec alg ana sta


38 .95455 50 .59091 50 .60227 46 .68182 42 .30682

cov ( scor )

mec vec alg ana sta


mec 305 .7680 127 .22257 101 .57941 106 .27273 117 .40491
vec 127 .2226 172 .84222 85 .15726 94 .67294 99 .01202
alg 101 .5794 85 .15726 112 .88597 112 .11338 121 .87056
ana 106 .2727 94 .67294 112 .11338 220 .38036 155 .53553
sta 117 .4049 99 .01202 121 .87056 155 .53553 297 .75536

37/73
Aplicación a datos multivariantes
I Se calculan los autovalores y autovectores de la matriz de
covarianzas.
I La matriz 5 × 5 de covarianzas tiene 5 autovalores positivos en orden
b1 ≥ λ
decreciente: λ b2 ≥ λb3 ≥ λb4 ≥ λ
b5 y a cada uno de ellos le
corresponde un autovector diferente.
round ( eigen ( cov ( scor )) $ values ,3) # Autovalores

[1] 686 .990 202 .111 103 .747 84 .630 32 .153

round ( eigen ( cov ( scor )) $ vectors ,3) # Autovectores

[ ,1] [ ,2] [ ,3] [ ,4] [ ,5]


[1 ,] -0 .505 0 .749 -0 .300 0 .296 -0 .079
[2 ,] -0 .368 0 .207 0 .416 -0 .783 -0 .189
[3 ,] -0 .346 -0 .076 0 .145 -0 .003 0 .924
[4 ,] -0 .451 -0 .301 0 .597 0 .518 -0 .286
[5 ,] -0 .535 -0 .548 -0 .600 -0 .176 -0 .151

38/73
Aplicación a datos multivariantes

I Los autovalores y autovectores de la matriz de covarianzas son


importantes para explicar la estructura multivariante de los datos.

I Se observa que las calificaciones en los exámenes están altamente


correlacionados entre sı́: un estudiante con calificaciones altas en
mecánica suele tenerlas altas también en cálculo vectorial.

I Un modelo posible para las medidas correlacionadas serı́a

xi = Q i v

para i = 1, . . . , 88

39/73
Aplicación a datos multivariantes

I Donde Qi es un número que representa la capacidad del estudiante i


mientras que v es un vector de valores fijos para todos los
estudiantes

I Qi se puede interpretar como el coeficiente intelectual (IQ) del


estudiante i-ésimo.

I Si el modelo anterior fuese cierto, entonces solo el primer autovalor


λ
b1 serı́a positivo y el resto de autovalores serı́an igual a 0.

I También, v serı́a igual al primer autovector b


v1 .

40/73
Aplicación a datos multivariantes

I Se define el ratio del mayor autovalor con respecto al total θ,


b

λb1
θb = P5
i=1 λi
b

I Ası́ el modelo anterior es equivalente a θb = 1.

I Aunque, en la práctica, no se espera que sea exactamente igual a 1

I En el caso de las calificaciones, la estimación de θb es

686,990
θb = = 0,619
686,990 + 202,111 + 103,747 + 84,630 + 32,153

41/73
Aplicación a datos multivariantes

I En muchas circunstancias es interesante tener un valor alto de θb


porque eso indica un alto poder explicativo del modelo.

I El valor de θb mide el porcentaje de varianza explicada por el primer


componente.

I Cuanto más cerca estén los puntos respecto al eje del componente
principal, mayor será el valor de θ.
b

I ¿Qué precisión tiene θ?


b ¿Cuál es el error estándar de θ?
b

I Esta serı́a una aplicación directa del bootstrap en este caso.

42/73
Aplicación a datos multivariantes

I La complejidad del cálculo de θb no resulta relevante, en tanto que se


pueda calcular θb∗ para cualquier muestra bootstrap.

I En este caso, una muestra bootstrap es una matriz remuestreada X∗


de tamaño 88×5.

I Las filas xi∗ de X∗ proceden de una m.a.s de tamaño 88 de las filas


de la matriz de datos original.

x1∗ = xi1 , x2∗ = xi2 , . . . , x88



= xi88

I De este modo, algunas filas de X aparecerán varias veces y otras


ninguna en la matriz remuestreada X∗ .

43/73
Aplicación a datos multivariantes

I Una vez generada X∗ se calcula la matriz de covarianzas G∗ de la


manera habitual y luego se calculan los autovalores correspondientes.

I Se calcula la réplica bootstrap de θb

b∗
λ
θb∗ = P5 1
b∗
j=1 λj

I Y se aplica el algoritmo general bootstrap para calcular el error


estándar.

44/73
Aplicación a datos multivariantes

library ( bootstrap )

# En componente principales svd es numericamente mas


# estable que la descom posicion por autovectores y
# autovalores , pero para aplicar bootstrap
# esta utima es mas rapida

autovals = eigen ( var ( scor ) , symmetric = TRUE ,


only.values = TRUE ) $ values
( teta = autovals [1] / sum ( autovals ))

[1] 0 .619115

45/73
Aplicación a datos multivariantes

theta = function ( ind ){


vals = eigen ( var ( scor [ ind ,]) , symmetric = TRUE ,
only.values = TRUE ) $ values
vals [1] / sum ( vals )
}

scor.boot = bootstrap (1: dim ( scor )[1] , 500 , theta )


sd ( scor.boot $ thetastar ) # error estandar del bootstrap

[1] 0 .04570752

library ( ggplot2 )
qplot ( scor.boot $ thetastar , geom = " histogram " , binwidth =0 .02 ,
fill = I ( " lightgreen " ) , xlab = " Samples " ) +
geom _ vline ( xintercept = teta , col = " red " )

46/73
47/73
Aplicación a datos multivariantes

library ( bootstrap )
data ( scor )
X = scor

eigenTeta = function ( X ) {
ee = eigen ( cov ( X ))[[ " values " ]]
ee [1] / sum ( ee )
}

ind = 1: dim ( X )[[1]]

eigendist = replicate (5000 ,


eigenTeta ( X [ sample ( ind , replace = TRUE ) ,]))

library ( latticeExtra )
densityplot ( eigendist , plot.points = FALSE ,
xlab = expression ( theta ))

48/73
49/73
Aplicación a datos multivariantes

I El autovector b
v1 que corresponde al mayor autovalor se le denomina
primer componente principal de G
I Supongamos que se trata de de resumir el rendimiento de los
estudiantes mediante un solo número, en lugar de con 5 notas.
I Se puede demostrar que la mejor combinación lineal de las 5 notas es
5
X
yi = vb1k xik
k=1

es decir, una combinación lineal donde los componentes b


v1 equivalen
a los pesos de las notas originales.
I Esta combinación lineal es óptima en el sentido de que captura la
mayor parte de la variabilidad de las 5 puntuaciones originales de
entre todos los posibles v.

50/73
Aplicación a datos multivariantes

I La segunda combinación lineal


5
X
zi = vb2k xik
k=1

es el segundo componente principal b


v2 es decir, el segundo
autovector de G.

I El primer componente se puede asociar a la media de puntuaciones


de un estudiante, mientras que el segundo parece asociarse más bien
a la relación que hay entre exámenes con libro abierto frente a
cerrado.

51/73
Aplicación a datos multivariantes

I Tanto b
v1 como bv2 son estadı́sticos del mismo modo que lo es θ,
b y de
este modo se puede aplicar el bootstrap para calcular su variabilidad.

library ( bootstrap )
data ( scor )
X = scor

eigenVec = function ( X ) {
ee = eigen ( cov ( X ))[[ " vectors " ]]
return ( cbind ( ee [ ,1] , ee [ ,2]))
}

52/73
Aplicación a datos multivariantes

ind = 1: dim ( X )[[1]]


eigendist = replicate (500 ,
eigenVec ( X [ sample ( ind , replace = TRUE ) ,]))

apply ( eigendist [1:5 ,1 ,] , 1 , sd )

[1] 0 .2525616 0 .1828719 0 .1704597 0 .2257041 0 .2645642

apply ( eigendist [1:5 ,2 ,] , 1 , sd )

[1] 0 .50333116 0 .20130403 0 .07855624 0 .23908044 0 .41575549

boxplot ( cbind ( eigendist [1 ,1 ,] , eigendist [2 ,1 ,] ,


eigendist [3 ,1 ,] , eigendist [4 ,1 ,] , eigendist [5 ,1 ,]) ,
main = " Componente 1 " , col = " lightblue " )

53/73
54/73
Cuando puede fallar el bootstrap
I Consideramos el siguiente problema:
X se distribuye como una distribución uniforme en (0, θ).
El estimador de máxima verosimilitud para θ es el máx(Xi )

I Tenemos una muestra de 50 observaciones.

I Comparamos el estimador bootstrap no paramétrico de θ con


respecto al estimador paramétrico del mismo.

N = 50
X = runif ( N )
( thetaHat = max ( X ))

[1] 0 .990335

55/73
Cuando puede fallar el bootstrap

standardBoot = replicate (500 ,


max ( sample (X , N , replace = TRUE )))

paramBoot = replicate (500 ,


max ( runif (N , min =0 , max = thetaHat )))

library ( latticeExtra )

densityplot (∼paramBoot + standardBoot , xlab = " " ,


auto.key = list ( columns =2 , size =2 , between =1 ,
col = c ( " violet " , " orange " )) ,
par.settings = list ( superpo se.line =
list ( col = c ( " violet " , " orange " ))) ,
col = c ( " violet " , " orange " ))

56/73
57/73
Estructuras de datos generales

I Hasta ahora se ha considerado una estructura simple de los datos: el


modelo unimuestral donde una distribución de probabilidad
desconocida F genera los datos X mediante muestreo aleatorio.

I Pero algunos datos xi pueden ser bastante complejos, como


vectores, mapas o imágenes.

I Estructuras complejas de datos aparecen en modelos como series


temporales, análisis de varianza, modelos de regresión, datos
censurados o muestreo estratificado.

I Pero el método bootstrap se puede adaptar a estructuras de datos


generales.

58/73
Problemas unimuestrales

I El esquema del método bootstrap para problemas unimuestrales se


basa en la existencia de dos mundos paralelos.

I Por un lado está el mundo real con una distribución desconocida F


de la que se toma una muestra aleatoria y se calcula un estadı́stico a
partir de x digamos θb = s(x). Despúes se trata de estudiar su
comportamiento: errores, intervalos de confianza, etc.

I Por otro lado está el mundo bootstrap de modo que la población se


reduce a la muestra original y a partir de la distribución empı́rica Fb
se obtienen la muestras bootstrap x∗ .

I A partir de ella se calcula el estadı́stico de interés θb∗ = s(x∗ ) y se


estudia su comportamiento.

59/73
60/73
Problemas unimuestrales

I La doble flecha del esquema indica el cálculo de Fb a partir de F .

I Conceptualmente este es el paso fundamental del bootstrap y el


resto de pasos se definen por analogı́a.

I El procedimiento bootstrap para estructuras más complejas es


inmediato una vez que se sabe como realizar el proceso de la doble
flecha, es decir cómo estimar el mecanismo probabilı́stico a partir de
los datos.

I Se usa la notación P → x para indicar que un modelo de


probabilidad desconocido P ha generado el conjunto de datos x.

61/73
Problemas de dos muestras

I En el caso del problema de inferencia de dos muestras, el modelo de


probabilidad se puede considerar como P = (F , G) donde F es la
distribución de probabilidad del primer grupo y G la del segundo
grupo.

I Se obtienen dos muestras aleatorias independientes x = (z, y) de


modo que la aplicación P → x se decribe como
F → z e independientemente G → y

I En este caso se toman las respectivas funciones de distribución


empı́ricas y el estimador natural de P se construye como Pb = (Fb , G)
b
∗ ∗ ∗
y se obtiene una muestra bootstrap x = (z , y ) como
Fb → z∗ e independientemente G b → y∗

62/73
Ejemplo de los ratones
I Se toma el ejemplo de las diferencias de medias entre ratones según
son tratamiento o control

Trata = c (94 ,197 ,16 ,38 ,99 ,141 ,23)


Cont = c (52 ,104 ,146 ,10 ,51 ,30 ,40 ,27 ,46)

mean ( Trata ) - mean ( Cont )

[1] 30 .63492

B = 1000

sd ( replicate (B , mean ( sample ( Trata , replace = TRUE )) -


mean ( sample ( Cont , replace = TRUE ))))

[1] 26 .19047

63/73
Ejemplo de los ratones

library ( simpleboot )

b = two.boot ( Trata , Cont , mean , R = B )

sd ( b $ t )

[1] 27 .27626

hist (b , col = " lightblue " )

64/73
Ejemplo de los ratones

# Programa original de la libreria bootstrap

B = 1000
library ( bootstrap )
mouse.boot.c = bootstrap ( mouse.c , B , mean )
mouse.boot.t = bootstrap ( mouse.t , B , mean )

m ou se .b o ot .d if f =
mouse.boot.t $ thetastar - mouse.boot.c $ thetastar

sd ( m ou s e. bo ot . di ff )

[1] 26 .58504

66/73
Ejemplo de los ratones

Trata = c (94 ,197 ,16 ,38 ,99 ,141 ,23)


Cont = c (52 ,104 ,146 ,10 ,51 ,30 ,40 ,27 ,46)
B = 1000
n = length ( Trata )
Losratones = c ( Trata , Cont )

( t.obs = mean ( Trata ) - mean ( Cont ))

[1] 30 .63492

67/73
Ejemplo de los ratones

library ( boot )
t.fun = function ( data ,i , n ){
bobo = data [ i ]
mean ( bobo [1: n ]) - mean ( bobo [ - c (1: n )])
}

( mouse.boot = boot ( Losratones , t.fun , R =1000 , n = n ))

ORDINARY NONPARAMETRIC BOOTSTRAP

Call :
boot ( data = Losratones , statistic = t.fun , R = 1000 , n = n )

Bootstrap Statistics :
original bias std. error
t1 * 30 .63492 -29 .97352 27 .83245

68/73
Estructuras de datos generales

I El esquema siguiente se aplica a estructuras generales P → x

I En el mundo real se tiene un distribución desconocida P que da


lugar al conjunto de datos x

I El paso principal es el indicado por ⇒ que da lugar un estimador P


b
de la distribución original P.
De este modo P b → x∗ es equivalente a P → x

I Y ası́ x∗ → θb∗ = s(x∗ ) es la misma función que x → θb = s(x).

I Generalmente la generación de muestras bootstrap P b → x∗ requiere


menos tiempo de computación que el cálculo de θb = s(x∗ ).

69/73
70/73
Otro ejemplo

I Se toma el ejemplo de los datos correspondientes a 20 observaciones


donde se mide el efecto de unas pastillas para dormir (el incremento
de horas de sueño en relación a mediciones control).

data ( sleep )

# test t Student
with ( sleep , t.test ( extra∼group ) $ statistic )

t
-1 .860813

71/73
Otro ejemplo
scores = sleep $ extra
R = 1000
t.valores = numeric ( R )

scoresG1 = subset ( scores , sleep $ group ==1)


scoresG2 = subset ( scores , sleep $ group ==2)

for ( i in 1: R ) {
grupo1 = sample ( scoresG1 , size =10 , replace = T )
grupo2 = sample ( scoresG2 , size =10 , replace = T )
t.valores [ i ] = t.test ( grupo1 , grupo2 ) $ statistic
}

sd ( t.valores )

[1] 1 .103293

ggplot2 :: qplot ( t.valores , geom = " histogram " , binwidth =1 ,


fill = I ( " lightgreen " ) , col = I ( " red " ))

72/73
73/73

También podría gustarte