1. El Teorema del Límite Central es uno de los más importantes en la inferencia estadística y habla sobre la convergencia de los estimadores como la proporción muestral a la distribución normal. Algunos autores afirman que esta aproximación es bastante buena a partir del umbral n>30.

* Realice una simulación en la cual genere una población de N=1000 (Lote) y además que el porcentaje de individuos (plantas) enfermas sea del 50%.

set.seed(845)
plantas <- rep(c(F,T), 500)
plantas <- sample(plantas, 1000)
table(plantas)
## plantas
## FALSE  TRUE 
##   500   500

* Genere una función que permita obtener una muestra aleatoria de la población y calcule el estimador de la proporción muestral para un tamaño de muestra dado n

prop_muestra <- function(data, n){
  x <- sample(data, n)
  return(sum(x)/n)
}

prop_muestra(plantas, 49)
## [1] 0.5102041

* Repita el escenario anterior (b) 500 veces y analice los resultados en cuanto al comportamiento de los 500 estimadores. ¿Qué tan simétricos son los datos?, ¿Son sesgados y qué pasa en cuanto a variabilidad?

x <- c()
for (i in 1:500) {
  x[i] <- prop_muestra(plantas, 49)
}

hist(x, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)

data.frame(media = mean(x), 
           desviacion = sd(x))
##       media desviacion
## 1 0.4952245 0.07007709

Lo que se aprecia en el histograma es que la distribución de los estimadores es simétrica, parece tener una distribución normal alrededor de la media poblacional. Es consistente puesto que la media tiende a ser igual a la media poblacional a medida que aumenta la muestra y es insesgada porque la media de la muestra es igual (casi) a la media poblacional

* Realice los ejercicios completos b y c para tamaños de muestra n=5, 10, 15, 20, 30, 50, 60, 100, 200, 500. Y compare los resultados de los estimadores en cuanto a la normalidad. Investigue y utilice pruebas de bondad y ajuste (shapiro wilks) y métodos gráficos (grafico qq de normalidad).

# Muestra n = 5
x_5 <- c()
for (i in 1:500) {
  x_5[i] <- prop_muestra(plantas, 5)
}

hist(x_5, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_5), col = "red")

n_5 <- data.frame(media = mean(x_5), 
           desviacion = sd(x_5))
qqPlot(x_5)

## [1] 11 14
shapiro.test(x_5)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_5
## W = 0.93113, p-value = 2.102e-14
# Muestra n = 10
x_10 <- c()
for (i in 1:500) {
  x_10[i] <- prop_muestra(plantas, 10)
}

hist(x_10, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_10), col = "red")

n_10 <- data.frame(media = mean(x_10), 
           desviacion = sd(x_10))
qqPlot(x_10)

## [1] 43 50
shapiro.test(x_10)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_10
## W = 0.96284, p-value = 6.335e-10
# Muestra n = 15
x_15 <- c()
for (i in 1:500) {
  x_15[i] <- prop_muestra(plantas, 15)
}

hist(x_15, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_15), col = "red")

n_15 <- data.frame(media = mean(x_15), 
           desviacion = sd(x_15))
qqPlot(x_15)

## [1]  40 241
shapiro.test(x_15)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_15
## W = 0.97486, p-value = 1.431e-07
# Muestra n = 20
x_20 <- c()
for (i in 1:500) {
  x_20[i] <- prop_muestra(plantas, 20)
}

hist(x_20, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_20), col = "red")

n_20 <- data.frame(media = mean(x_20), 
           desviacion = sd(x_20))
qqPlot(x_20)

## [1] 190 292
shapiro.test(x_20)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_20
## W = 0.98004, p-value = 2.352e-06
# Muestra n = 30
x_30 <- c()
for (i in 1:500) {
  x_30[i] <- prop_muestra(plantas, 30)
}

hist(x_30, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_30), col = "red")

n_30 <- data.frame(media = mean(x_30), 
           desviacion = sd(x_30))
qqPlot(x_30)

## [1] 395 336
shapiro.test(x_30)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_30
## W = 0.98422, p-value = 2.998e-05
# Muestra n = 40
x_40 <- c()
for (i in 1:500) {
  x_40[i] <- prop_muestra(plantas, 40)
}

hist(x_40, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_40), col = "red")

n_40 <- data.frame(media = mean(x_40), 
           desviacion = sd(x_40))
qqPlot(x_40)

## [1] 61 32
shapiro.test(x_40)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_40
## W = 0.98652, p-value = 0.0001391
# Muestra n = 50
x_50 <- c()
for (i in 1:500) {
  x_50[i] <- prop_muestra(plantas, 50)
}

hist(x_50, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_50), col = "red")

n_50 <- data.frame(media = mean(x_50), 
           desviacion = sd(x_50))
qqPlot(x_50)

## [1] 205 416
shapiro.test(x_50)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_50
## W = 0.99167, p-value = 0.006615
# Muestra n = 60
x_60 <- c()
for (i in 1:500) {
  x_60[i] <- prop_muestra(plantas, 60)
}

hist(x_60, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_60), col = "red")

n_60 <- data.frame(media = mean(x_60), 
           desviacion = sd(x_60))
qqPlot(x_60)

## [1] 398 468
shapiro.test(x_60)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_60
## W = 0.99174, p-value = 0.007018
# Muestra n = 100
x_100 <- c()
for (i in 1:500) {
  x_100[i] <- prop_muestra(plantas, 100)
}

hist(x_100, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_100), col = "red")

n_100 <- data.frame(media = mean(x_100), 
           desviacion = sd(x_100))
qqPlot(x_100)

## [1] 428 424
shapiro.test(x_100)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_100
## W = 0.99352, p-value = 0.03062
# Muestra n = 200
x_200 <- c()
for (i in 1:500) {
  x_200[i] <- prop_muestra(plantas, 200)
}

hist(x_200, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_200), col = "red")

n_200 <- data.frame(media = mean(x_200), 
           desviacion = sd(x_200))
qqPlot(x_200)

## [1]  55 394
shapiro.test(x_200)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_200
## W = 0.99607, p-value = 0.2519
# Muestra n = 500
x_500 <- c()
for (i in 1:500) {
  x_500[i] <- prop_muestra(plantas, 500)
}

hist(x_500, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_500), col = "red")

n_500 <- data.frame(media = mean(x_500), 
           desviacion = sd(x_500))
qqPlot(x_500)

## [1] 114 251
shapiro.test(x_500)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_500
## W = 0.9954, p-value = 0.1467
rbind(n_5, n_10, n_15, n_20, n_30, n_40, n_50, n_60, n_100, n_200, n_500)
##        media desviacion
## 1  0.5080000 0.23316245
## 2  0.5000000 0.16362624
## 3  0.4950667 0.12641581
## 4  0.5000000 0.11115174
## 5  0.4997333 0.09061958
## 6  0.5028500 0.07393619
## 7  0.4989200 0.06948734
## 8  0.4986000 0.06258582
## 9  0.4990000 0.04682933
## 10 0.5006700 0.03069014
## 11 0.4998080 0.01610870

Se observa que a medida que el tamaño de las muestras crece, la distribución de los estimadores se hace cada vez más parecida a la distribución normal. Además la media se acerca más a la media poblacional y la varianza disminuye, probando la aplicación del Teorema de Límite Central. Muestra simetría, consistencia e insesgamiento. Además, para apoyar lo anterior, los test de shapiro comienzan a reflejar normalidad a partir de una muestra mayor a 30 y los graficos qq también evidencian una tendencia hacia la normalidad de la distribución a medida que aumenta el tamaño de la muestra.

* Repita toda la simulación (puntos a – d) pero ahora con lotes con 10% y 90% de plantas enfermas. Concluya todo el ejercicio.

plantas <- rep(c(F,T), c(100, 900))
plantas <- sample(plantas, 1000)
table(plantas)
## plantas
## FALSE  TRUE 
##   100   900

* Genere una función que permita obtener una muestra aleatoria de la población y calcule el estimador de la proporción muestral para un tamaño de muestra dado n

prop_muestra <- function(data, n){
  x <- sample(data, n)
  return(sum(x)/n)
}

prop_muestra(plantas, 49)
## [1] 0.8979592

* Repita el escenario anterior (b) 500 veces y analice los resultados en cuanto al comportamiento de los 500 estimadores. ¿Qué tan simétricos son los datos?, ¿Son sesgados y qué pasa en cuanto a variabilidad?

x <- c()
for (i in 1:500) {
  x[i] <- prop_muestra(plantas, 49)
}

hist(x, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)

data.frame(media = mean(x), 
           desviacion = sd(x))
##       media desviacion
## 1 0.8951429  0.0423368

Lo que se aprecia en el histograma es que la distribución de las medias de los estimadores es simétrica, parece tener una distribución normal alrededor de la media poblacional. Es consistente puesto que la media tiende a ser igual a la media poblacional a medida que aumenta la muestra y es insesgada porque la media de la muestra es igual (casi) a la media poblacional

* Realice los ejercicios completos b y c para tamaños de muestra n=5, 10, 15, 20, 30, 50, 60, 100, 200, 500. Y compare los resultados de los estimadores en cuanto a la normalidad. Investigue y utilice pruebas de bondad y ajuste (shapiro wilks) y métodos gráficos (grafico qq de normalidad).

# Muestra n = 5
x_5 <- c()
for (i in 1:500) {
  x_5[i] <- prop_muestra(plantas, 5)
}

hist(x_5, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_5), col = "red")

n_5 <- data.frame(media = mean(x_5), 
           desviacion = sd(x_5))
qqPlot(x_5)

## [1] 149 218
shapiro.test(x_5)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_5
## W = 0.71345, p-value < 2.2e-16
# Muestra n = 10
x_10 <- c()
for (i in 1:500) {
  x_10[i] <- prop_muestra(plantas, 10)
}

hist(x_10, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_10), col = "red")

n_10 <- data.frame(media = mean(x_10), 
           desviacion = sd(x_10))
qqPlot(x_10)

## [1]  17 209
shapiro.test(x_10)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_10
## W = 0.84109, p-value < 2.2e-16
# Muestra n = 15
x_15 <- c()
for (i in 1:500) {
  x_15[i] <- prop_muestra(plantas, 15)
}

hist(x_15, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_15), col = "red")

n_15 <- data.frame(media = mean(x_15), 
           desviacion = sd(x_15))
qqPlot(x_15)

## [1] 450   2
shapiro.test(x_15)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_15
## W = 0.89827, p-value < 2.2e-16
# Muestra n = 20
x_20 <- c()
for (i in 1:500) {
  x_20[i] <- prop_muestra(plantas, 20)
}

hist(x_20, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_20), col = "red")

n_20 <- data.frame(media = mean(x_20), 
           desviacion = sd(x_20))
qqPlot(x_20)

## [1] 489 228
shapiro.test(x_20)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_20
## W = 0.92785, p-value = 8.828e-15
# Muestra n = 30
x_30 <- c()
for (i in 1:500) {
  x_30[i] <- prop_muestra(plantas, 30)
}

hist(x_30, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_30), col = "red")

n_30 <- data.frame(media = mean(x_30), 
           desviacion = sd(x_30))
qqPlot(x_30)

## [1]  72 177
shapiro.test(x_30)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_30
## W = 0.94925, p-value = 4.528e-12
# Muestra n = 40
x_40 <- c()
for (i in 1:500) {
  x_40[i] <- prop_muestra(plantas, 40)
}

hist(x_40, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_40), col = "red")

n_40 <- data.frame(media = mean(x_40), 
           desviacion = sd(x_40))
qqPlot(x_40)

## [1] 231 112
shapiro.test(x_40)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_40
## W = 0.96807, p-value = 5.742e-09
# Muestra n = 50
x_50 <- c()
for (i in 1:500) {
  x_50[i] <- prop_muestra(plantas, 50)
}

hist(x_50, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_50), col = "red")

n_50 <- data.frame(media = mean(x_50), 
           desviacion = sd(x_50))
qqPlot(x_50)

## [1]  21 345
shapiro.test(x_50)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_50
## W = 0.96635, p-value = 2.723e-09
# Muestra n = 60
x_60 <- c()
for (i in 1:500) {
  x_60[i] <- prop_muestra(plantas, 60)
}

hist(x_60, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_60), col = "red")

n_60 <- data.frame(media = mean(x_60), 
           desviacion = sd(x_60))
qqPlot(x_60)

## [1] 110 159
shapiro.test(x_60)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_60
## W = 0.97672, p-value = 3.766e-07
# Muestra n = 100
x_100 <- c()
for (i in 1:500) {
  x_100[i] <- prop_muestra(plantas, 100)
}

hist(x_100, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_100), col = "red")

n_100 <- data.frame(media = mean(x_100), 
           desviacion = sd(x_100))
qqPlot(x_100)

## [1] 48 47
shapiro.test(x_100)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_100
## W = 0.98324, p-value = 1.606e-05
# Muestra n = 200
x_200 <- c()
for (i in 1:500) {
  x_200[i] <- prop_muestra(plantas, 200)
}

hist(x_200, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_200), col = "red")

n_200 <- data.frame(media = mean(x_200), 
           desviacion = sd(x_200))
qqPlot(x_200)

## [1] 271 428
shapiro.test(x_200)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_200
## W = 0.99269, p-value = 0.01527
# Muestra n = 500
x_500 <- c()
for (i in 1:500) {
  x_500[i] <- prop_muestra(plantas, 500)
}

hist(x_500, main = "Histograma de Estimadores", xlim = range(0,1), breaks = 20)
abline(v = mean(x_500), col = "red")

n_500 <- data.frame(media = mean(x_500), 
           desviacion = sd(x_500))
qqPlot(x_500)

## [1] 242 265
shapiro.test(x_500)
## 
##  Shapiro-Wilk normality test
## 
## data:  x_500
## W = 0.99266, p-value = 0.01489
rbind(n_5, n_10, n_15, n_20, n_30, n_40, n_50, n_60, n_100, n_200, n_500)
##        media desviacion
## 1  0.8988000 0.13249025
## 2  0.9022000 0.09291084
## 3  0.8986667 0.08067849
## 4  0.8995000 0.06666061
## 5  0.8966000 0.05658188
## 6  0.8990500 0.04479800
## 7  0.8967200 0.04117076
## 8  0.8991333 0.03673862
## 9  0.8993000 0.02970985
## 10 0.8995600 0.01990474
## 11 0.9000760 0.01002933

Lo que se aprecia en el histograma es que la distribución de las medias de los estimadores es simétrica, parece tener una distribución normal alrededor de la media poblacional. Es consistente puesto que la media tiende a ser igual a la media poblacional a medida que aumenta la muestra y es insesgada porque la media de la muestra es igual (casi) a la media poblacional. Esto prueba que cuando el tamaño de la muestra es grande, la distribución de las variable aleatoria correspondiente a las medias de los estimadores será normal sin importar la distribución subyacente de la población.