Problema

Estimación del valor de π 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 π .

Estimación

El problema se dividirá en problemas más pequeños por lo que para realizar la estimación del valor de π, se seguirán los pasos descritos a continuación:

Validación de punto

Inicialmente se construye una función que realiza lo siguiente:

  • Permite validar si el punto en cuestión se encuentra dentro del círculo si su distancia desde el centro (0.5,0.5) es menor a 0.5 usando (Xi−0.5)^2 + (Yi−0.5)^2 que es el cuadrado de la distancia, y al determinar si es menor que 0.25.
  • Una vez se produce la validación, la función devuelve 1 si cumple la condición y 0 si no.
validacion_puntos <- function(x, y){
  validacion <- (x-0.5)^2 + (y-0.5)^2
  if (validacion < 0.25){
    return(1)
  }
  else{
    return(0)
  }
}

Conteo de puntos

Se construye una función que realiza lo siguiente:

  • Cuenta la cantidad de puntos que están dentro del círculo.
  • Inicializa el contador conteo en 0
  • Esta función llama a la función validación de punto para que valide si el punto generado con runif está en el círculo. Si es así, aumenta en 1 el contador.
  • Repite el proceso con la estructura de control for hasta que complete el tamaño de la muestra que recibe como parámetro.
  • Imprime la estimación de π.
conteo_puntos <- function(n){
  conteo <- 0
  puntos_x <- numeric(n)
  puntos_y <- numeric(n)
  resultados <- numeric(n)
  for (i in 1:n) {
    a <- runif(1)
    b <- runif(1)
    puntos_x[i] <- a
    puntos_y[i] <- b
    resultados[i] <- validacion_puntos(a, b)
    if(validacion_puntos(a, b) == 1){
      conteo <- conteo + 1
    }
  }
  # Graficamos los puntos
  plot(puntos_x, puntos_y, col = ifelse(resultados, 'green', 'black'),
       asp = 1, xlab = "X", ylab = "Y", main = paste("Distribución de puntos
para n =", n))
  legend("topright", c("Dentro del círculo", "Fuera del círculo"),
         fill = c("green", "black"))
  cat("La cantidad de puntos para una muestra de ", n, " puntos es ", conteo, "\n")
  cat("La estimación de π es ", 4 * conteo / n, "\n")
  
}

Iteración por muestra

Se construye una función que realiza lo siguiente:

  • Se establece una semilla
  • Se crea un vector con los tamaños de muestra
  • Se usa una estructura de control for para iterar cada tamaño de muestra que pasa como parámetro de la función conteo puntos.
set.seed(1234)
muestras <- c(1000, 10000, 100000)
iteracion <- function(){
  for (j in 1:length(muestras)) {
    conteo_puntos(muestras[j])
  }
}

Llamado de función

Finalmente se llama a la función iteración que a su vez activa la función conteo puntos que a su vez activa la función validación puntos

iteracion()

## La cantidad de puntos para una muestra de  1000  puntos es  789 
## La estimación de π es  3.156

## La cantidad de puntos para una muestra de  10000  puntos es  7866 
## La estimación de π es  3.1464

## La cantidad de puntos para una muestra de  1e+05  puntos es  78571 
## La estimación de π es  3.14284