Problema 1

Estimación del valor de \(\pi\)

La siguiente figura sugiere como estimar el valor de π con una simulación. En la figura, un circuito con un área igual a π/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 π/4. Por tanto, se puede estimar el valor de π/4 al contar el número de puntos dentro del círculo, para obtener la estimación de π/4 . De este último resultado se encontrar una aproximación para el valor de π.

Pasos sugeridos:

A.Genere n coordenadas x: X1, . . . , Xn. 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).

n <- 1000
x = runif(n,0,1)

B.Genere 1000 coordenadas y : Y1,…,Yn, utilizando nuevamente la distribución uniforme con valor mínimo de 0 y valor máximo de 1.

# Definir el número de puntos a generar

n <- 1000

x = runif(n,0,1)

C.Cada punto (Xi,Yi) se encuentra dentro del círculo si su distancia desde el centro (0.5,0.5) es menor a 0.5. Para cada par (Xi,Yi) determine si la distancia desde el centro es menor a 0.5. Esto último se puede realizar al calcular el valor (Xi−0.5)2+(Yi−0.5)2, que es el cuadrado de la distancia, y al determinar si es menor que 0.25.

n <- 1000

# Generar coordenadas aleatorias para x e y
x <- runif(n, 0, 1)
y <- runif(n, 0, 1)

# Calcular la distancia de cada punto al centro (0.5, 0.5)
distancia <- (x - 0.5)^2 + (y - 0.5)^2

# Contar cuántos puntos están dentro del círculo (distancia < radio^2)
contador <- as.numeric(distancia < 0.25)

D.¿Cuántos de los puntos están dentro del círculo? ¿Cuál es su estimación de π?

# Contar cuántos puntos están dentro del círculo (distancia < radio^2)
contador <- as.numeric(distancia < 0.25)

cat("Suma de puntos dentro del circulo : ", sum(contador))
Suma de puntos dentro del circulo :  795
# Calcular el área aproximada del círculo usando la proporción de puntos dentro del círculo respecto al total de puntos generados
area_aprox <- sum(contador) / n * 4

#¿Cuál es su estimación de π?



cat("Estimación de pi :", sum(contador)/n*4)
Estimación de pi : 3.18

Conclusión:

Se concluye que entre mas grande el tamaño de la muestra sea , el valor se aproxima a el pi.

n= 100
x = runif(n,0,1)
y = runif(n,0,1)
distancia = (x-0.5)^2+(y-0.5)^2
contador= as.numeric(distancia<0.25)
cat("Numero de puntos ",sum(contador))
Numero de puntos  73
cat("Estimación de pi ",sum(contador)/n*4)
Estimación de pi  2.92
# Definir el número de puntos a generar
n <- 100

# Generar coordenadas aleatorias para x e y
x <- runif(n, 0, 1)
y <- runif(n, 0, 1)

# Calcular la distancia de cada punto al centro (0.5, 0.5)
distancia <- (x - 0.5)^2 + (y - 0.5)^2

# Contar cuántos puntos están dentro del círculo (distancia < radio^2)
contador <- as.numeric(distancia < 0.25)

# Calcular el área aproximada del círculo usando la proporción de puntos dentro del círculo
area_aprox <- sum(contador) / n * 4

# Crear un gráfico de dispersión
plot(x, y, col = ifelse(distancia < 0.25, "#9fc3e7", "#E0E5E9"), pch = 20,
     main = "Simulación de puntos aleatorios en un círculo",
     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 = "#0051a1", lwd = 2)

n= 1000
x = runif(n,0,1)
y = runif(n,0,1)
distancia = (x-0.5)^2+(y-0.5)^2
contador= as.numeric(distancia<0.25)
cat("Numero de puntos ",sum(contador))
Numero de puntos  790
cat("Estimación de pi ",sum(contador)/n*4)
Estimación de pi  3.16
# Definir el número de puntos a generar
n <- 1000

# Generar coordenadas aleatorias para x e y
x <- runif(n, 0, 1)
y <- runif(n, 0, 1)

# Calcular la distancia de cada punto al centro (0.5, 0.5)
distancia <- (x - 0.5)^2 + (y - 0.5)^2

# Contar cuántos puntos están dentro del círculo (distancia < radio^2)
contador <- as.numeric(distancia < 0.25)

# Calcular el área aproximada del círculo usando la proporción de puntos dentro del círculo
area_aprox <- sum(contador) / n * 4

# Crear un gráfico de dispersión
plot(x, y, col = ifelse(distancia < 0.25, "#9fc3e7", "#E0E5E9"), pch = 20,
     main = "Simulación de puntos aleatorios en un círculo",
     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 = "#0051a1", lwd = 2)

n= 10000
x = runif(n,0,1)
y = runif(n,0,1)
distancia = (x-0.5)^2+(y-0.5)^2
contador= as.numeric(distancia<0.25)
cat("Numero de puntos ",sum(contador))
Numero de puntos  7844
cat("Estimación de pi ",sum(contador)/n*4)
Estimación de pi  3.1376
# Definir el número de puntos a generar
n <- 10000

# Generar coordenadas aleatorias para x e y
x <- runif(n, 0, 1)
y <- runif(n, 0, 1)

# Calcular la distancia de cada punto al centro (0.5, 0.5)
distancia <- (x - 0.5)^2 + (y - 0.5)^2

# Contar cuántos puntos están dentro del círculo (distancia < radio^2)
contador <- as.numeric(distancia < 0.25)

# Calcular el área aproximada del círculo usando la proporción de puntos dentro del círculo
area_aprox <- sum(contador) / n * 4

# Crear un gráfico de dispersión
plot(x, y, col = ifelse(distancia < 0.25, "#9fc3e7", "#E0E5E9"), pch = 20,
     main = "Simulación de puntos aleatorios en un círculo",
     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 = "#0051a1", lwd = 2)

n= 100000
x = runif(n,0,1)
y = runif(n,0,1)
distancia = (x-0.5)^2+(y-0.5)^2
contador= as.numeric(distancia<0.25)
cat("Numero de puntos ",sum(contador))
Numero de puntos  78549
cat("Estimación de pi ",sum(contador)/n*4)
Estimación de pi  3.14196
# Definir el número de puntos a generar
n <- 100000

# Generar coordenadas aleatorias para x e y
x <- runif(n, 0, 1)
y <- runif(n, 0, 1)

# Calcular la distancia de cada punto al centro (0.5, 0.5)
distancia <- (x - 0.5)^2 + (y - 0.5)^2

# Contar cuántos puntos están dentro del círculo (distancia < radio^2)
contador <- as.numeric(distancia < 0.25)

# Calcular el área aproximada del círculo usando la proporción de puntos dentro del círculo
area_aprox <- sum(contador) / n * 4

# Crear un gráfico de dispersión
plot(x, y, col = ifelse(distancia < 0.25, "#9fc3e7", "#E0E5E9"), pch = 20,
     main = "Simulación de puntos aleatorios en un círculo",
     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 = "#0051a1", lwd = 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%.

poblacion = c(rep(1,500),rep(0,500))
barplot(table(poblacion),las=1, col= c("#e0e5e9","#9fc3e7"))

B.Genere una función que permita:

  • Obtener una muestra aleatoria de la población y
#Genere una función que permita obtener una muestra aleatoria de la población y 


muestra = function(x,n){sample(x,n,replace = TRUE)}
muestra = function(x,n){sample(x,n,replace = TRUE)}


  • Calcule el estimador de la proporción muestral \(\widehat{p}\) para un tamaño de muestra dado n.

    • tamaño de las muestras : n = 5
    • número de muestras : m = 10
n=5
m=10

y= matrix(muestra(poblacion,n*m), ncol=n)

#y
  • Calcule el estimador de la proporción muestral para un tamaño de muestra dado n.
#calcule el estimador de la proporción muestral para un tamaño de muestra dado n.

p=apply(y,1,sum)/n

hist(p, las=1, main = "proporción muestral", col = "#9fc3e7")

C. Repita el escenario anterior (b) n=500 veces y analice los resultados en cuanto al comportamiento de los 500 resultados del estimador \(\widehat{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.

poblacion = c(rep(1,500),rep(0,500))

n = 500

muestra = function(x,n){sample(x,n,replace = TRUE)}
n=500
m=500

y= matrix(muestra(poblacion,n*m), ncol=m)

#y

p=apply(y,1,sum)/n

hist(p,main = "Proporción, n=500",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")

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

# n=5

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=5
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p5=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p5,main = "Proporción, n=5",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p5,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

En la grafica tanto como en la prueba de shapiro ,como el valor de p es mayor a 0.05 no se rechaza la hipotesis de normalidad, se asume que la variable es normal, y la grafica indica que las estimaciones de la proporcion para un n=5 son normales

shapiro.test(p5)

    Shapiro-Wilk normality test

data:  p5
W = 0.93083, p-value = 1.936e-14



#n=10

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=10
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p6=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p6,main = "Proporción, n=10",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p6,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p6)

    Shapiro-Wilk normality test

data:  p6
W = 0.96142, p-value = 3.601e-10
#n=15

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=15
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p7=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p7,main = "Proporción, n=15",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p7,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p7)

    Shapiro-Wilk normality test

data:  p7
W = 0.97871, p-value = 1.109e-06
#n=20

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=20
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p8=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p8,main = "Proporción, n=20",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p8,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p8)

    Shapiro-Wilk normality test

data:  p8
W = 0.98127, p-value = 4.837e-06
#n=30


poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=30
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p9=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p9,main = "Proporción, n=30",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p9,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p9)

    Shapiro-Wilk normality test

data:  p9
W = 0.98649, p-value = 0.0001368
#n=50

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=50
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p10=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p10,main = "Proporción, n=50",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p10,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p10)

    Shapiro-Wilk normality test

data:  p10
W = 0.98929, p-value = 0.001029
#n=60

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=60
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p11=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p11,main = "Proporción, n=60",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p11,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p11)

    Shapiro-Wilk normality test

data:  p11
W = 0.99262, p-value = 0.01439
#n=100

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=100
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p12=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p12,main = "Proporción, n=100",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p12,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p12)

    Shapiro-Wilk normality test

data:  p12
W = 0.99185, p-value = 0.007637
#n=200

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=200
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p13=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p13,main = "Proporción, n=200",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p13,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p13)

    Shapiro-Wilk normality test

data:  p13
W = 0.99643, p-value = 0.3306
#n=500

poblacion = c(rep(1,500),rep(0,500))
muestra = function(x,n){sample(x,n,replace = TRUE)}
n=500
m=500

y5= matrix(muestra(poblacion,n*m), ncol=n)


p13=apply(y5,1,sum)/n

par(mfrow = c(1, 2))
hist(p13,main = "Proporción, n=500",ylab = "Frecuencia", xlab = "p", col="#9fc3e7")
qqnorm(p13,main = " ", ylab = "Percenties Muestra", xlab = "Percentiles Normal", col="#0051a1")

shapiro.test(p13)

    Shapiro-Wilk normality test

data:  p13
W = 0.99629, p-value = 0.298
#como el valor p es mayor a 0.05 no se rechaza la hipotesis de normalidad, se asume que la variable es normal, que la estimacion es normal.
shapiro.test(p)

    Shapiro-Wilk normality test

data:  p
W = 0.9963, p-value = 0.2986
#tanto la prueba de shapiro como el grafico indica que las estimacones de la proporcion para un n= 500 son normales.
qqnorm(p)