Para la resolución de este problema se siguió la siguiente aproximación:
Se sabe que el radio del círculo es 0.5 (1/2) y, por tanto, su área es π*(1/2)^2 que en su igualdad es π/4, la cual a su vez es la probabilidad de que un punto caiga dentro del círculo al ser generado aleatoriamente.
Teniendo lo anterior en cuenta se define que -> (Número de puntos dentro del círculo)/(total puntos generados) ≈ π/4.
En R:
Se creo una función para hacer la simulación con diferentes tamaños de muestra y un seed (semilla de aleatoriedad) predeterminada, permitiendo así la comprobación externa del experimento.
Esta función contiene además la fórmula para calcular la distancia de los puntos respecto al centro y así definir cuántos se encuentran dentro del círculo.
Lo anterior genera una sumatorias de TRUEs en casos que estén dentro del círculo lo cuál luego es el valor que se usa para calcular el pi estimado.
Para mi ejercicio, realicé cinco simulaciones cuyos resultados y detalles se muestran en la Tabla 1. Se observa que la distancia con el Pi real no es necesariamente menor en cuanto la muestra sea mayor, lo cual está relacionado con la aleatoriedad de las muestras usadas.
Se realizaron dos visualizaciones con las dos simulaciones más pequeñas (muestras n=1000 y n=10000) para constatar la naturaleza de las mismas. No se realizaron visualizaciones para las otras pues se considera que no se iba a ver, al ojo humano, mucha diferencia.
# Función para realizar la simulación de Monte Carlo y obtener los puntos generados
simulacion_pi_con_puntos <- function(n, seed) {
set.seed(seed) # Fijar la semilla
x <- runif(n, min = 0, max = 1)
y <- runif(n, min = 0, max = 1)
# Calcular el cuadrado de la distancia desde el centro (0.5, 0.5)
dist_squared <- (x - 0.5)^2 + (y - 0.5)^2
# Contar cuántos puntos están dentro del círculo
dentro_del_circulo <- sum(dist_squared <= 0.25)
# Estimar π
pi_estimado <- (dentro_del_circulo / n) * 4
# Retornar la estimación y los puntos generados
return(list(pi_estimado = pi_estimado, x = x, y = y, dentro_del_circulo = dist_squared <= 0.25))
}
# Escenarios y seeds
escenarios <- c(1000, 10000, 100000, 500000, 1000000)
seeds <- c(42, 123, 456, 789, 101112)
pi_real <- pi # Valor real de pi
# Desactivar la notación científica para la salida
options(scipen = 999)
# Crear un dataframe para almacenar los resultados
resultados <- data.frame(
Tamaño_Muestra = integer(),
Seed = integer(),
Pi_Estimado = numeric(),
Error_Estimado = numeric()
)
# Bucle para realizar la simulación en cada escenario con diferentes seeds
for (i in 1:length(escenarios)) {
# Realizar la simulación
simulacion <- simulacion_pi_con_puntos(escenarios[i], seeds[i])
pi_est <- simulacion$pi_estimado
# Calcular el error estimado
error_estimado <- abs(pi_est - pi_real)
# Agregar los resultados al dataframe
resultados <- rbind(resultados, data.frame(
Tamaño_Muestra = escenarios[i],
Seed = seeds[i],
Pi_Estimado = pi_est,
Error_Estimado = error_estimado
))
# Graficar las primeras dos simulaciones
if (i <= 2) {
plot(simulacion$x, simulacion$y, col = ifelse(simulacion$dentro_del_circulo, "blue", "red"), pch = 16,
main = paste("Simulación con n =", escenarios[i], "y Seed =", seeds[i]),
xlab = "X", ylab = "Y", asp = 1)
# Dibujar el círculo
symbols(0.5, 0.5, circles = 0.5, add = TRUE, inches = FALSE)
}
}
## Tamaño_Muestra Seed Pi_Estimado Error_Estimado
## 1 1000 42 3.104000 0.037592654
## 2 10000 123 3.157600 0.016007346
## 3 100000 456 3.129520 0.012072654
## 4 500000 789 3.139288 0.002304654
## 5 1000000 101112 3.137536 0.004056654