ACTIVIDAD 2

Problema 3. Teorema del Límite Central

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.

A continuación se describen los siguientes pasos para su verificación:

A. Realice una simulación en la cual genere una población de n=1000 (Lote), donde el porcentaje de individuos (supongamos plantas) enfermas sea del 50%.

B. Genere una función que permita: -Obtener una muestra aleatoria de la población y -Calcule el estimador de la proporción muestral pˆ para un tamaño de muestra dado n.

C. Repita el escenario anterior (b) n=500 veces y analice los resultados en cuanto al comportamiento de los 500 resultados del estimador pˆ . ¿Qué tan simétricos o sesgados son los resultados obtenidos? y ¿qué se puede observar en cuanto a la variabilidad?. Realice en su informe un comentario sobre los resultados obtenidos.

D. Repita los puntos b y c para tamaños de muestra n=5 , 10 , 15 , 20 , 30 , 50 , 60 , 100 , 200 , 500 . Compare los resultados obtenidos para los diferentes tamaños de muestra en cuanto a la normalidad. Utilice pruebas de bondad y ajuste (shapiro wilks :shspiro.test()) y métodos gráficos (gráfico de normalidad: qqnorm()). Comente en su informe los resultados obtenidos

E. Repita toda la simulación (puntos a – d), pero ahora para lotes con 10% de plantas enfermas y de nuevo para lotes con un 90% de plantas enfermas. Concluya sobre los resultados del ejercicio.

A Simulación de población n=1000

Para esta actividad se usaron las librerias: dplyr, ggplot2 y magrittr y e1071.

library(dplyr)

## Warning: package 'dplyr' was built under R version 4.2.3

## 
## Attaching package: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(ggplot2)

## Warning: package 'ggplot2' was built under R version 4.2.3

library(magrittr)
library(e1071)

## Warning: package 'e1071' was built under R version 4.2.3

Antes de crear la simulación se generó semilla para su reproductibilidad, posteriormente se creó el vector n con valor igual a 1000 correspondiente al numero total de plantas y el vector de probabilidad de plantas enfermas correspondiente al 50%.

# crear la semilla 
set.seed(123)
# Número total de plantas
n <- 1000 
# Probabilidad de que una planta esté enferma
prob_enferma <- 0.5  

B Generación de la función

Se creó un vector de tamaño n, para cada elemento del vector en donde un valor de 1 indica una planta enferma y un valor de 0 indica una planta sana, tambien se generó el dataframe con los campos id_planta y estado_enfermo.

# Generar la población
estado_plantas <- rbinom(n, 1, prob_enferma)

# Crear dataframe
plantas <- data.frame(id_planta = 1:n, estado_enfermo = estado_plantas)

Posteriormente se generó la muestra aleatoria de esta población y se calculó el estimador de la proporción muestral para un tamaño de muestra dado n=100.

 # Valor de la muestra
n_muestra <- 100 

muestra_plantas <- plantas %>% sample_n(n_muestra)

Se calculó el estimador de la proporción muestral pˆ usando la media, en donde se observa un valor cercano a la proporción de plantas enfermas de 0.53.

estimador <- mean(muestra_plantas$estado_enfermo)
estimador

## [1] 0.53

C 500 simulaciones

# Número de veces que se repite la simulación
n_simulaciones <- 500  
# Función para obtener una muestra  
calcular_estimador <- function() {
  muestra_plantas <- sample_n(plantas, n_muestra)
  return(mean(muestra_plantas$estado_enfermo))
}

# Repetir la simulación n_simulaciones veces y calcular los estimadores
estimadores <- replicate(n_simulaciones, calcular_estimador())

# Analizar los resultados con un resumen estadístico
summary(estimadores)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  0.3600  0.4600  0.4900  0.4939  0.5225  0.6200

hist(estimadores, breaks = 30, main = "Distribución de los estimadores de p^", xlab = "Estimador de p^", ylab = "Frecuencia")

¿Qué tan simétricos o sesgados son los resultados obtenidos? y ¿qué se puede observar en cuanto a la variabilidad?.

De aqui se observa que la mediana y la media son casi iguales, cercanas al valor real de 50% de plantas enfermas, tambien que la mayoria de estimaciones se encuentra con valores entre el 45% y 55%, con una minima de 36% y una maxima de 62%.

En cuanto a la desviación de los estimadores se tiene que es de 0.046 lo cual nos indica que en promedio los estimadores de la proporción muestral se desvían en un 4.6%, lo cual denota que existe poca variabilidad entre ellos.

En cuanto al sesgo se obtuvo un resultado de -0.10, muy cercano a 0, lo cual nos indica que los estimadores no estan sesgados, sino que son más bien simetricos.

# Cálculo de la desviación estándar
desviacion_estandar <- sd(estimadores)

# Cálculo del sesgo (skewness)
sesgo <- skewness(estimadores)

# Mostrar resultados
cat("Desviación estándar de los estimadores:", desviacion_estandar, "\n")

## Desviación estándar de los estimadores: 0.04646357

cat("Sesgo de los estimadores:", sesgo, "\n")

## Sesgo de los estimadores: -0.1059264

D Repetición de los puntos B y C para diferentes tamaños de muestra

Al realizar la simulación para los diferentes tamaños de muestra se tiene que si bien la media es cercana al valor real de 50% desde una muestra incluso menor a 30, esta tiene alta variación, sin embargo en la medida en que la muestra crece, la variación disminuye como es de esperarse, este resultado puede observarse de mejor forma en el grafico de cajas que nos enseña como la varibilidad en los datos disminuye cada vez que crece la muestra hasta el punto en que la caja se encuentra casi plana y los bigotes se hacen muy cortos.

# Tamaños de muestra a evaluar
tamanos_muestra <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)


# Función para realizar las simulaciones y calcular el estimador de p^ para diferentes tamaños de muestra
simular_estimaciones <- function(tamano) {
  estimaciones <- replicate(n_simulaciones, {
    muestra_plantas <- sample_n(plantas, tamano)
    mean(muestra_plantas$estado_enfermo)
  })
  data.frame(TamanoMuestra = tamano, EstimadorP = estimaciones)
}

# Aplicar la simulación para cada tamaño de muestra y combinar los resultados
resultados <- do.call(rbind, lapply(tamanos_muestra, simular_estimaciones))

# Crear un resumen estadístico para cada tamaño de muestra
resumen_estadistico <- resultados %>%
  group_by(TamanoMuestra) %>%
  summarise(Media = mean(EstimadorP),
            DesviacionStd = sd(EstimadorP),
            Min = min(EstimadorP),
            Max = max(EstimadorP),
            Q1 = quantile(EstimadorP, 0.25),
            Q3 = quantile(EstimadorP, 0.75))

# Mostrar el resumen estadístico
print(resumen_estadistico)

## # A tibble: 10 × 7
##    TamanoMuestra Media DesviacionStd   Min   Max    Q1    Q3
##            <dbl> <dbl>         <dbl> <dbl> <dbl> <dbl> <dbl>
##  1             5 0.503        0.226  0     1     0.4   0.6  
##  2            10 0.500        0.157  0     0.9   0.4   0.6  
##  3            15 0.51         0.122  0.2   0.8   0.4   0.6  
##  4            20 0.490        0.108  0.2   0.85  0.4   0.55 
##  5            30 0.501        0.0861 0.233 0.767 0.433 0.567
##  6            50 0.495        0.0666 0.3   0.7   0.46  0.54 
##  7            60 0.494        0.0600 0.317 0.65  0.45  0.533
##  8           100 0.493        0.0479 0.36  0.62  0.46  0.52 
##  9           200 0.493        0.0309 0.41  0.585 0.475 0.51 
## 10           500 0.493        0.0151 0.442 0.54  0.482 0.502

# Opcional: Crear un gráfico boxplot para visualizar la distribución de los estimadores para cada tamaño de muestra
ggplot(resultados, aes(x = factor(TamanoMuestra), y = EstimadorP)) +
  geom_boxplot() +
  xlab("Tamaño de la muestra") +
  ylab("Estimador de p^") +
  ggtitle("Distribución de los estimadores de p^ por tamaño de muestra")

Se compararon los resultados obtenidos para los diferentes tamaños de muestra usando shapiro wilks y grafico de normalidad

# Aplicar la prueba de Shapiro-Wilk por tamaño de muestra
pruebas_shapiro <- resultados %>%
  group_by(TamanoMuestra) %>%
  summarise(p_valor_shapiro = shapiro.test(EstimadorP)$p.value)

# Mostrar los p-valores de la prueba
print(pruebas_shapiro)

## # A tibble: 10 × 2
##    TamanoMuestra p_valor_shapiro
##            <dbl>           <dbl>
##  1             5        1.39e-14
##  2            10        5.83e-10
##  3            15        2.05e- 8
##  4            20        6.62e- 7
##  5            30        1.02e- 4
##  6            50        2.09e- 3
##  7            60        2.97e- 3
##  8           100        1.77e- 2
##  9           200        2.70e- 2
## 10           500        1.44e- 1

Analisis de resultados

Para los tamaños de muestra de 5 a 200 todos los valores son menores que 0.05, lo que nos indica que para estos tamaños de muestra los datos no siguen una distribución normal, para el tamaño de muestra 500 por su parte el valor es 0.1440891 mayor que 0.05 indicando que para este tamaño de muestral los estimadores estan siguiendo una distribución normal, se observa que a medida que el tamaño de muestra aumenta, gracias al Teorema del Límite Central la distribución de los estimadores tiende a aproximarse a una normal.

library(dplyr)


# Lista de tamaños de muestra a visualizar
tamanos_muestra <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)

# Bucle para generar y visualizar un gráfico Q-Q por cada tamaño de muestra
for (tamano in tamanos_muestra) {
  # Filtrar los datos para el tamaño de muestra actual
  datos_ejemplo <- filter(resultados, TamanoMuestra == tamano)
  
  # Realizar gráfico Q-Q para el tamaño de muestra seleccionado
  qqnorm(datos_ejemplo$EstimadorP, main = paste("Q-Q Plot, tamaño de muestra =", tamano))
  qqline(datos_ejemplo$EstimadorP, col = "red")
  
  Sys.sleep(1)  
}

El grafico de normalidad nos muestra como los puntos pasan de estar desalineados con una muestra n=5 y empiezan a alinearse a partir de los 20 puntos, con algunos puntos de las colas desviados, se observa como para la muestra n=500 los puntos se encuentran alineados con la línea de referencia, incluyendo las colas, indicando normalidad en la distribución de los datos.

E Repetición de la simulación para lotes con 10% de plantas enfermas

# Establecer la semilla para reproducibilidad
set.seed(123)

# Función para realizar simulaciones y análisis para el 10% de plantas enfermas
simular_y_analizar <- function(prob_enferma) {
  n <- 1000  # Tamaño de la población
  tamanos_muestra <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)  # Tamaños de muestra a evaluar
  n_simulaciones <- 500  # Número de simulaciones
  
  # Generar la población con el porcentaje especificado de plantas enfermas
  estado_plantas <- rbinom(n, 1, prob_enferma)
  plantas <- data.frame(id_planta = 1:n, estado_enfermo = estado_plantas)
  
  # Función para realizar las simulaciones y calcular el estimador de p^ para diferentes tamaños de muestra
  simular_estimaciones <- function(tamano) {
    estimaciones <- replicate(n_simulaciones, {
      muestra_plantas <- sample_n(plantas, tamano)
      mean(muestra_plantas$estado_enfermo)
    })
    data.frame(TamanoMuestra = tamano, EstimadorP = estimaciones)
  }
  
  # Aplicar la simulación para cada tamaño de muestra y combinar los resultados
  resultados <- do.call(rbind, lapply(tamanos_muestra, simular_estimaciones))
  
  # Realizar prueba de Shapiro-Wilk solo para tamaños de muestra adecuados
  resultados %>%
    group_by(TamanoMuestra) %>%
    summarise(Media = mean(EstimadorP),
              DesviacionStd = sd(EstimadorP),
              ShapiroP = ifelse(n() >= 3 & n() <= 5000, shapiro.test(EstimadorP)$p.value, NA)) %>%
    print()
  
  # Opcional: Crear un gráfico boxplot para visualizar la distribución de los estimadores para cada tamaño de muestra
  ggplot(resultados, aes(x = factor(TamanoMuestra), y = EstimadorP)) +
    geom_boxplot() +
    xlab("Tamano de la muestra") +
    ylab("Estimador de p^") +
    ggtitle(paste("Distribucion de los estimadores de p^ para prob_enferma =", prob_enferma))
}

# Simulación con 10% de plantas enfermas
simular_y_analizar(0.1)

## # A tibble: 10 × 4
##    TamanoMuestra  Media DesviacionStd ShapiroP
##            <dbl>  <dbl>         <dbl>    <dbl>
##  1             5 0.092        0.135   1.30e-29
##  2            10 0.084        0.0860  6.93e-24
##  3            15 0.0925       0.0712  3.78e-18
##  4            20 0.095        0.0630  2.05e-16
##  5            30 0.0969       0.0519  3.25e-11
##  6            50 0.0934       0.0398  2.31e- 9
##  7            60 0.0923       0.0367  1.60e- 7
##  8           100 0.0933       0.0277  1.19e- 4
##  9           200 0.0911       0.0178  3.42e- 3
## 10           500 0.0921       0.00960 8.31e- 2

# Ajustes previos: Asegúrate de tener el dataframe 'resultados' correctamente definido
# resultados <- tu_dataframe_con_resultados

# Lista de tamaños de muestra
tamanos_muestra <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)

# Loop para generar y visualizar un gráfico Q-Q por cada tamaño de muestra
for (tamano in tamanos_muestra) {
  # Filtrar los datos para el tamaño de muestra actual
  datos_ejemplo <- filter(resultados, TamanoMuestra == tamano)
  
  # Crear y mostrar el gráfico Q-Q
  qqnorm(datos_ejemplo$EstimadorP, main = paste("Q-Q Plot, Size =", tamano))
  qqline(datos_ejemplo$EstimadorP, col = "red")
  
  # Imprimir el mensaje para proceder manualmente en entornos interactivos
  # Remover o ajustar según sea necesario para tu flujo de trabajo
  cat("Displayed Q-Q Plot for sample size:", tamano, "\n")
  # Pausa para visualizar el gráfico; ajustar según necesidad o remover si no deseado
  # Sys.sleep(time) # time en segundos, ajusta a tu preferencia
  readline(prompt = "Press [Enter] to continue to the next plot...")
}

## Displayed Q-Q Plot for sample size: 5 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 10 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 15 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 20 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 30 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 50 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 60 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 100 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 200 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for sample size: 500 
## Press [Enter] to continue to the next plot...

E Repetición de la simulación para lotes con 90% de plantas enfermas

# Establecer la semilla para reproducibilidad
set.seed(123)

# Función ajustada para realizar simulaciones y análisis con un 90% de plantas enfermas
simular_y_analizar_90 <- function(prob_enferma) {
  n <- 1000  # Tamaño de la población
  tamanos_muestra_90 <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)  # Tamaños de muestra a evaluar
  n_simulaciones <- 500  # Número de simulaciones
  
  # Generar la población con el 90% de plantas enfermas
  estado_plantas <- rbinom(n, 1, prob_enferma)
  plantas <- data.frame(id_planta = 1:n, estado_enfermo = estado_plantas)
  
  # Función para realizar las simulaciones y calcular el estimador de p^ para diferentes tamaños de muestra
  simular_estimaciones <- function(tamano) {
    estimaciones <- replicate(n_simulaciones, {
      muestra_plantas <- sample_n(plantas, tamano)
      mean(muestra_plantas$estado_enfermo)
    })
    data.frame(TamanoMuestra = tamano, EstimadorP = estimaciones)
  }
  
  # Aplicar la simulación para cada tamaño de muestra y combinar los resultados
  resultados_90 <- do.call(rbind, lapply(tamanos_muestra, simular_estimaciones))
  
  # Realizar prueba de Shapiro-Wilk solo para tamaños de muestra adecuados
  resultados_90 %>%
    group_by(TamanoMuestra) %>%
    summarise(Media = mean(EstimadorP),
              DesviacionStd = sd(EstimadorP),
              ShapiroP = ifelse(n() >= 3 & n() <= 5000, shapiro.test(EstimadorP)$p.value, NA)) %>%
    print()
  
  # Crear un gráfico boxplot para visualizar la distribución de los estimadores para cada tamaño de muestra
  ggplot(resultados_90, aes(x = factor(TamanoMuestra), y = EstimadorP)) +
    geom_boxplot() +
    xlab("Tamano de la muestra") +
    ylab("Estimador de p^") +
    ggtitle(paste("Distribucion de los estimadores de p^ para prob_enferma =", prob_enferma))
}

# Simulación con 90% de plantas enfermas
simular_y_analizar_90(0.9)

## # A tibble: 10 × 4
##    TamanoMuestra Media DesviacionStd ShapiroP
##            <dbl> <dbl>         <dbl>    <dbl>
##  1             5 0.908       0.135   1.30e-29
##  2            10 0.916       0.0860  6.93e-24
##  3            15 0.907       0.0712  3.78e-18
##  4            20 0.905       0.0630  2.05e-16
##  5            30 0.903       0.0519  3.25e-11
##  6            50 0.907       0.0398  2.31e- 9
##  7            60 0.908       0.0367  1.60e- 7
##  8           100 0.907       0.0277  1.19e- 4
##  9           200 0.909       0.0178  3.42e- 3
## 10           500 0.908       0.00960 8.31e- 2

# Asegúrate de que el dataframe 'resultados_90' está correctamente definido
# Aquí deberías tener ya el dataframe 'resultados_90' listo y disponible

# Lista de tamaños de muestra
tamanos_muestra <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)

  # Aplicar la simulación para cada tamaño de muestra y combinar los resultados
  resultados_90 <- do.call(rbind, lapply(tamanos_muestra, simular_estimaciones))
  

# Loop para generar y visualizar un gráfico Q-Q por cada tamaño de muestra para 'resultados_90'
for (tamano in tamanos_muestra) {
  # Filtrar los datos para el tamaño de muestra actual de 'resultados_90'
  datos_ejemplo_90 <- filter(resultados_90, TamanoMuestra == tamano)
  
  # Crear y mostrar el gráfico Q-Q para 'resultados_90'
  qqnorm(datos_ejemplo_90$EstimadorP, main = paste("Q-Q Plot, 90% Sick, Size =", tamano))
  qqline(datos_ejemplo_90$EstimadorP, col = "red")
  
  # Mensaje para proceder manualmente en entornos interactivos
  cat("Displayed Q-Q Plot for 90% sick, sample size:", tamano, "\n")
  readline(prompt="Press [Enter] to continue to the next plot...")
}

## Displayed Q-Q Plot for 90% sick, sample size: 5 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 10 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 15 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 20 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 30 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 50 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 60 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 100 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 200 
## Press [Enter] to continue to the next plot...

## Displayed Q-Q Plot for 90% sick, sample size: 500 
## Press [Enter] to continue to the next plot...

Análisis de los resultados para 10% y 90% de plantas enfermas

De aqui se observa que la mediana y la media al igual que para las plantas con 50% son cercanas cercanas al valor real, tambien que la mayoria de estimaciones se encuentra con valores entre el 0% y 20% para el 10% de probabilidad de plantas enfermas y el 80% y 100% para el 90% de probabilidad de plantas enfermas.

En cuanto a la desviación de los estimadores en ambos casos tambien es cercana a 0 lo cual nos indica existe poca variabilidad entre ellos y tambien que el tamaño de la muestra tiene una relación directa con la variación.

En cuanto al sesgo se observa que para este tipo de proporciones en donde la mayoria de población esta enferma, los valores se encuentran más sesgados cuando las muestras son menores de acuerdo con el grafico de normalidad,sin embargo y en la medida en que se incrementa la muestra, el comportamiento tiende hacia la normalidad.

Conclusiones

La actividad desarrollada ha permitido probar las propiedades del teorema de limite central en cuanto se han usado diferentes tamaños de muestras en poblaciones de plantas totalmente distintas, con 50%,10% y 90% de plantas enfermas, en donde se ha logrado probar que a partir de n=>30 se logran alcanzar las propiedades de estimadores esperadas, permitiendo realizar inferencias sobre la población a partir de muestras aleatorias.