Introducción

En el ámbito de la estadística, el proceso de simulación emerge como una herramienta poderosa que permite abordar relaciones complejas y estimar valores difíciles de calcular directamente. Esta técnica, que involucra la generación de datos sintéticos mediante modelos matemáticos y la posterior exploración de múltiples escenarios, ofrece una vía para comprender la dinámica subyacente de diversos fenómenos y tomar decisiones.

En este trabajo, se desarrollan cuatro actividades con el fin de explorar el potencial de la simulación estadística apoyada en el software R. Estas actividades buscan principalmente comprender cómo la simulación puede aplicarse en diversos contextos.

Problema 1. Estimación del valor de \(\pi\)

A continuación, se presenta un procedimiento para calcular el valor de \(\pi\) mediante estimaciones:

Tomemos un círculo inscrito dentro de un cuadrado de lado igual a 1. El radio del círculo es \(\frac{1}{2}\), por lo tanto, el área del círculo es \(\frac{\pi}{4}\). Si ubicamos puntos aleatorios dentro del cuadrado, la probabilidad de que un punto de estos esté dentro del círculo es igual a la porción de área que el círculo abarca dentro del cuadrado, es decir, la probabilidad es \(\frac{\pi}{4}\).

Ubiquemos 1000 puntos aleatorios dentro del cuadrado de lado 1 y determinemos cuantos puntos de estos están dentro del círculo inscrito.

set.seed(123)
n <- 1000 
x <- runif(n)
y <- runif(n)


plot(c(0, 1), c(0, 1), type="n", asp=1, xlab="", ylab="", )
rect(0, 0, 1, 1, border="blue") 
theta <- seq(0, 2*pi, length.out=100)
lines(0.5 + 0.5*cos(theta), 0.5 + 0.5*sin(theta), col="red") 


points(x, y, col="black", pch=20)

Para determinar si un punto está dentro del círculo o no, utilizamos el concepto de distancia. Si ubicamos el cuadrado en un plano cartesiano y tomamos vértices en \((0 , 0)\),\((0 , 1)\),\((1 , 1)\) y \((1 , 0)\) , el círculo inscrito tendría centro\((0.5 , 0.5)\); entonces, un punto dentro del cuadrado, está también dentro de la circunferencia si su distancia hasta el centro del círculo es menor a 0.5.

.

distancias <- sqrt((x - 0.5)^2 + (y - 0.5)^2)
radio <- 0.5
puntos_dentro <- sum(distancias < radio)
cat("Número de puntos dentro de la circunferencia:", puntos_dentro)
## Número de puntos dentro de la circunferencia: 800

Con ayuda de r, se obtuvo que el número de puntos dentro del círculo es 800, por lo tanto, usando el enfoque clásico de la probabilidad, se obtiene que la probabilidad de que un punto esté dentro del círculo inscrito es igual a \(\frac{800}{1000}\), simplificando esto es \(\frac{4}{5}\), pero , desde el principio sabemos que la probabilidad de que un punto que esté dentro del cuadrado también esté dentro del círculo es \(\frac{\pi}{4}\), lo que indica que \(\frac{4}{5} = \frac{\pi}{4}\), de donde se obtiene que \[ \pi = \frac{16}{5} = 3.2\]

observemos que el error es mayor a 0.05 para la simulación con 1000 puntos, con el fin de mejorar el resultado ahora trabajaremos con 10000 puntos y con 100000 puntos.

Simulación usando 10000 puntos

set.seed(123)
n <- 10000 
x <- runif(n)
y <- runif(n)



distancias <- sqrt((x - 0.5)^2 + (y - 0.5)^2)
radio <- 0.5
puntos_dentro <- sum(distancias < radio)

cat("Número de puntos dentro de la circunferencia:", puntos_dentro)
## Número de puntos dentro de la circunferencia: 7894
calcular_pi <- function(a, b) {
  cociente <- 4*(a / b)  
  return(cociente)  
}


resultado <- calcular_pi(puntos_dentro, n)
cat("El valor de pi es :",resultado) 
## El valor de pi es : 3.1576

Simulación usando 100000 puntos

set.seed(123)
n <- 100000 
x <- runif(n)
y <- runif(n)


distancias <- sqrt((x - 0.5)^2 + (y - 0.5)^2)
radio <- 0.5
puntos_dentro <- sum(distancias < radio)

cat("Número de puntos dentro de la circunferencia:", puntos_dentro)
## Número de puntos dentro de la circunferencia: 78658
calcular_pi <- function(a, b) {
  cociente <- 4*(a / b)  
  return(cociente)  
}


resultado <- calcular_pi(puntos_dentro, n)
cat("El valor de pi es :",resultado) 
## El valor de pi es : 3.14632

Resultados

A continuación se presentan los resultados obtenidos en las simulaciones

Tabla de resultados simulación

Puntos Estimación Error
1000 3.2 aprox 0.5
10000 3.1576 aprox 0.016
100000 3.14632 aprox 0,0047

Conclusiones

Se evidencia que a medida que aumenta el tamaño de la muestra de números aleatorios generados, el error relativo disminuye. Esto indica que la estimación se acerca cada vez más al valor real de π.

Los resultados obtenidos demuestran que el método utilizado es una aproximación válida para π. La estrategia de generar puntos aleatorios dentro del cuadrado inscrito en el círculo y calcular la fracción de puntos que caen dentro del círculo proporciona una estimación razonablemente precisa del valor de π.

Limitaciones:

La precisión de la estimación está limitada por el número de puntos generados aleatoriamente. Aunque aumentar el número de puntos mejora la precisión de la estimación, hay un límite en cuánto puede mejorar la precisión simplemente aumentando la muestra.

A medida que aumenta el número de puntos, se requiere más tiempo y recursos computacionales para el procesamiento. Por lo tanto, hay una compensación entre la precisión deseada y el costo computacional asociado con la generación y procesamiento de un gran número de puntos. Es importante encontrar un equilibrio entre la precisión deseada y los recursos disponibles.

Problema 2. Propiedades de los estimadores.

Los estimadores son variables aleatorias que se basan en los datos recolectados en una muestra y siguen una distribución conocida llamada distribución muestral. Estos estimadores se destacan por sus propiedades deseables, como la insesgadez, eficiencia y consistencia. La simulación desempeña un rol crucial en la comprensión y validación de estas propiedades.

A continuación, se presenta un problema diseñado para resaltar las características principales de un conjunto de estimadores propuestos para la estimación de un parámetro asociado a un modelo de probabilidad.

Primero tomemos una muestra aleatoria de tamaño \(n=4\) , ( \(X_{1}, X_{2}, X_{3}, X_{4}\)), cuya población la conforma una distribución exponencial con parámetro \(\theta =1\)

Ahora definamos los siguientes estimadores:

\[ \hat{\theta}_{1} = \frac{X_{1}+X_{2}}{6} + \frac{X_{3}+X_{4}}{3} \]

\[ \hat{\theta}_{2} = \frac{X_{1}+2X_{2}+3X_{3}+4X_{4}}{5} \]

\[\hat{\theta}_{3}=\frac{X_{1}+X_{2}+X_{3}+X_{4}}{4} \]

\[ \hat{\theta_4 }=\frac{min(X_1, X_2,X_3, X_4)+max(X_1, X_2,X_3, X_4) }{2} \]

estimadores=function(theta){
  x=rexp(4,rate = 1)
  theta1=((x[1]+x[2])/6)+((x[3]+x[4])/3)
  theta2=((x[1]+2*x[2]+3*x[3]+4*x[4])/5)
  theta3=((x[1]+x[2]+x[3]+x[4])/4)
  theta4=(min(x) + max(x))/2
  theta=c(theta1,theta2,theta3,theta4)
  return(theta)
}

Caso 1, muestra de n=20

Con la muestra aleatoria de tamaño \(n=4\) , ( \(X_{1}, X_{2}, X_{3}, X_{4}\)), cuya población la conforma una distribución exponencial con parámetro \(\theta =1\)

resultados <- t(sapply(1:20, function(x) estimadores()))


resultados <- as.data.frame(resultados)
colnames(resultados) <- c("Estimador 1", "Estimador 2", "Estimador 3", "Estimador 4")

Histogramas

par(mfrow = c(2, 2)) # Dividir el área de la trama en 2x2
library(ggplot2)
for (j in 1:4) {
  hist(resultados[, j], main = paste("Estimador", j), xlab = "Valor estimado")
}

library(ggplot2)


boxplot(resultados, las=1, main="Comparación estimadores con n=20")  # gráfico de comparación   
abline(h=1,  col="red")     

apply(resultados,2,mean)   
## Estimador 1 Estimador 2 Estimador 3 Estimador 4 
##    0.934147    1.911135    1.007984    1.331279
apply(resultados,2,var) 
## Estimador 1 Estimador 2 Estimador 3 Estimador 4 
##   0.3164758   1.4492040   0.4020187   1.0802144
apply(resultados,2,sd)
## Estimador 1 Estimador 2 Estimador 3 Estimador 4 
##   0.5625618   1.2038289   0.6340494   1.0393336

Caso 2, muestra de n=50

Con la muestra aleatoria de tamaño \(n=4\) , ( \(X_{1}, X_{2}, X_{3}, X_{4}\)), cuya población la conforma una distribución exponencial con parámetro \(\theta =1\)

estimadores=function(theta){
  x=rexp(4,rate = 1)
  theta1=((x[1]+x[2])/6)+((x[3]+x[4])/3)
  theta2=((x[1]+2*x[2]+3*x[3]+4*x[4])/5)
  theta3=((x[1]+x[2]+x[3]+x[4])/4)
  theta4=(min(x) + max(x))/2
  theta=c(theta1,theta2,theta3,theta4)
  return(theta)
}
resultados <- t(sapply(1:50, function(x) estimadores()))


resultados <- as.data.frame(resultados)
colnames(resultados) <- c("theta1", "theta2", "theta3", "theta4")

Histogramas

par(mfrow = c(2, 2)) # Dividir el área de la trama en 2x2
library(ggplot2)
for (j in 1:4) {
  hist(resultados[, j], main = paste("Estimador", j), xlab = "Valor estimado")
}

library(ggplot2)


boxplot(resultados, las=1, main="Comparación estimadores con n=50")  # gráfico de comparación   
abline(h=1,  col="red")     

apply(resultados,2,mean)   
##    theta1    theta2    theta3    theta4 
## 0.9942197 1.9986551 0.9972690 1.1887274
apply(resultados,2,var) 
##    theta1    theta2    theta3    theta4 
## 0.2630406 1.1820158 0.2299353 0.3595139
apply(resultados,2,sd)
##    theta1    theta2    theta3    theta4 
## 0.5128748 1.0872055 0.4795157 0.5995948

Caso 3, muestra de n=100

Con la muestra aleatoria de tamaño \(n=4\) , ( \(X_{1}, X_{2}, X_{3}, X_{4}\)), cuya población la conforma una distribución exponencial con parámetro \(\theta =1\)

estimadores=function(theta){
  x=rexp(4,rate = 1)
  theta1=((x[1]+x[2])/6)+((x[3]+x[4])/3)
  theta2=((x[1]+2*x[2]+3*x[3]+4*x[4])/5)
  theta3=((x[1]+x[2]+x[3]+x[4])/4)
  theta4=(min(x) + max(x))/2
  theta=c(theta1,theta2,theta3,theta4)
  return(theta)
}
resultados <- t(sapply(1:100, function(x) estimadores()))


resultados <- as.data.frame(resultados)
colnames(resultados) <- c("theta1", "theta2", "theta3", "theta4")

Histogramas

par(mfrow = c(2, 2)) # Dividir el área de la trama en 2x2
library(ggplot2)
for (j in 1:4) {
  hist(resultados[, j], main = paste("Estimador", j), xlab = "Valor estimado")
}

library(ggplot2)


boxplot(resultados, las=1, main="Comparación estimadores con n=100")  # gráfico de comparación   
abline(h=1,  col="red")     

apply(resultados,2,mean)   
##    theta1    theta2    theta3    theta4 
## 0.9287932 1.8588598 0.9534211 1.1129014
apply(resultados,2,var) 
##    theta1    theta2    theta3    theta4 
## 0.2541739 1.1789399 0.2387615 0.3778919
apply(resultados,2,sd)
##    theta1    theta2    theta3    theta4 
## 0.5041566 1.0857900 0.4886323 0.6147291

Caso 4, muestra de n=1000

Con la muestra aleatoria de tamaño \(n=4\) , ( \(X_{1}, X_{2}, X_{3}, X_{4}\)), cuya población la conforma una distribución exponencial con parámetro \(\theta =1\)

estimadores=function(theta){
  x=rexp(4,rate = 1)
  theta1=((x[1]+x[2])/6)+((x[3]+x[4])/3)
  theta2=((x[1]+2*x[2]+3*x[3]+4*x[4])/5)
  theta3=((x[1]+x[2]+x[3]+x[4])/4)
  theta4=(min(x) + max(x))/2
  theta=c(theta1,theta2,theta3,theta4)
  return(theta)
}
resultados <- t(sapply(1:1000, function(x) estimadores()))


resultados <- as.data.frame(resultados)
colnames(resultados) <- c("theta1", "theta2", "theta3", "theta4")

Histogramas

par(mfrow = c(2, 2)) # Dividir el área de la trama en 2x2
library(ggplot2)
for (j in 1:4) {
  hist(resultados[, j], main = paste("Estimador", j), xlab = "Valor estimado")
}

library(ggplot2)


boxplot(resultados, las=1, main="Comparación estimadores con n=1000")  # gráfico de comparación   
abline(h=1,  col="red") 

apply(resultados,2,mean)   
##    theta1    theta2    theta3    theta4 
## 0.9812648 1.9588880 0.9885603 1.1510213
apply(resultados,2,var) 
##    theta1    theta2    theta3    theta4 
## 0.2737777 1.1777004 0.2538531 0.3950889
apply(resultados,2,sd)
##    theta1    theta2    theta3    theta4 
## 0.5232377 1.0852191 0.5038384 0.6285610

Resultados

A continuación se presentan los resultados obtenidos en las simulaciones

Tabla de resultados simulación

** Muestra n=20**

Estimador Media Varianza
\(\hat{\theta}_{1}\) 0.934147 0.3164758
\(\hat{\theta}_{2}\) 1.911135 1.4492040
\(\hat{\theta}_{3}\) 1.007984 0.4020187
\(\hat{\theta}_{4}\) 1.331279 1.0802144

De acuerdo con los resultados se observa que \(\hat{\theta}_{3}\) tiene la media más cercana al valor esperado y \(\hat{\theta}_{1}\) tiene la menor varianza.

** Muestra n=50**

Estimador Media Varianza
\(\hat{\theta}_{1}\) 0.9942197 0.2630406
\(\hat{\theta}_{2}\) 1.9986551 1.1820158
\(\hat{\theta}_{3}\) 0.9972690 0.2299353
\(\hat{\theta}_{4}\) 1.1887274 0.3595139

De acuerdo con los resultados se observa que \(\hat{\theta}_{3}\) tiene la media más cercana al valor esperado y la menor varianza.

** Muestra n=100**

Estimador Media Varianza
\(\hat{\theta}_{1}\) 0.9287932 0.2541739
\(\hat{\theta}_{2}\) 1.8588598 1.1789399
\(\hat{\theta}_{3}\) 0.9534211 0.2387615
\(\hat{\theta}_{4}\) 1.1129014 0.3778919

De acuerdo con los resultados se observa que \(\hat{\theta}_{3}\) tiene la media más cercana al valor esperado y la menor varianza.

** Muestra n=1000**

Estimador Media Varianza
\(\hat{\theta}_{1}\) 0.9812648 0.2737777
\(\hat{\theta}_{2}\) 1.9588880 1.1777004
\(\hat{\theta}_{3}\) 0.9885603 0.2538531
\(\hat{\theta}_{4}\) 1.1510213 0.3950889

De acuerdo con los resultados se observa que \(\hat{\theta}_{3}\) tiene la media más cercana al valor esperado y la menor varianza.

Conclusiones

Los resultados de la simulación con diferentes tamaños de muestras sugieren que el estimador \(\hat{\theta}_{3}\) exhibe una notable característica de ser insesgado y eficiente, ya su media es muy cercana al valor esperado, su varianza es la menor entre los estimadores y disminuye a medida que aumenta el tamaño de la muestra.

Por otro lado, los estimadores \(\hat{\theta}_{1}\) y \(\hat{\theta}_{4}\) , aunque muestran un ligero sesgo mínimo, tienden a converger hacia el valor esperado a medida que el tamaño de la muestra aumenta, lo que sugiere que son consistentes.

Sin embargo, el estimador \(\hat{\theta}_{2}\) no cumple con las propiedades deseables de un estimador, ya que no muestra signos de convergencia hacia el valor esperado a medida que aumenta el tamaño de la muestra y su varianza es alta comparada con la varianza de los otros estimadores.

Además, al analizar la forma de los diagramas de cajas asociados con cada estimador, se observa que la distribución de los datos puede influir en la eficacia y confiabilidad de los estimadores. Esta información resalta la importancia de considerar no solo el sesgo y la eficiencia, sino también la forma de los datos al evaluar la idoneidad de los estimadores.

Problema 3. Teorema del límite central

El Teorema del Límite Central representa uno de los principios fundamentales en la inferencia estadística. Este teorema establece que, a medida que el tamaño de la muestra aumenta, la distribución de la media muestral tiende a converger hacia una distribución normal. Esta convergencia se produce independientemente de la forma de la distribución original. Según varios autores, esta aproximación es considerada confiable cuando el tamaño muestral supera un umbral típico, generalmente establecido en \(n>30\)

A continuación se realiza una simulación con el fin de verificar el teorema del límite central:

Paso 1

Generar una población de 1000 individuos, donde el porcentaje de individuos enfermos es del 50%.

set.seed(123)
total<- 1000 
enfermos<- floor(total*0.5) 
sanos <- total - enfermos 
Poblacion <- c(rep(1, enfermos), rep(0, sanos)) 
lote_total <- sample(Poblacion)

Paso 2

Generar una función que permita obtener una muestra aleatoria de la población de tamaño n=50 y calcular el estimador de la proporción muestral \(\hat{p}\) para un tamaño de muestra dado n.

muestra_n<- function(Poblacion, n){
  
muestra<- sample(Poblacion, size = n)  
proporcion_muestral <-sum(muestra==1)/ n 
return(proporcion_muestral)

}

tamaño<- 50
proporcion<- muestra_n(lote_total, tamaño)
cat("El estimador p^ de la proporción muestral es:", proporcion, "\n")
## El estimador p^ de la proporción muestral es: 0.48

Paso 3

Repitir el paso anterior \(500\) veces y analizar los resultados del estimador \(\hat{p}\) para el total de casos.

tamaño <- 50  
repeticiones <- 500  

proporciones <- numeric(repeticiones)  

for (i in 1:repeticiones) {
  proporciones[i] <- muestra_n(lote_total, tamaño)
}
summary(proporciones)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  0.2400  0.4600  0.5000  0.5011  0.5400  0.6800
hist(proporciones,  main = "n=500", freq=FALSE)

library(e1071)

library(moments)
## 
## Attaching package: 'moments'
## The following objects are masked from 'package:e1071':
## 
##     kurtosis, moment, skewness
media <- mean(proporciones)
mediana <- median(proporciones)
desviacion <- sd(proporciones)
varianza <- var(proporciones)
asimetria <- skewness(proporciones)
curtosis <- kurtosis(proporciones)
data.frame(media, mediana, varianza, desviacion, asimetria, curtosis)
shapiro.test(proporciones)
## 
##  Shapiro-Wilk normality test
## 
## data:  proporciones
## W = 0.98921, p-value = 0.0009728

Análisis de los datos obtenidos

  1. Simetría y sesgo: La media y la mediana son cercanas entre sí (0.50112 y 0.5 respectivamente), lo que sugiere una distribución relativamente simétrica. La medida de simetría es ligeramente negativa (-0.1437757), lo que indica una ligera asimetría hacia la izquierda, aunque esta asimetría no es significativa. La curtosis es 3.339249, lo que sugiere una distribución leptocúrtica, es decir, con colas más pesadas y picos más altos que la distribución normal.
  2. Variabilidad: La varianza (0.004503753) y la desviación estándar (0.06711001) son relativamente pequeñas, lo que sugiere que los datos tienen poca dispersión alrededor de la media. Los resultados son consistentes y tienen poca variabilidad entre los estimadores.
  3. Normalidad: El resultado del test de Shapiro-Wilk muestra un valor p muy bajo (0.0009728), lo que indica que los datos no se ajustan bien a una distribución normal.

En resumen, los resultados sugieren que los datos tienen una distribución relativamente simétrica, con una ligera asimetría hacia la izquierda, y una variabilidad baja. Sin embargo, los datos no se ajustan bien a una distribución normal, lo que puede tener implicaciones en la elección de los métodos estadísticos apropiados para analizar estos datos.

Paso 4

Repita los pasos 2 y 3 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()).

# Define la función para realizar la simulación de proporciones muestrales
simulacion_proporciones <- function(tamaños, repeticiones, lote_total) {
  resultados <- matrix(nrow = repeticiones, ncol = length(tamaños))  # Matriz para almacenar los resultados
  
  for (i in seq_along(tamaños)) {
    tamaño <- tamaños[i]
    proporciones <- replicate(repeticiones, muestra_n(lote_total, tamaño))
    resultados[, i] <- proporciones
  }
  
  colnames(resultados) <- tamaños  # Etiquetar las columnas con los tamaños de muestra
  return(resultados)
}

# Función para realizar pruebas de normalidad y gráficos QQ
pruebas_normalidad <- function(resultados_simulacion) {
  for (i in seq_along(colnames(resultados_simulacion))) {
    tamaño <- colnames(resultados_simulacion)[i]
    resultados_p_n <- resultados_simulacion[, i]
    
    # Test de Shapiro-Wilk
    shapiro_result <- shapiro.test(resultados_p_n)
    cat("Para una muestra de tamaño:", tamaño, "\n")
    cat("Test de Shapiro-Wilk:\n")
    print(shapiro_result)
 
   
    # Gráfico QQ-Normal
    qqnorm(resultados_p_n, main = paste("QQ-Normal para tamaño:", tamaño))
    qqline(resultados_p_n)
    
    
  }
}

# Parámetros
tamaño_n <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)
repeticiones <- 500  

# Llamada a la función para realizar la simulación
resultados_simulacion <- simulacion_proporciones(tamaño_n, repeticiones, lote_total)

# Realizar pruebas de normalidad y gráficos QQ
pruebas_normalidad(resultados_simulacion)
## Para una muestra de tamaño: 5 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.93117, p-value = 2.12e-14

## Para una muestra de tamaño: 10 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.95783, p-value = 9.13e-11

## Para una muestra de tamaño: 15 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.97469, p-value = 1.315e-07

## Para una muestra de tamaño: 20 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.97258, p-value = 4.627e-08

## Para una muestra de tamaño: 30 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.98643, p-value = 0.0001311

## Para una muestra de tamaño: 50 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99213, p-value = 0.009596

## Para una muestra de tamaño: 60 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99164, p-value = 0.006483

## Para una muestra de tamaño: 100 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.9935, p-value = 0.02996

## Para una muestra de tamaño: 200 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99165, p-value = 0.006522

## Para una muestra de tamaño: 500 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99417, p-value = 0.05251

Paso 5

Repita toda la simulación, 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.

set.seed(123)
total<- 1000 
enfermos<- floor(total*0.1) 
sanos <- total - enfermos 
Poblacion <- c(rep(1, enfermos), rep(0, sanos)) 
lote_total <- sample(Poblacion) 
# Define la función para realizar la simulación de proporciones muestrales
simulacion_proporciones <- function(tamaños, repeticiones, lote_total) {
  resultados <- matrix(nrow = repeticiones, ncol = length(tamaños))  # Matriz para almacenar los resultados
  
  for (i in seq_along(tamaños)) {
    tamaño <- tamaños[i]
    proporciones <- replicate(repeticiones, muestra_n(lote_total, tamaño))
    resultados[, i] <- proporciones
  }
  
  colnames(resultados) <- tamaños  # Etiquetar las columnas con los tamaños de muestra
  return(resultados)
}

# Función para realizar pruebas de normalidad y gráficos QQ
pruebas_normalidad <- function(resultados_simulacion) {
  for (i in seq_along(colnames(resultados_simulacion))) {
    tamaño <- colnames(resultados_simulacion)[i]
    resultados_p_n <- resultados_simulacion[, i]
    
    # Test de Shapiro-Wilk
    shapiro_result <- shapiro.test(resultados_p_n)
    cat("Para una muestra de tamaño:", tamaño, "\n")
    cat("Test de Shapiro-Wilk:\n")
    print(shapiro_result)
    
    # Gráfico QQ-Normal
    qqnorm(resultados_p_n, main = paste("QQ-Normal para tamaño:", tamaño))
    qqline(resultados_p_n)
   
  }
}

# Parámetros
tamaño_n <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)
repeticiones <- 500  

# Llamada a la función para realizar la simulación
resultados_simulacion <- simulacion_proporciones(tamaño_n, repeticiones, lote_total)

# Realizar pruebas de normalidad y gráficos QQ
pruebas_normalidad(resultados_simulacion)
## Para una muestra de tamaño: 5 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.73821, p-value < 2.2e-16

## Para una muestra de tamaño: 10 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.83393, p-value < 2.2e-16

## Para una muestra de tamaño: 15 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.90171, p-value < 2.2e-16

## Para una muestra de tamaño: 20 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.92715, p-value = 7.372e-15

## Para una muestra de tamaño: 30 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.95115, p-value = 8.537e-12

## Para una muestra de tamaño: 50 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.97657, p-value = 3.476e-07

## Para una muestra de tamaño: 60 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.97736, p-value = 5.305e-07

## Para una muestra de tamaño: 100 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.98171, p-value = 6.292e-06

## Para una muestra de tamaño: 200 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99128, p-value = 0.004844

## Para una muestra de tamaño: 500 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99089, p-value = 0.003559

set.seed(123)
total<- 1000 
enfermos<- floor(total*0.9) 
sanos <- total - enfermos 
Poblacion <- c(rep(1, enfermos), rep(0, sanos)) 
lote_total <- sample(Poblacion) 
# Define la función para realizar la simulación de proporciones muestrales
simulacion_proporciones <- function(tamaños, repeticiones, lote_total) {
  resultados <- matrix(nrow = repeticiones, ncol = length(tamaños))  # Matriz para almacenar los resultados
  
  for (i in seq_along(tamaños)) {
    tamaño <- tamaños[i]
    proporciones <- replicate(repeticiones, muestra_n(lote_total, tamaño))
    resultados[, i] <- proporciones
  }
  
  colnames(resultados) <- tamaños  # Etiquetar las columnas con los tamaños de muestra
  return(resultados)
}

# Función para realizar pruebas de normalidad y gráficos QQ
pruebas_normalidad <- function(resultados_simulacion) {
  for (i in seq_along(colnames(resultados_simulacion))) {
    tamaño <- colnames(resultados_simulacion)[i]
    resultados_p_n <- resultados_simulacion[, i]
    
    # Test de Shapiro-Wilk
    shapiro_result <- shapiro.test(resultados_p_n)
    cat("Para una muestra de tamaño:", tamaño, "\n")
    cat("Test de Shapiro-Wilk:\n")
    print(shapiro_result)
    
    # Gráfico QQ-Normal
    qqnorm(resultados_p_n, main = paste("QQ-Normal para tamaño:", tamaño))
    qqline(resultados_p_n)
    
    
  }
}

# Parámetros
tamaño_n <- c(5, 10, 15, 20, 30, 50, 60, 100, 200, 500)
repeticiones <- 500  

# Llamada a la función para realizar la simulación
resultados_simulacion <- simulacion_proporciones(tamaño_n, repeticiones, lote_total)

# Realizar pruebas de normalidad y gráficos QQ
pruebas_normalidad(resultados_simulacion)
## Para una muestra de tamaño: 5 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.6933, p-value < 2.2e-16

## Para una muestra de tamaño: 10 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.82451, p-value < 2.2e-16

## Para una muestra de tamaño: 15 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.89293, p-value < 2.2e-16

## Para una muestra de tamaño: 20 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.92519, p-value = 4.464e-15

## Para una muestra de tamaño: 30 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.95614, p-value = 4.908e-11

## Para una muestra de tamaño: 50 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.96953, p-value = 1.106e-08

## Para una muestra de tamaño: 60 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.97705, p-value = 4.478e-07

## Para una muestra de tamaño: 100 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.98631, p-value = 0.0001207

## Para una muestra de tamaño: 200 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99242, p-value = 0.01224

## Para una muestra de tamaño: 500 
## Test de Shapiro-Wilk:
## 
##  Shapiro-Wilk normality test
## 
## data:  resultados_p_n
## W = 0.99458, p-value = 0.07447

Resultados

Comparación de normalidad para los diferentes tamaños de muestras

El Test de Shapiro-Wilk es una prueba estadística utilizada para determinar si un conjunto de datos sigue una distribución normal. La hipótesis nula \(H_{0}\) de la prueba de Shapiro-Wilk es que los datos siguen una distribución normal. La hipótesis alternativa \(H_{1}\) es que los datos no siguen una distribución normal.

Los siguientes fueron los resultados obtenidos para una población de 1000 con el 50% de la población enferma

Tamaño <- c(5,10,15, 20,30,50,60,100,200,500)
w <- c(0.93117,0.95783, 0.97469 , 0.97258, 0.98643 , 0.99213, 0.99164 ,0.9935, 0.99165, 0.99417)
p <-c(2.12e-14, 9.13e-11 , 1.315e-07,  4.627e-08, 0.0001311 , 0.009596, 0.006483, 0.02996 , 0.006522 ,0.05251 )
Normalidad<- c("no","no","no","no","no","no","no","no","no","si")
resultados1<-data.frame(Tamaño, w, p, Normalidad)
resultados1

Los siguientes fueron los resultados obtenidos para una población de 1000 con el 10% de la población enferma

Tamaño <- c(5,10,15, 20,30,50,60,100,200,500)
w <- c(0.73821,0.83393, 0.90171 ,0.92715,0.95115 ,  0.97657, 0.97736,0.98171, 0.99128, 0.99089 )
p <-c(2.2e-16 ,2.2e-16 , 2.2e-16 ,7.372e-15 , 8.537e-12,3.476e-07 , 5.305e-07 , 6.292e-06 , 0.004844 ,0.003559 )
Normalidad<- c("no","no","no","no","no","no","no","no","no","no")
resultados2<-data.frame(Tamaño, w, p, Normalidad)
resultados2

Los siguientes fueron los resultados obtenidos para una población de 1000 con el 90% de la población enferma

Tamaño <- c(5,10,15, 20,30,50,60,100,200,500)
w <- c(0.6933, 0.82451, 0.89293, 0.92519, 0.95614, 0.96953 , 0.97705 , 0.98631, 0.99242, 0.99458 )
p <-c(2.2e-16, 2.2e-16, 2.2e-16 ,4.464e-15, 4.908e-11 , 1.106e-08, 4.478e-07, 0.0001207 , 0.01224, 0.07447 )
Normalidad<- c("no","no","no","no","no","no","no","no","no","si")
resultados3<-data.frame(Tamaño, w, p, Normalidad)
resultados3

Conclusiones

En resumen, al analizar los resultados de la simulación, se observa que los p-valores de la prueba de Shapiro-Wilk y la apariencia de los gráficos Q-Q cambian significativamente a medida que aumenta el tamaño de la muestra. Este comportamiento refleja una tendencia hacia la normalidad, lo cual es consistente con el Teorema del Límite Central.

Al considerar lotes con diferentes porcentajes de plantas enfermas (10%, 50% y 90%), los resultados del análisis muestran que la normalidad de la distribución del estimador de la proporción muestral está influenciada por la variabilidad en la población. Tanto una alta como una baja proporción de plantas enfermas tienden a alejar la distribución de la normalidad, mientras que una proporción más equitativa se acerca más a una distribución normal.

Para tamaños de muestra pequeños, como 5, 10 y 15, los valores p de las pruebas de normalidad son extremadamente bajos, lo que indica una falta de ajuste a una distribución normal. A medida que el tamaño de muestra aumenta, la evidencia en contra de la normalidad disminuye, aunque para tamaños de muestra más grandes, como 500, aún puede persistir una falta de ajuste perfecto a la normalidad.

Problema 4. Estimación boostrap

Cuando se extrae una muestra de una población que no es normal y se requiere estimar un intervalo de confianza se pueden utilizar los métodos de estimación bootstrap. Esta metodología supone que se puede reconstruir la población objeto de estudio mediante un muestreo con reemplazo de la muestra que se tiene. Existen varias versiones del método. Una presentación básica del método se describe a continuación:

El artículo de In-use Emissions from Heavy Duty Dissel Vehicles (J.Yanowitz, 2001) presenta las mediciones de eficiencia de combustible en millas/galón de una muestra de siete camiones. Los datos obtenidos son los siguientes: 7.69, 4.97, 4.56, 6.49, 4.34, 6.24 y 4.45. Se supone que es una muestra aleatoria de camiones y que se desea construir un intervalo de confianza del 95 % para la media de la eficiencia de combustible de esta población. No se tiene información de la distribución de los datos. El método bootstrap permite construir intervalos de confianza del 95 % - Para ilustrar el método suponga que coloca los valores de la muestra en una caja y extrae uno al azar. Este correspondería al primer valor de la muestra bootstrap X1. Después de anotado el valor se regresa X1 a la caja y se extrae el valor X2, regresandolo nuevamente. Este procedimiento se repite hasta completar una muestra de amaño n, X1, X2, X2, Xn, conformando la muestra bootstrap.

Es necesario extraer un gran número de muestras (suponga k = 1000). Para cada una de las muestra bootstrap obtenidas se calcula la media Xi, obteniéndose un valor para cada muestra. El intervalo de confianza queda conformado por los percentiles P2.5 y P97.5. Existen dos métodos para estimarlo:

  • Método 1 (P2.5;P97.5)
  • Método 2 (2X−P97.5;2X−P2.5)

Construya el intervalo de confianza por los dos métodos y compare los resultados obtenidos. Comente los resultados. Confiaría en estas estimaciones

# Cargar información inicial
datos <- c(7.69, 4.97, 4.56, 6.49, 4.34, 6.24, 4.45)
k <- 1000

# Vector para almacenar las medias de las muestras bootstrap
medias_bootstrap <- vector("double", length = k)

# Realizar el procedimiento bootstrap
for (i in 1:k) {
  muestra_bootstrap <- sample(datos, replace = TRUE)
  medias_bootstrap[i] <- mean(muestra_bootstrap)
}

# Cálculo de los percentiles
percentil_025 <- quantile(medias_bootstrap, 0.025)
percentil_975 <- quantile(medias_bootstrap, 0.975)

# Cálculo del intervalo de confianza
intervalo_confianza_1 <- c(percentil_025, percentil_975)
intervalo_confianza_2 <- c(2 * mean(datos) - percentil_975, 2 * mean(datos) - percentil_025)

# Impresión de los resultados
cat("Media de la muestra original:", mean(datos), "\n")
## Media de la muestra original: 5.534286
cat("Percentil 2.5%:", percentil_025, "\n")
## Percentil 2.5%: 4.694179
cat("Percentil 97.5%:", percentil_975, "\n")
## Percentil 97.5%: 6.417143
cat("Intervalo de confianza metodo 1:", intervalo_confianza_1, "\n")
## Intervalo de confianza metodo 1: 4.694179 6.417143
cat("Intervalo de confianza metodo 2:", intervalo_confianza_2, "\n")
## Intervalo de confianza metodo 2: 4.651429 6.374393

Comentarios

  • Los intervalos de confianza son diferentes para ambos métodos. Esto se debe a que el método 2 utiliza la media de la muestra original como referencia, mientras que el método 1 no.
  • A pesar de que se evidencia que el intervalo de confianza del método 1 es más grande que el del método 2, la diferencia entre los intervalos no es muy grande, lo que indica que ambos métodos son relativamente robustos.

Confianza en las estimaciones

  • La confianza en las estimaciones depende de varios factores, como el tamaño de la muestra, la representatividad de la muestra y la distribución de los datos. En este caso, la muestra es pequeña y no se tiene información sobre la distribución de los datos. Por lo tanto, la confianza en las estimaciones es limitada. Sin embargo, el método bootstrap proporciona una forma de estimar la incertidumbre en las estimaciones, incluso cuando no se conoce la distribución de los datos.

5. Referencias