Pasos sugeridos:

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).

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

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.

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

Estimación con n = 1,000

n=1000
x=runif(n,min = 0, max = 1)
y=runif(n,min = 0, max = 1)

d=function(a){(a[1]-0.5)^2 + (a[2]-0.5)^2}
xy=data.frame(x,y)
dist=apply (xy,1,d)
cont=as.numeric(dist<=0.25)

plot(x, y, pch = 20,cex=0.5, xlim = c(0, 1), ylim = c(0, 1), asp = 1, xlab = "X", ylab = "Y",
     symbols(0.5, 0.5, circles = 0.5, add = TRUE, fg = "red", bg = "transparent"))

cat("Puntos dentro del círculo: ", sum(cont))
Puntos dentro del círculo:  769
cat("Estimación de pi: ", 4*(sum(cont)/n))
Estimación de pi:  3.076
cat("Error absoluto: ", abs((4*(sum(cont)/n))-pi))
Error absoluto:  0.06559265



Estimación con n = 10,000

n=10000
x=runif(n,min = 0, max = 1)
y=runif(n,min = 0, max = 1)

d=function(a){(a[1]-0.5)^2 + (a[2]-0.5)^2}
xy=data.frame(x,y)
dist=apply (xy,1,d)
cont=as.numeric(dist<=0.25)

plot(x, y, pch = 20, cex=0.5, xlim = c(0, 1), ylim = c(0, 1), asp = 1, xlab = "X", ylab = "Y",
     symbols(0.5, 0.5, circles = 0.5, add = TRUE, fg = "red", bg = "transparent"))

cat("Puntos dentro del círculo: ", sum(cont))
Puntos dentro del círculo:  7849
cat("Estimación de pi: ", 4*(sum(cont)/n))
Estimación de pi:  3.1396
cat("Error absoluto: ", abs((4*(sum(cont)/n))-pi))
Error absoluto:  0.001992654



Estimación con n = 100,000

n=100000
x=runif(n,min = 0, max = 1)
y=runif(n,min = 0, max = 1)

d=function(a){(a[1]-0.5)^2 + (a[2]-0.5)^2}
xy=data.frame(x,y)
dist=apply (xy,1,d)
cont=as.numeric(dist<=0.25)

plot(x, y, pch = 20, cex=0.5, xlim = c(0, 1), ylim = c(0, 1), asp = 1, xlab = "X", ylab = "Y",
     symbols(0.5, 0.5, circles = 0.5, add = TRUE, fg = "red", bg = "transparent"))

cat("Puntos dentro del círculo: ", sum(cont))
Puntos dentro del círculo:  78451
cat("Estimación de pi: ", 4*(sum(cont)/n))
Estimación de pi:  3.13804
cat("Error absoluto: ", abs((4*(sum(cont)/n))-pi))
Error absoluto:  0.003552654



Conclusión

Mediante este ejercicio se pudo realizar la estimación de un parámetro (pi) usando técnicas estadísticas y probabilísticas.

Para la estimación del valor de Pi, se logra observar que a medida que se usa un n más elevado, el error tiende ser más bajo y el valor estimado es más cercano al valor real, aunque en ciertas simulaciones esto no se cumple para el n = 100,000, ya que el error en casos aislados puede ser mayor que en una simulación con n = 10,000.