Última atualização: 02 de julho, 2025 às 18:39._

1 Estimação por simulação do valor de \(\pi\)

1.1 Área em um círculo e \(\pi\)

Suponha que esquecemos a fórmula que calcula a área de um círculo. Nem tudo está perdido! Podemos facilmente utilizar simulação computacional para calcular a área do círculo. Por conveniência, vamos considerar determinar a área do círculo unitário de raio 1 unidade e chamar esta constante desconhecida de \(\pi\).

Sabemos que a área do quadrado de vértices \((-1, 1)\), \((1,1)\), \((1, -1)\) e \((-1, -1)\), é 4 uma vez que o tamanho de cada lado é 2. Este quadrado é circunscrito ao círculo.

Se tomarmos um ponto ao acaso do quadrado, podemos calcular a probabilidade de ele ser interior ao círculo unitário. Um ponto aleatório \(Z\) é determinado por coordenadas \((X, Y)\) que possuem distribuições uniformes independentes no intervalo \([-1, 1]\). A probabilidade do ponto estar no círculo é a razão entre a área do círculo e do quadrado, ou seja: \(P(Z \mbox{ está no círculo}) = \frac{\mbox{área do círculo}}{\mbox{área quadrado}} = \frac{\pi}{4}\). Sabendo disto podemos estimar o valor de \(\pi\) (desconhecido) multiplicando por 4 a probabilidade de \(Z\) estar no círculo.

O que resta para obter a probabilidade é gerar um grande número de pontos aleatórios no quadrado e verificar a proporção que também está dentro do círculo. Multiplicando a proporção por 4 temos uma estimativa de \(\pi\).

O objetivo aqui é obter a estimativa para diferentes números de simulações como, por exemplo, 100, 1000, 10000. Avalie a aproximações, já que na realidade conhecemos o valor verdadeiro de \(\pi\).

Dica: você vai precisar saber se um ponto \((x,y)\) aleatório interior do quadrado está também dentro do círculo. O círculo é definido pela origem \((0, 0)\) e pelo raio (que estamos considerando 1). Portanto pode-se calcular a distância do ponto aleatório até o centro do círculo usando Pitágoras e verificar se a distância é inferior ao valor do raio.

1.2 Agulha de Buffon

O problema da “agulha de Buffon” é um dos problemas mais antigos da teoria de probabilidades. Foi proposto pelo naturalista francês George Louis Leclerc posteriormente “Comte de Buffon” e é considerado marco inicial em probabilidade geométrica.

Considere um piso laminado construído com tábuas de mesma largura \(d\). Seja ainda uma agulha de tamanho \(l\), sendo \(l\) menor que a largura \(d\) das tábuas do assoalho. A pergunta é: se uma agulha é lançada sobre este assoalho, qual a probabilidade de que cruze alguma linha que separa duas das tábuas?

Figura: Ilustração do experimento de Comte de Buffon.

O problema pode ser expresso da seguinte forma.
Seja uma região com retas paralelas separadas pela distância \(d\). Seja um segmento de reta de comprimento \(l\) posicionado ao acaso nesta região. Qual a probabilidade de que o segmento intercepte a região?

A solução deste problema leva a: \[ P = P[\text{Interceptar}] = \frac{2l}{\pi d} = \frac{2r}{\pi},\] em que \(r = \frac{l}{d}\).

Os valores de \(l\) e \(d\) serão fornecidos e \(\pi\) é uma constante conhecida. Mas, em analogia com outros problemas de expressões conhecidas a uma constante, vamos considerar que \(\pi\) é desconhecido e deve ser estimado por um experimento de “Monte Carlo” ou seja, um ensaio aleatório.

Um experimento “real” consistiria em usar um piso com estas características e realizar nele diversos lançamentos de alguma agulha (ou simular). Alternativamente, o experimento poderia ser feito usando o computador. O experimento deve fornecer uma estimativa \(\hat{P}\) da proporção de segmentos que interceptam alguma reta e desta forma, \[ \hat{\pi} = \frac{2l}{\hat{P} d} = \frac{2r}{\hat{P}}.\]

Escreva uma função para realizar tal experimento um certo número de vezes e encontrar uma estimativa de \(\pi\).

Dica 1: A figura a seguir ilustra um “lançamento da agulha”. Considere o ângulo \(\theta\) entre a agulha e a linha que divide as tábuas. Para que ocorra a intersecção, a distância (\(y\)) do ponto médio da agulha até a divisão entre as tábuas deve ser maior que a da distância vertical da ponta da agulha até a linha mais próxima, ou seja, \(y \leq \frac{l}{2} \text{sin}(\theta)\).

Dica 2: Por simetria, pode-se considerar o lançamento apenas na metade superior da tábua e a distância até a linha superior, ou seja, centro da agulha entre \([0, d/2]\). Analogamente, para o ângulo de posição da agulha pode-se considerar apenas valores no intervalo \([0, \pi/2]\).

Ilustração:

Código inicial para lançamento de agulha e verificação se intercepta linha:

d <- 2; l <- 1
(y <- runif(1, min = 0, max = d/2))
## [1] 0.1813371
(th <- runif(1, min = 0, max = pi/2))
## [1] 0.4877831
y < (l/2)*sin(th)
## [1] TRUE

Sendo a linguagem vetorial, para 10 lançamentos fazemos:

d <- 2; l <- 1
y <- runif(10, min = 0, max = d/2)
th <- runif(10, min = 0, max = pi/2)
y < (l/2)*sin(th)
##  [1]  TRUE FALSE FALSE FALSE FALSE  TRUE  TRUE  TRUE FALSE FALSE

Um tema interessante é a variância da estimativa obtida. Tem-se que

\[ \begin{align*} X&: \text{número de intersecções em N ensaios (simulações)}\\ X &\sim \text{Bin}(N, P)\\ \text{Var}(\hat{P}) &= \frac{P(1-P)}{N} \end{align*} \] Como \(\hat{\pi}\) está relacionado com \(1/\hat{P}\), vamos optar por examinar a variância de \(\hat{\phi} = 1/\hat{\pi}\). Tem-se que: \[ \phi = \frac{1}{\pi} = \frac{P}{2r}. \] Desta forma, como \(P = 2r\phi\) \[ \text{Var}(\hat{\phi}) = \text{Var}\left(\frac{\hat{P}}{2r}\right) = \frac{1}{4r^2}\frac{2r\phi(1-2r\phi)}{N} = \frac{2r\phi - 4r^2\phi^2}{4r^2N} = \frac{\phi^2}{N}\left(\frac{1}{2r\phi} - 1\right). \]

Supomos aqui que \(r \leq 1\) (ou seja, \(l \leq d\)) e a variância é mínima para \(r = 1\). Portanto escolher \(r = 1\) (\(l = d\)) é a forma mais eficiente de conduzir a simulação no sentido de obter a menor variância para \(1/\hat{\pi}\).

As simulações podem ser representadas visualmente pelos gráficos a seguir.

A partir disto voce pode escrever uma função para simular o lançamento da agulha, obter a estimativa de \(\pi\) e explorar as questões e possibilidades deste exemplo!

1.3 Agulha de Buffon (2)

Examine, no problema da agulha de Buffon, as propriedades do estimador (algorítmico) de \(\pi\).

1.4 Agulha de Buffon (3)

Para rodar o exemplo voce deve ter escolhido valores de \(d\) e \(l\), ou da razão entre estes valores. Verificar se há alguma escolha mais adequada.

1.5 Comparando os experimentos para estimação de \(\pi\)

Desenhe e rode um experimento computacional para comparar a estimação de \(\pi\) pelos dois experimentos de simulação desta lista, o da da razão entre as areas do quadrado e círculo cirscunscrito e o experimento da agulha de Buffon.

Será que um tipo de experimento é preferível ao outro para estimar \(\pi\)? Examine a questão utilizando resultados de simulação.

Apêndice: Propriedades de estimadores

Seja \(T\) um estimador de um parâmetro \(\theta\).

  1. O estimador é dito não-viciado se \({\rm E}[T] = \theta\).
  2. A variância do estimador é \({\rm E}[(T - \text{E}[T])^2]\) e o erro padrão é a raiz quadrada deste valor.
  3. O erro quadrático médio (MSE de sigla em inglês) é \({\rm MSE}(T) = {\rm E}[(T - \theta)^2] = {\rm Var[T]} + (\text{E}(T) - \theta)^2\).
  4. Um estimador \(T_2\) é dito mais eficiente que \(T_1\) se \(\text{Var}(T_1) > \text{Var}(T_2)\) para todos valores de \(\theta\).
  5. \(T\) é dito consistente se \(\underset{n \to \infty}{\lim} \text{Pr}(|T_n - \theta| > \epsilon) = 0\).