Capítulo 55 Usando números aleatórios para o cálculo de Integrais
Seja \(g(x)\) uma função e suponha que se quer calcular \(\theta\) onde: \[\theta = \int_{0}^{1}g(x)\, dx.\] Um método eficiente para realizar esse cálculo é o método de Monte Carlo, que utiliza números aleatórios para aproximar o valor do integral. Vamos ver como isso funciona em detalhes.
Aproximação do Integral Usando Variáveis Aleatórias
Considere \(U\), variável aleatória que segue uma distribuição uniforme no intervalo (0,1). Sabemos que o valor esperado de \(g(U)\), ou seja, \(E[g(U)]\), pode ser usado para calcular o integral de \(g(x)\) no intervalo (0,1). Isto é, podemos reescrever: \[\theta = \int_{0}^{1}g(x)\,dx = E[g(U)].\] Agora, se \(U_1,\ldots,U_k\) são variáveis aleatórias independentes e uniformemente distribuídas no intervalo (0,1), então as variáveis aleatórias \(g(U_1),g(U_2),\ldots,g(U_k)\) também são independentes e identicamente distribuídas, com média \(\theta\). Portanto, pela Lei Forte dos Grandes Números, segue que, com probabilidade 1, \[\sum_{i=1}^{k}\frac{g(U_i)}{k} \longrightarrow E[g(U)]=\theta, \quad \text{quando} \quad k\rightarrow \infty.\] Em outras palavras, à medida que o número \(k\) de variáveis aleatórias \(U_i\) aumenta, a média das avaliações de \(g\) nestes \(U_i\)’s se aproxima do valor do integral \(\theta\).
Portanto, para aproximar \(\theta\), podemos gerar uma grande quantidade de números aleatórios \(u_1, u_2, \ldots, u_k\) no intervalo (0,1) e calcular a média dos valores de \(g(u_i)\): \[\theta \approx \frac{1}{k}\sum_{i=1}^{k}g(u_i).\] Esta abordagem é o que chamamos de método de Monte Carlo para o cálculo de integrais.
55.1 Cálculo de \(\int_{a}^{b}g(x)\, dx\)
Agora, suponha que queremos calcular a integral de \(g(x)\) em um intervalo diferente, digamos [a, b]. Podemos resolver este problema realizando uma substituição de variáveis.
Definimos uma nova variável \(y\) como \[y=\frac{x-a}{b-a}, \quad \text{de modo que} \quad y = \frac{dx}{(b-a)}.\] Isso transforma o integral original em: \[\theta = \int_{0}^{1}g(a+[b-a]y)(b-a)\,dy = \int_{0}^{1}h(y)\,dy,\] onde \(h(y) = (b-a)g(a+[b-a]y)\). Assim, podemos usar novamente o método de Monte Carlo para aproximar essa integral. Basta gerar números aleatórios uniformemente distribuídos no intervalo (0,1), calcular a função \(h(y)\) nesses pontos, e tomar a média. Ou seja, podemos aproximar: \[ \theta \approx \frac{1}{k} \sum_{i=1}^{k} h(y_i), \] onde \(y_1, y_2, \ldots, y_k\) são números aleatórios uniformemente distribuídos em (0,1), e \(h(y) = (b - a)g(a + (b - a)y)\).
Exemplo prático: Vamos considerar um exemplo simples para ilustrar o uso do método de Monte Carlo. Suponha que queremos calcular a integral de \(g(x) = x^2\) no intervalo [0, 1]: \[\theta = \int_{0}^{1} x^2\, dx.\]
Sabemos que a resposta exata é \(\theta = \frac{1}{3}\), mas vamos usar o método de Monte Carlo para aproximar este valor.
Geramos \(n = 10.000\) números aleatórios \(u_i\) no intervalo [0, 1], calculamos \(g(u_i) = u_i^2\) para cada um desses valores e tomamos a média:
\[\theta \approx \frac{1}{10000} \sum_{i=1}^{10000} u_i^2.\]
Ao rodar este processo no R, esperamos obter uma aproximação próxima de \(\frac{1}{3}\).
# Função a ser integrada
f <- function(x) x^2
# Quantidade de números aleatórios
n <- 10000
# Números aleatórios no intervalo [0,1]
u <- runif(n = n, min = 0, max = 1)
# Valor aproximado do integral
sum(f(u))/n## [1] 0.3293839## 0.3333333 with absolute error < 3.7e-15Vamos calcular agora o \[\int_{1}^{2} x^2\, dx.\]
# Limites de integração
a <- 1
b <- 2
# Função a ser integrada
f <- function(x) x^2
# Quantidade de números aleatórios
n <- 10000
# Números aleatórios no intervalo [0,1]
u <- runif(n = n, min = 0, max = 1)
# Valor aproximado do integral
sum(f(a+(b-a)*u)*(b-a))/n## [1] 2.33584## 2.333333 with absolute error < 2.6e-14## [1] 2.31112355.2 Cálculo de \(\int_{0}^{\infty}g(x)\, dx\)
Se o interesse é calcular \(\int_{0}^{\infty}g(x)\, dx\), pode-se aplicar a substituição \[y=\frac{1}{(x+1)}, \quad \text{de modo que} \quad dy=\frac{-dx}{(x+1)^2}=-y^2\,dx,\]
para obter a identidade:
\[\theta = \int_{0}^{1}h(y)\,dy,\] onde \(h(y)=\frac{g(\frac{1}{y}-1)}{y^2}\).
55.3 Cálculo de integrais multidimensionais
O método também pode ser utilizado para o caso de integrais multidimensionais. Suponha que \(g\) é uma função com argumento n-dimensional e se quer calcular: \[\theta = \int_{0}^{1}\int_{0}^{1}\cdots \int_{0}^{1} g(x_1, \ldots,x_n)\, dx_1 \, dx_2\cdots dx_n.\] Tem-se que \(\theta\) pode ser expresso como \(\theta = E[g(U_1,\ldots,U_n)]\), onde \(U_1,\ldots,U_n\) são variáveis aleatórias independentes uniformes no intervalo (0,1). Portanto, para calcular \(\theta\) é preciso gerar \(k\) cojuntos independentes, cada um consistindo de \(n\) variáveis uniformes (0,1) independentes:
\[\begin{align*} U_{1}^{1},&\ldots,U_{n}^{1} \\ U_{1}^{2},&\ldots,U_{n}^{2} \\ & \;\;\vdots \\ U_{1}^{k},&\ldots,U_{n}^{k}. \end{align*}\]
Assim, as variáveis \(g(U_{1}^{i},\ldots, U_{n}^{i})\), \(i=1,\ldots,k\) são todas independentes e identicamente distribuídas com média \(\mu\). Pode-se estimar \(\theta\) por \(\sum_{i=1}^{k}g(U_{1}^{i},\ldots, U_{n}^{i})/k\).
55.4 Exercícios
Use simulação para aproximar as seguintes integrais. Compare sua estimativa com a resposta exata, se conhecida.
1. \(\int_{0}^{1}exp(e^x)\, dx\).
2. \(\int_{0}^{1}(1-x^2)^{3/2}\, dx\).
3. \(\int_{-2}^{2}e^{x+x^2}\, dx\).
4. \(\int_{0}^{\pi/2} cos(x)\, dx\).
5. \(\int_{0}^{\infty}x(1+x^2)^{-2}\, dx\).
6. \(\int_{-\infty}^{\infty}e^{-x^2}\, dx\).
7. \(\int_{0}^{1}\int_{0}^{1}e^{(x+y)^2}\,dy\,dx\).
8. \(\int_0^1 \int_0^1 x^2 + y^2 \, dx \, dy\).