Essa resposta está baseada em "Monte Carlo Methods - Malvin H. Kalos and Paula A. Whitlock" seção 4.4.
O método das variáveis antiéticas explora a queda da variância que ocorre quando variáveis aleatórias são correlacionadas negativamente.
Suponha a função a seguir:
\[ G=\int_{1}^{0} g(x)dx \]
Logo, G também pode ser escrito da seguinte maneira:
\[ G=\int_{1}^{0} \frac{1}{2}[g(x) + g(1 -x)dx \]
G pode ser avaliado via Monte Carlo selecionando x uniformemente entre 0 e 1, o valor de g(x) e g(1-x) é determinado e o estimador de G é formado:
\[ G_{N}=1/N \sum_{i=1}^{N}1/2[g(x_i )+g(1-x_i )] \]
Para funções lineares, esse metódo reduzirá a variância. Por exemplo considere a função exponencial \( G = \int_{0}^{1} e^{x}dx \) , a variância associada a uma avaliação Monte Carlo direta é de 0,243 e utilizando o estimador Gn com variáves antiéticas a variância é reduzida para 0,0035.
Essa técnica também pode ser combinada com outras técnicas para reduzir a variância, um exemplo é utilizar importance sampling, a função G se torna:
\[ G = \int_{0}^{1} \frac{1/2[g(x) - g(1-x)]}{f(x)}dx \]
Um f(x) aproximado para o exemplo é: \( f(x) = \frac{24}{25}[1 + \frac{1}{2}(x - \frac{1}{2})^2] \) . E a variância cai para 0,0000012.
Abaixo segue todo o código utilizado para conferir os exemplos.
import numpy as np
import math
import matplotlib.pyplot as plt
def var_strgh_mc(n):
x = np.random.uniform(0,1,n)
a = []
for i in range(n):
h_x = math.exp(x[i])
a = np.append(a,h_x)
GN = sum(a)/n
var = np.var(a)
return var, GN
def est_monte_carlo(n):
x = np.random.uniform(0,1,n)
a = []
for i in range(n):
h_x = 0.5*(math.exp(x[i]) + math.exp(1-x[i]))
#print(h_x)
a = np.append(a,h_x)
GN = sum(a)/n
var = np.var(a)
return var, GN
def ImpSamp_AnthVar(n):
x = np.random.uniform(0,1,n)
a = []
for i in range(n):
f_x = (24/25)*(1 + 0.5*((x[i] - 0.5)**(2)))
h_x = 0.5*((math.exp(x[i]) + math.exp(1-x[i]))/f_x)
a = np.append(a,h_x)
GN = sum(a)/n
var = np.var(a)
return var, GN
if __name__ == '__main__':
numeroAmostras=1000
n=1000
vetor_var = np.empty([numeroAmostras])
vetor_GN = np.empty([numeroAmostras])
vetor_var1 = np.empty([numeroAmostras])
vetor_GN1 = np.empty([numeroAmostras])
vetor_var2 = np.empty([numeroAmostras])
vetor_GN2 = np.empty([numeroAmostras])
for i in range(0,numeroAmostras):
[var, GN] = var_strgh_mc(n)
[var1, GN1] = est_monte_carlo(n)
[var2 , GN2] = ImpSamp_AnthVar(n)
vetor_var[i] = var
vetor_var1[i] = var1
vetor_var2[i] = var2
vetor_GN[i] = GN
vetor_GN1[i] = GN1
vetor_GN2[i] = GN2
plt.figure(0)
plt.hist(vetor_var)
plt.figure(1)
plt.hist(vetor_var1)
plt.figure(2)
plt.hist(vetor_var2)
print('Variância função exponencial:',var)
print('Variância com Variáveis Antiéticas:',var1)
print('Variância com Var Antiéticas com Importance Sampling:', var2)
var_strgh_mc(n)
est_monte_carlo(n)
ImpSamp_AnthVar(n)