Metodos y Simulación Estadistica
Actividad 2. Ejercicios de Simulación
Gustavo Moreno Collazos y Nicolás García Peñaloza
2024-03-12
Introducción
“El proceso de simulación constituye una herramienta poderosa para la estadística que se pueden emplear para entender relaciones complejas y estimar valores difíciles de calcular directamente. Para entenderlo utilizaremos se plantean los siguientes problemas:”
En el presente documento se encuentra la solución propuesta para los problemas de simulación estadística correspondientes a la actividad número 2 del curso de simulación y métodos estadística de la maestría en Ciencia de Datos de la Pontificia Universidad Javeriana.
Problema 1
Estimación del valor de \(\pi\) La siguiente figura sugiere como estimar el valor de \(\pi\) con una simulación. En la figura, un circuito con un área igual a \(\pi/4\), está inscrito en un cuadrado cuya área es igual a 1. Se elige de forma aleatoria n puntos dentro del cuadrado . La probabilidad de que un punto esté dentro del círculo es igual a la fracción del área del cuadrado que abarca a este, la cual es \(\pi/4\) Por tanto, se puede estimar el valor de \(\pi/4\) al contar el número de puntos dentro del círculo, para obtener la estimación de \(\pi/4\) De este último resultado se encontrar una aproximación para el valor de \(\pi\)
figura 1
Pasos sugeridos:
- Genere \(n\) coordenadas \(x: X_{1}, . . . , X_{n}\). Utilice la distribución uniforme con valor mínimo de \(0\) y valor máximo de \(1\) La distribución uniforme genera variables aleatorias que tienen la misma probabilidad de venir de cualquier parte del intervalo \((0,1)\).
# Definir la función para determinar si un punto está dentro del círculo
esta_dentro_del_circulo = function(coordenada_x, coordenada_y) {
distancia_desde_centro = (coordenada_x - 0.5)^2 + (coordenada_y - 0.5)^2
return(distancia_desde_centro < 0.25)
}Ejercicio para una simulación de n = 1.000
Genere \(1000\) coordenadas \(y: Y_{1},...,Y_{n}\), utilizando nuevamente la distribución uniforme con valor mínimo de \(0\) y valor máximo de \(1\).
Cada punto \((X_{i},Y_{i})\) se encuentra dentro del círculo si su distancia desde el centro \((0.5,0.5)\) es menor a \(0.5\). Para cada par \((X_{i},Y_{i})\) determine si la distancia desde el centro es menor a \(0.5\). Esto último se puede realizar al calcular el valor \((X_{i}−0.5)^{2}+(Y_{i}−0.5)^{2}\), que es el cuadrado de la distancia, y al determinar si es menor que \(0.25\).
¿Cuántos de los puntos están dentro del círculo? ¿Cuál es su estimación de \(\pi\)?
n = 1000
# Generar n coordenadas de X y Y
set.seed(2023)
coordenada_x = runif(n, min = 0, max = 1)
coordenada_y = runif(n, min = 0, max = 1)
# Determinar si los puntos están dentro del círculo
dentro_del_circulo = esta_dentro_del_circulo(coordenada_x, coordenada_y)
# Contar cuántos puntos están dentro del círculo
puntos_dentro = sum(dentro_del_circulo)
# Calcular la estimación de π
estimacion_pi = 4 * (puntos_dentro / n)
cat("Número de puntos dentro del círculo para n =", n, ":", puntos_dentro, "\n")## Número de puntos dentro del círculo para n = 1000 : 801cat("Estimación de π para n =", n, ":", estimacion_pi, "\n")## Estimación de π para n = 1000 : 3.204# Crear un gráfico n = 1000
plot(coordenada_x, coordenada_y, col = ifelse(dentro_del_circulo, "#63B8FF", "#EEE9E9"), pch = 20,
main = "Simulación de puntos aleatorios n = 1.000",
xlab = "Coordenada X", ylab = "Coordenada Y")
# Dibujar el círculo sin borrar los puntos existentes
theta = seq(0, 2 * pi, length.out = 100)
circ_x = 0.5 + 0.5 * cos(theta)
circ_y = 0.5 + 0.5 * sin(theta)
lines(circ_x, circ_y, col = "#0000CD", lwd = 2)
Ejercicio para una simulación de n = 10.000
n = 10000
set.seed(2030)
coordenada_x = runif(n, min = 0, max = 1)
coordenada_y = runif(n, min = 0, max = 1)
dentro_del_circulo = esta_dentro_del_circulo(coordenada_x, coordenada_y)
puntos_dentro = sum(dentro_del_circulo)
estimacion_pi = 4 * (puntos_dentro / n)
print(paste("Número de puntos dentro del círculo para n =", n, ":", puntos_dentro))## [1] "Número de puntos dentro del círculo para n = 10000 : 7792"print(paste("Estimación de π para n =", n, ":", estimacion_pi))## [1] "Estimación de π para n = 10000 : 3.1168"# Crear un gráfico n = 10.000
plot(coordenada_x, coordenada_y, col = ifelse(dentro_del_circulo, "#63B8FF", "#EEE9E9"), pch = 20,
main = "Simulación de puntos aleatorios n = 10.000",
xlab = "Coordenada X", ylab = "Coordenada Y")
# Dibujar el círculo
theta = seq(0, 2 * pi, length.out = 100)
circ_x = 0.5 + 0.5 * cos(theta)
circ_y = 0.5 + 0.5 * sin(theta)
lines(circ_x, circ_y, col = "#0000CD", lwd = 2)
Ejercicio para una simulación de n = 100.000
n = 100000
set.seed(2030)
coordenada_x = runif(n, min = 0, max = 1)
coordenada_y = runif(n, min = 0, max = 1)
dentro_del_circulo = esta_dentro_del_circulo(coordenada_x, coordenada_y)
puntos_dentro = sum(dentro_del_circulo)
estimacion_pi = 4 * (puntos_dentro / n)
print(paste("Número de puntos dentro del círculo para n =", n, ":", puntos_dentro))## [1] "Número de puntos dentro del círculo para n = 1e+05 : 78624"print(paste("Estimación de π para n =", n, ":", estimacion_pi))## [1] "Estimación de π para n = 1e+05 : 3.14496"# Crear un gráfico n = 100.000
plot(coordenada_x, coordenada_y, col = ifelse(dentro_del_circulo, "#63B8FF", "#EEE9E9"), pch = 20,
main = "Simulación de puntos aleatorios n = 100.000",
xlab = "Coordenada X", ylab = "Coordenada Y")
# Dibujar el círculo
theta = seq(0, 2 * pi, length.out = 100)
circ_x = 0.5 + 0.5 * cos(theta)
circ_y = 0.5 + 0.5 * sin(theta)
lines(circ_x, circ_y, col = "#0000CD", lwd = 2)
Conclusión:
Al realizar un mayor numero de simulaciones, aumentando el tamaño de la muestra (\(n\)), se obtienen estimaciones más precisas y confiables para el fenómeno del valor esperado, en este caso \(\pi\) (3.14159). Esto se debe a que se reduce el error de muestreo y se mejora la representatividad de la muestra.
Problema 2
Propiedades de los estimadores La simulación ayuda a entender y validad las propiedades de los estimadores estadísticos como son. insesgadez, eficiencia y la consistencia principalmente. El siguiente problema permite evidenciar las principales características de un grupo de estimadores propuestos para la estimación de un parámetro asociado a un modelo de probabilidad.
Sean \(X_{1}, X_{2}, X_{3}\) y \(X_{4}\), una muestra aleatoria de tamaño \(n=4\) cuya población la conforma una distribución exponencial con parámetro desconocido. Determine las características de cada uno de los siguientes estimadores propuestos:
\(\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}\)
- Genere una muestras de \(n=20\), \(50\), \(100\) y \(1000\) para cada uno de los estimadores planteados.
Funcion_estimadores = function(estimador, muestra, lambda) {
# Crear matriz de estimadores
mx = matrix(rexp(estimador * muestra, rate = lambda), nrow = estimador, ncol = muestra, byrow = TRUE)
# Seleccionar subconjunto de estimadores
m_estimadores = mx[1:estimador, ]
# Calcular estimadores para cada muestra
t1 = (m_estimadores[, 1] + m_estimadores[, 2] / 6) + (m_estimadores[, 3] + m_estimadores[, 4] / 3)
t2 = (m_estimadores[, 1] + (2 * m_estimadores[, 2]) + (3 * m_estimadores[, 3]) + (4 * m_estimadores[, 4])) / 5
t3 = ((m_estimadores[, 1] + m_estimadores[, 2] + m_estimadores[, 3] + m_estimadores[, 4]) / 4)
t4 = (min(m_estimadores[, 1], m_estimadores[, 2], m_estimadores[, 3], m_estimadores[, 4]) + max(m_estimadores[, 1], m_estimadores[, 2], m_estimadores[, 3], m_estimadores[, 4])) / 2
# Graficar boxplot
boxplot(m_estimadores, xlab = "Estimadores", ylab = "Valores")
abline(h = lambda, col = "red")
# Calcular propiedades de los estimadores
m_medias = apply(m_estimadores, 2, mean)
m_varianza = apply(m_estimadores, 2, var)
sesgo = lambda - m_medias
# Imprimir propiedades de los estimadores
valores = data.frame(Estimador = paste("Estimador", 1:4), Media = m_medias, Varianza = m_varianza, Sesgo = sesgo)
print(valores)
}En cada caso evalue las propiedades de insesgadez, eficiencia y consistencia
Suponga un valor para el parámetro \(\theta\)
Introducción.
A continuación, se evaluarán las propiedades de insesgadez, eficiencia y consistencia de los estimadores. - Entendemos por insesgadez a la diferencia entre la media de un estimador y el valor del parámetro. Un estimador es insesgado si su sesgo es cercano a cero. Sesgo \(E[\theta_{1}]\)=\(\theta\) - La eficiencia por su parte, está relacionada con la varianza del estimador. El cual considera que entre menor sea su varianza, es mas eficiente. \(V[\hat{\theta_{2}}]\)<\(V[\hat{\theta_{2}}]\) - Y la consistencia está relacionada con el comportamiento del estimador a medida que la muestra tiende al infinito. Es decir, un estimador es consistente si converge al verdadero valor del parámetro a medida que el tamaño de la muestra aumenta.
Simulación con n = 20
Respecto a la primera muestra, con n = 20, podemos observar que tanto el estimador 1 como el estimador 4, tienen sesgos muy cercanos a cero. Por lo que ambos pueden considerarse insesgados. Sin embargo, el estimador 4 posee el valor más cercano a cero (-0.0669). En temas de eficiencia el estimador 4 tiene una varianza (0.638) más baja que los demás estimadores, lo que indica ser el más eficiente en términos de dispersión de los datos alrededor de la media.
# Llamar a la función con diferentes parámetros
lambda = 1
mue = 4
sim = 20
set.seed(1015)
Funcion_estimadores(sim, mue, lambda)
## Estimador Media Varianza Sesgo
## 1 Estimador 1 1.099279 0.8542797 -0.09927873
## 2 Estimador 2 1.349316 2.8253037 -0.34931568
## 3 Estimador 3 1.322809 3.1608728 -0.32280896
## 4 Estimador 4 1.066961 0.6389195 -0.06696082Simulación con n = 50
Para la siguiente muestra, con n = 50, podemos observar ahora que tanto el estimador 1, 3 y 4 tienen sesgos muy cercanos a cero, considerándolos como insesgados; se conserva la constante del estimador 4 como el más insesgado (0.0508). En temas de eficiencia el estimador 4 conserva la menor varianza posible entre los cuatro estimadores (0.47) seguido del estimador 1 (0.91).
lambda = 1
mue = 4
sim = 50
set.seed(1015)
Funcion_estimadores(sim, mue, lambda)
## Estimador Media Varianza Sesgo
## 1 Estimador 1 1.0965000 0.9183035 -0.09650002
## 2 Estimador 2 1.1346336 1.6719717 -0.13463364
## 3 Estimador 3 1.0892831 1.6662160 -0.08928305
## 4 Estimador 4 0.9491779 0.4712841 0.05082210Simulación con n = 100
Para la siguiente muestra, con n = 100, podemos observar la consistencia del estimador 4 al permanecer como el estimador de mayor insesgadez entre los 4 estimadores evaluados. De igual forma, su eficiencia es la mejor entre los demás con un valor de 0.73 seguido del estimador 2 con un valor de 1.05.
lambda = 1
mue = 4
sim = 100
set.seed(1015)
Funcion_estimadores(sim, mue, lambda)
## Estimador Media Varianza Sesgo
## 1 Estimador 1 1.1214879 1.4072443 -0.121487883
## 2 Estimador 2 0.9679395 1.0553541 0.032060531
## 3 Estimador 3 1.0250448 1.2348286 -0.025044757
## 4 Estimador 4 0.9918490 0.7367481 0.008150997Simulación con n = 1000
Para la siguiente muestra, con n = 1000, se puede observar que ahora el estimador 1 es el más insesgado con un valor de 0.0019 en comparación con el estimar 4 quien le sigue con un valor de 0.012. Por su parte, la mejor confianza es para el estimador 2 cuya varianza es la inferior a los otros 3 estimadores con un valor de 0.86.
lambda = 1
mue = 4
sim = 1000
set.seed(1015)
Funcion_estimadores(sim, mue, lambda)
## Estimador Media Varianza Sesgo
## 1 Estimador 1 1.0019591 1.0430985 -0.001959072
## 2 Estimador 2 0.9653345 0.8652995 0.034665511
## 3 Estimador 3 0.9575673 0.9374961 0.042432724
## 4 Estimador 4 0.9875813 0.9374443 0.012418745Conclusión:
En conclusión, dado el aumento de los valores de la muestra podemos observar cierta constancia en las propiedades de los estimadores como insesgadez y eficiencia. Sin embargo, a pesar de que la estimación exhibió variaciones que no se alejaron significativamente de su valor inicial, se observó que la varianza podría experimentar cambios en respuesta a un aumento en los valores atípicos y tienden a perder consistencia.
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:
- 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\%\).
set.seed(1997)
LOTE = rbinom(1000 , 1 , .5)Genere una función que permita: Obtener una muestra aleatoria de la población y Calcule el estimador de la proporción muestral \(\hat{p}\) para un tamaño de muestra dado \(n\).
Repita el escenario anterior (b) \(n=500\) veces y analice los resultados en cuanto al comportamiento de los \(500\) resultados del estimador \(\hat{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.
muestra_aleatoria_aleatoria = function(vector) {
tamano_muestra = sample(1:length(vector), 1)
muestra = sample(vector, size = tamano_muestra, replace = FALSE)
return(muestra)
}
set.seed(1997)
muestra_1 = muestra_aleatoria_aleatoria(LOTE) ## Se tomo el 17% de la muestra
set.seed(1928)
LOTE_500 = rbinom(500 , 1 , .5)
set.seed(1928)
muestra_2 = muestra_aleatoria_aleatoria(LOTE_500) ## Se tomo el 26% de la muestra
prop.table(table(muestra_2))## muestra_2
## 0 1
## 0.5419847 0.45801530.5402-0.458 ;rm(LOTE ,LOTE_500 , muestra_1 , muestra_2 )## [1] 0.0822Los resultados varian en la unidad de elección 1 Mientras en la primera simulación se obtiene un 54% de los datos hacia plantas enfermas dato la seleción del \(70\%\) de la muestra. en el segundo ejercicio este resultado es el \(45\%\) de los datos una diferencia de 0.08 puntos porcentuales.
- 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 :shapiro.test()) y métodos gráficos (gráfico de normalidad: qqnorm()). Comente en su informe los resultados obtenidos
set.seed(1943)
LOTE_5 = rbinom(5 , 1 , .5)
LOTE_10 = rbinom(10 , 1 , .5)
LOTE_15 = rbinom(15 , 1 , .5)
LOTE_20 = rbinom(20 , 1 , .5)
LOTE_30 = rbinom(30 , 1 , .5)
LOTE_50 = rbinom(50 , 1 , .5)
LOTE_60 = rbinom(60 , 1 , .5)
LOTE_100 = rbinom(100 , 1 , .5)
LOTE_200 = rbinom(200 , 1 , .5)
LOTE_500 = rbinom(500 , 1 , .5)
#
set.seed(1943)
muestra_1 = muestra_aleatoria_aleatoria(LOTE_5)
muestra_2 = muestra_aleatoria_aleatoria(LOTE_10)
muestra_3 = muestra_aleatoria_aleatoria(LOTE_15)
muestra_4 = muestra_aleatoria_aleatoria(LOTE_20)
muestra_5 = muestra_aleatoria_aleatoria(LOTE_30)
muestra_6 = muestra_aleatoria_aleatoria(LOTE_50)
muestra_7 = muestra_aleatoria_aleatoria(LOTE_60)
muestra_8 = muestra_aleatoria_aleatoria(LOTE_100)
muestra_9 = muestra_aleatoria_aleatoria(LOTE_200)
muestra_10 = muestra_aleatoria_aleatoria(LOTE_500)
# Realizar la prueba de Shapiro para cada muestra en la lista
muestras = list(muestra_2, muestra_3, muestra_4, muestra_5, muestra_6, muestra_7, muestra_8, muestra_9, muestra_10)
for (muestra in muestras) {
print(shapiro.test(muestra))
print(qqnorm(muestra))
}##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.75, p-value < 2.2e-16
## $x
## [1] -0.8694238 0.0000000 0.8694238
##
## $y
## [1] 0 0 1
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.6457, p-value = 0.0001673
## $x
## [1] 0.0000000 -1.7688250 0.1940281 -1.1983797 0.3957253 -0.8694238
## [7] 0.6151411 -0.6151411 0.8694238 -0.3957253 1.1983797 -0.1940281
## [13] 1.7688250
##
## $y
## [1] 1 0 1 0 1 0 1 0 1 0 1 0 1
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.63748, p-value = 7.387e-06
## $x
## [1] 0.18911843 0.31863936 0.45376219 -1.95996398 0.59776013 -1.43953147
## [7] 0.75541503 -1.15034938 -0.93458929 0.93458929 -0.75541503 1.15034938
## [13] -0.59776013 -0.45376219 -0.31863936 -0.18911843 1.43953147 -0.06270678
## [19] 0.06270678 1.95996398
##
## $y
## [1] 1 1 1 0 1 0 1 0 0 1 0 1 0 0 0 0 1 0 0 1
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.62622, p-value = 4.249e-07
## $x
## [1] -2.08535557 0.28221615 -1.59321882 0.38032564 0.48224821 -1.32495769
## [7] -1.12814365 -0.96742157 0.58945580 0.70392179 -0.82846465 -0.70392179
## [13] -0.58945580 -0.48224821 0.82846465 -0.38032564 -0.28221615 -0.18675612
## [19] 0.96742157 1.12814365 1.32495769 -0.09297185 0.00000000 0.09297185
## [25] 1.59321882 2.08535557 0.18675612
##
## $y
## [1] 0 1 0 1 1 0 0 0 1 1 0 0 0 0 1 0 0 0 1 1 1 0 0 0 1 1 0
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.61691, p-value = 1.651e-07
## $x
## [1] 0.35293399 0.44658820 -2.11438077 -1.62836141 -1.36448875 0.54434091
## [7] -1.17154617 0.64760358 -1.01449875 -0.87916772 0.75829256 -0.75829256
## [13] 0.87916772 -0.64760358 1.01449875 -0.54434091 -0.44658820 -0.35293399
## [19] -0.26228277 1.17154617 -0.17374106 1.36448875 -0.08654337 0.00000000
## [25] 1.62836141 0.08654337 0.17374106 0.26228277 2.11438077
##
## $y
## [1] 1 1 0 0 0 1 0 1 0 0 1 0 1 0 1 0 0 0 0 1 0 1 0 0 1 0 0 0 1
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.63682, p-value = 2.616e-09
## $x
## [1] -2.28654795 0.05573169 0.11163715 -1.83391464 -1.59321882 -1.42017907
## [7] 0.16789400 -1.28155157 0.22468772 -1.16394958 -1.06056224 0.28221615
## [13] 0.34069483 -0.96742157 -0.88199821 0.40036338 0.46149369 -0.80257189
## [19] 0.52440051 -0.72791329 -0.65710857 0.58945580 0.65710857 -0.58945580
## [25] 0.72791329 0.80257189 0.88199821 0.96742157 -0.52440051 -0.46149369
## [31] 1.06056224 1.16394958 1.28155157 1.42017907 -0.40036338 -0.34069483
## [37] -0.28221615 1.59321882 -0.22468772 1.83391464 -0.16789400 -0.11163715
## [43] 2.28654795 -0.05573169 0.00000000
##
## $y
## [1] 0 1 1 0 0 0 1 0 1 0 0 1 1 0 0 1 1 0 1 0 0 1 1 0 1 1 1 1 0 0 1 1 1 1 0 0 0 1
## [39] 0 1 0 0 1 0 0
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.63647, p-value = 5.953e-08
## $x
## [1] -2.17792307 -0.03687053 0.03687053 -1.70478090 -1.44999891 0.11081291
## [7] 0.18536702 -1.26496924 0.26096740 0.33809166 0.41728414 0.49918784
## [13] 0.58458986 -1.11533736 0.67448975 0.77020808 0.87356913 -0.98723099
## [19] 0.98723099 -0.87356913 -0.77020808 1.11533736 -0.67448975 1.26496924
## [25] 1.44999891 -0.58458986 1.70478090 -0.49918784 -0.41728414 2.17792307
## [31] -0.33809166 -0.26096740 -0.18536702 -0.11081291
##
## $y
## [1] 0 1 1 0 0 1 1 0 1 1 1 1 1 0 1 1 1 0 1 0 0 1 0 1 1 0 1 0 0 1 0 0 0 0
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.63638, p-value = 3.414e-13
## $x
## [1] -2.51912447 0.00000000 0.02949402 0.05901372 -2.10496772 -1.88950996
## [7] 0.08858496 -1.73719396 0.11823387 0.14798711 0.17787196 0.20791656
## [13] 0.23815006 0.26860288 -1.61698114 -1.51642562 0.29930691 -1.42921942
## [19] -1.35170224 -1.28155157 0.33029579 0.36160523 -1.21719663 -1.15752281
## [25] 0.39327329 -1.10170893 0.42534083 0.45785193 0.49085443 -1.04913140
## [31] -0.99930480 -0.95184328 0.52440051 0.55854749 0.59335862 0.62890422
## [37] 0.66526290 -0.90643457 -0.86282209 -0.82079209 -0.78016439 -0.74078547
## [43] 0.70252320 0.74078547 -0.70252320 0.78016439 -0.66526290 -0.62890422
## [49] -0.59335862 -0.55854749 0.82079209 0.86282209 0.90643457 -0.52440051
## [55] -0.49085443 -0.45785193 -0.42534083 0.95184328 0.99930480 1.04913140
## [61] 1.10170893 -0.39327329 -0.36160523 -0.33029579 -0.29930691 1.15752281
## [67] -0.26860288 1.21719663 1.28155157 1.35170224 -0.23815006 1.42921942
## [73] -0.20791656 1.51642562 -0.17787196 -0.14798711 -0.11823387 1.61698114
## [79] 1.73719396 1.88950996 2.10496772 -0.08858496 -0.05901372 2.51912447
## [85] -0.02949402
##
## $y
## [1] 0 1 1 1 0 0 1 0 1 1 1 1 1 1 0 0 1 0 0 0 1 1 0 0 1 0 1 1 1 0 0 0 1 1 1 1 1 0
## [39] 0 0 0 0 1 1 0 1 0 0 0 0 1 1 1 0 0 0 0 1 1 1 1 0 0 0 0 1 0 1 1 1 0 1 0 1 0 0
## [77] 0 1 1 1 1 0 0 1 0
##
##
## Shapiro-Wilk normality test
##
## data: muestra
## W = 0.62378, p-value < 2.2e-16
## $x
## [1] -0.228075394 -0.214833790 -2.797207657 -2.421327991 -0.201629750
## [6] -0.188460772 -0.175324397 -2.229612250 -2.095996807 -1.991811422
## [11] -1.905608404 -0.162218209 -1.831604880 -0.149139828 -1.766454744
## [16] -1.708040295 -0.136086908 -0.123057134 -0.110048219 -0.097057902
## [21] -1.654932729 -0.084083943 -0.071124123 -1.606120563 -0.058176239
## [26] -0.045238100 -1.560860364 -0.032307530 -0.019382360 -1.518589009
## [31] -0.006460427 0.006460427 -1.478869270 0.019382360 -1.441354576
## [36] 0.032307530 0.045238100 0.058176239 -1.405765350 -1.371872628
## [41] -1.339486424 0.071124123 0.084083943 0.097057902 0.110048219
## [46] -1.308447267 0.123057134 0.136086908 0.149139828 0.162218209
## [51] -1.278619931 0.175324397 -1.249888705 -1.222153765 0.188460772
## [56] 0.201629750 0.214833790 -1.195328362 -1.169336603 0.228075394
## [61] 0.241357111 0.254681544 0.268051349 -1.144111692 0.281469243
## [66] -1.119594506 0.294938007 0.308460491 -1.095732453 0.322039619
## [71] 0.335678392 0.349379897 0.363147310 0.376983904 0.390893054
## [76] -1.072478525 0.404878245 0.418943080 -1.049790524 -1.027630415
## [81] -1.005963788 -0.984759399 0.433091285 0.447326723 0.461653396
## [86] -0.963988789 -0.943625954 0.476075461 0.490597238 0.505223222
## [91] 0.519958094 -0.923647061 -0.904030212 -0.884755232 0.534806735
## [96] 0.549774243 0.564865945 0.580087415 -0.865803484 0.595444493
## [101] -0.847157720 0.610943307 -0.828801932 0.626590290 -0.810721240
## [106] 0.642392210 0.658356192 0.674489750 0.690800818 -0.792901778
## [111] 0.707297783 -0.775330604 0.723989527 0.740885469 -0.757995614
## [116] 0.757995614 0.775330604 0.792901778 0.810721240 0.828801932
## [121] -0.740885469 -0.723989527 -0.707297783 -0.690800818 -0.674489750
## [126] 0.847157720 -0.658356192 0.865803484 -0.642392210 0.884755232
## [131] -0.626590290 0.904030212 0.923647061 0.943625954 0.963988789
## [136] -0.610943307 -0.595444493 -0.580087415 0.984759399 -0.564865945
## [141] 1.005963788 1.027630415 -0.549774243 1.049790524 -0.534806735
## [146] -0.519958094 -0.505223222 -0.490597238 -0.476075461 1.072478525
## [151] -0.461653396 1.095732453 1.119594506 1.144111692 1.169336603
## [156] 1.195328362 -0.447326723 1.222153765 1.249888705 1.278619931
## [161] 1.308447267 1.339486424 -0.433091285 -0.418943080 -0.404878245
## [166] -0.390893054 -0.376983904 1.371872628 1.405765350 -0.363147310
## [171] 1.441354576 1.478869270 -0.349379897 1.518589009 1.560860364
## [176] -0.335678392 -0.322039619 1.606120563 1.654932729 1.708040295
## [181] 1.766454744 -0.308460491 -0.294938007 -0.281469243 1.831604880
## [186] 1.905608404 -0.268051349 1.991811422 2.095996807 2.229612250
## [191] -0.254681544 2.421327991 2.797207657 -0.241357111
##
## $y
## [1] 1 1 0 0 1 1 1 0 0 0 0 1 0 1 0 0 1 1 1 1 0 1 1 0 1 1 0 1 1 0 1 1 0 1 0 1 1
## [38] 1 0 0 0 1 1 1 1 0 1 1 1 1 0 1 0 0 1 1 1 0 0 1 1 1 1 0 1 0 1 1 0 1 1 1 1 1
## [75] 1 0 1 1 0 0 0 0 1 1 1 0 0 1 1 1 1 0 0 0 1 1 1 1 0 1 0 1 0 1 0 1 1 1 1 0 1
## [112] 0 1 1 0 1 1 1 1 1 0 0 0 0 0 1 0 1 0 1 0 1 1 1 1 0 0 0 1 0 1 1 0 1 0 0 0 0
## [149] 0 1 0 1 1 1 1 1 0 1 1 1 1 1 0 0 0 0 0 1 1 0 1 1 0 1 1 0 0 1 1 1 1 0 0 0 1
## [186] 1 0 1 1 1 0 1 1 0qqnorm(muestra_1)
- 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.
Las gráficas y la prueba de Shapiro evidencian normalidad en la muestra. Sin embargo, debido a problemas de escala, no se logra evidenciar claramente la tendencia por la que se mueven los puntos. No obstante, como se mencionó, las pruebas de normalidad muestran que la muestra sigue una distribución normal.
# 10%
set.seed(1933)
LOTE_5 = rbinom(5 , 1 , .1)
muestra_1 = muestra_aleatoria_aleatoria(LOTE_5) ## Se tomo el 17% de la muestra
prop.table(table(muestra_1))## muestra_1
## 0
## 1#shapiro.test(muestra_1)
qqnorm(muestra_1)
LOTE_10 = rbinom(10 , 1 , .1)
muestra_2 = muestra_aleatoria_aleatoria(LOTE_10) ## Se tomo el 17% de la muestra
prop.table(table(muestra_2))## muestra_2
## 0 1
## 0.9 0.1#shapiro.test(muestra_2)
qqnorm(muestra_2)
LOTE_15 = rbinom(15 , 1 , .1)
muestra_3 = muestra_aleatoria_aleatoria(LOTE_15) ## Se tomo el 17% de la muestra
prop.table(table(muestra_3))## muestra_3
## 0 1
## 0.875 0.125shapiro.test(muestra_3)##
## Shapiro-Wilk normality test
##
## data: muestra_3
## W = 0.4184, p-value = 1.047e-06qqnorm(muestra_3)
LOTE_20 = rbinom(20 , 1 , .1)
muestra_4 = muestra_aleatoria_aleatoria(LOTE_20) ## Se tomo el 17% de la muestra
prop.table(table(muestra_4))## muestra_4
## 0
## 1#shapiro.test(muestra_4)
qqnorm(muestra_4)
LOTE_30 = rbinom(30 , 1 , .1)
muestra_5 = muestra_aleatoria_aleatoria(LOTE_30) ## Se tomo el 17% de la muestra
prop.table(table(muestra_5))## muestra_5
## 0 1
## 0.5 0.5shapiro.test(muestra_5)##
## Shapiro-Wilk normality test
##
## data: muestra_5
## W = 0.72863, p-value = 0.02386qqnorm(muestra_5)
LOTE_50 = rbinom(50 , 1 , .1)
muestra_6 = muestra_aleatoria_aleatoria(LOTE_50) ## Se tomo el 17% de la muestra
prop.table(table(muestra_6))## muestra_6
## 0 1
## 0.9047619 0.0952381shapiro.test(muestra_6)##
## Shapiro-Wilk normality test
##
## data: muestra_6
## W = 0.34141, p-value = 9.399e-09qqnorm(muestra_6)
LOTE_60 = rbinom(60 , 1 , .1)
muestra_7 = muestra_aleatoria_aleatoria(LOTE_60) ## Se tomo el 17% de la muestra
prop.table(table(muestra_7))## muestra_7
## 0 1
## 0.91666667 0.08333333shapiro.test(muestra_7)##
## Shapiro-Wilk normality test
##
## data: muestra_7
## W = 0.32693, p-value = 1.207e-06qqnorm(muestra_7)
LOTE_100 = rbinom(100 , 1 , .1)
muestra_8 = muestra_aleatoria_aleatoria(LOTE_100) ## Se tomo el 17% de la muestra
prop.table(table(muestra_8))## muestra_8
## 0
## 1#shapiro.test(muestra_8)
qqnorm(muestra_8)
LOTE_200 = rbinom(200 , 1 , .1)
muestra_9 = muestra_aleatoria_aleatoria(LOTE_200) ## Se tomo el 17% de la muestra
prop.table(table(muestra_9))## muestra_9
## 0 1
## 0.8571429 0.1428571shapiro.test(muestra_9)##
## Shapiro-Wilk normality test
##
## data: muestra_9
## W = 0.41588, p-value < 2.2e-16qqnorm(muestra_9)
LOTE_500 = rbinom(500 , 1 , .1)
muestra_10 = muestra_aleatoria_aleatoria(LOTE_500) ## Se tomo el 17% de la muestra
prop.table(table(muestra_10))## muestra_10
## 0 1
## 0.8985507 0.1014493shapiro.test(muestra_10)##
## Shapiro-Wilk normality test
##
## data: muestra_10
## W = 0.34486, p-value < 2.2e-16qqnorm(muestra_10)
## Con 90%
set.seed(1944)
LOTE_5 = rbinom(5 , 1 , .1)
muestra_1 = muestra_aleatoria_aleatoria(LOTE_5) ## Se tomo el 17% de la muestra
prop.table(table(muestra_1))## muestra_1
## 0
## 1#shapiro.test(muestra_1)
qqnorm(muestra_1)
LOTE_10 = rbinom(10 , 1 , .1)
muestra_2 = muestra_aleatoria_aleatoria(LOTE_10) ## Se tomo el 17% de la muestra
prop.table(table(muestra_2))## muestra_2
## 0
## 1#shapiro.test(muestra_2)
qqnorm(muestra_2)
LOTE_15 = rbinom(15 , 1 , .1)
muestra_3 = muestra_aleatoria_aleatoria(LOTE_15) ## Se tomo el 17% de la muestra
prop.table(table(muestra_3))## muestra_3
## 0
## 1#shapiro.test(muestra_3)
qqnorm(muestra_3)
LOTE_20 = rbinom(20 , 1 , .1)
muestra_4 = muestra_aleatoria_aleatoria(LOTE_20) ## Se tomo el 17% de la muestra
prop.table(table(muestra_4))## muestra_4
## 0 1
## 0.8823529 0.1176471shapiro.test(muestra_4)##
## Shapiro-Wilk normality test
##
## data: muestra_4
## W = 0.38525, p-value = 1.586e-07qqnorm(muestra_4)
LOTE_30 = rbinom(30 , 1 , .1)
muestra_5 = muestra_aleatoria_aleatoria(LOTE_30) ## Se tomo el 17% de la muestra
prop.table(table(muestra_5))## muestra_5
## 0 1
## 0.95454545 0.04545455shapiro.test(muestra_5)##
## Shapiro-Wilk normality test
##
## data: muestra_5
## W = 0.22147, p-value = 7.417e-10qqnorm(muestra_5)
LOTE_50 = rbinom(50 , 1 , .1)
muestra_6 = muestra_aleatoria_aleatoria(LOTE_50) ## Se tomo el 17% de la muestra
prop.table(table(muestra_6))## muestra_6
## 0
## 1#shapiro.test(muestra_6)
qqnorm(muestra_6)
LOTE_60 = rbinom(60 , 1 , .1)
muestra_7 = muestra_aleatoria_aleatoria(LOTE_60) ## Se tomo el 17% de la muestra
prop.table(table(muestra_7))## muestra_7
## 0 1
## 0.8571429 0.1428571shapiro.test(muestra_7)##
## Shapiro-Wilk normality test
##
## data: muestra_7
## W = 0.41806, p-value = 1.207e-10qqnorm(muestra_7)
LOTE_100 = rbinom(100 , 1 , .1)
muestra_8 = muestra_aleatoria_aleatoria(LOTE_100) ## Se tomo el 17% de la muestra
prop.table(table(muestra_8))## muestra_8
## 0 1
## 0.8214286 0.1785714#shapiro.test(muestra_8)
qqnorm(muestra_8)
LOTE_200 = rbinom(200 , 1 , .1)
muestra_9 = muestra_aleatoria_aleatoria(LOTE_200) ## Se tomo el 17% de la muestra
prop.table(table(muestra_9))## muestra_9
## 0 1
## 0.90140845 0.09859155shapiro.test(muestra_9)##
## Shapiro-Wilk normality test
##
## data: muestra_9
## W = 0.34016, p-value = 3.067e-16qqnorm(muestra_9)
LOTE_500 = rbinom(500 , 1 , .1)
muestra_10 = muestra_aleatoria_aleatoria(LOTE_500) ## Se tomo el 17% de la muestra
prop.table(table(muestra_10))## muestra_10
## 0 1
## 0.8902954 0.1097046shapiro.test(muestra_10)##
## Shapiro-Wilk normality test
##
## data: muestra_10
## W = 0.36041, p-value < 2.2e-16qqnorm(muestra_10)
Conclusión:
Se ha evidenciado que dada la elección de una muestra aleatoria, el aumento del tamaño muestral converge hacia una distribución de tipo normal. Las pruebas estadísticas proporcionan evidencia de que el teorema central del límite mantiene su consistencia; sin embargo, este resultado no se cumple en las muestras 1, 2, 6 y 8 debido a una elección de muestra muy pequeña.
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 \(X^{∗}_{1}\). Después de anotado el valor se regresa X∗1 a la caja y se extrae el valor \(X^{∗}_{2}\), regresandolo nuevamente. Este procedimiento se repite hasta completar una muestra de tamaño \(n\), \(X^{∗}_{1}\), \(X^{∗}_{2}\),\(X^{∗}_{3}\),\(X^{∗}_{n}\), 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 {x^{*}_{i}} , 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étodos | Formula |
|---|---|
| Método 1 | \((P_{2.5};P_{97.5})\) |
| Método 2 | \((2\bar{x}-P_{97.5};2\bar{x}-P_{2.5})\) |
Construya el intervalo de confianza por los dos métodos y compare los resultados obtenidos. Comente los resultados. Confiaría en estas estimaciones?
# Estimar los intervalos de confianza para la media de una muestra de datos
x = c(7.69, 4.97, 6.49, 4.34, 6.24, 4.45) # Datos de la muestra
n = 1000
muest_bootstrap = matrix(NA, nrow = n, ncol = length(x))
set.seed(1015)
for (i in 1:n) {
muestra_bootstrap = sample(x, replace = TRUE)
muest_bootstrap[i, ] = muestra_bootstrap
}
medias_x_muestra = apply(muest_bootstrap, 1, mean)
metodo1 = quantile(medias_x_muestra, probs = c(0.025, 0.975)) #Se calcula el IC para el método 1
metodo1## 2.5% 97.5%
## 4.778167 6.678625metodo2 = c(2 * mean(medias_x_muestra) - metodo1[2], 2 * mean(medias_x_muestra) - metodo1[1]) #Se calcula el IC para el método 2
metodo2## 97.5% 2.5%
## 4.715962 6.616420# Calcular intervalo de confianza para la media.
media_metodos = c(mean(medias_x_muestra) - sd(medias_x_muestra), mean(medias_x_muestra) + sd(medias_x_muestra))
# Graficar el histograma con líneas para intervalos de confianza
hist(medias_x_muestra,
main = "Histograma con intervalos de confianza y medida",
xlab = "Media Bootstrap",
ylab = "Frecuencia", col = "grey")
abline(v = metodo1, col = "red", lwd = 2, lty = 2)
abline(v = metodo2, col = "blue", lwd = 2, lty = 2)
abline(v = mean(medias_x_muestra), col = "green", lwd = 2, lty = 2)
legend("topright", legend = c("Método 1", "Método 2", "Media"), col = c("red", "blue", "green"), lty = 2, lwd = 2)
Conclusion:
Las bandas de selección exhiben una significancia consistente y robusta frente a las muestras proporcionadas. Este hallazgo sugiere que el proceso de selección se realiza de manera efectiva y precisa en el conjunto de datos evaluado. El histograma evidencia las medias obtenidas con base en los métodos del boostrap. Las líneas paralelas roja y azul muestran los límites de los intervalos de confianza que se obtuvieron a partir de los dos métodos, las cuales dan resultados positivos; con una confianza del 95% el IC para el método 1 resulta entre 4.77 y 6.67 mientras que para el segundo método resulta entre 4.71 y 6.61 en relación con la media de aprox. 5.5.