Variáveis antitéticas

A precisão é cara

O erro padrão de uma estimativa Monte Carlo cai como 1 / sqrt(N): dez vezes mais precisão exige cem vezes mais trajetórias. Quando uma simulação é cara, uma alavanca melhor que mais trajetórias é menos ruído por trajetória. São os métodos de redução de variância, e o mais simples são as variáveis antitéticas. Elas usam duas vezes os mesmos números aleatórios: com um sorteio z também se roda a trajetória de -z.

import numpy as np

def terminal_price(s0, rate, sigma, t, z):
    return s0 * np.exp((rate - 0.5 * sigma ** 2) * t + sigma * np.sqrt(t) * z)

z = np.array([-1.0, 0.0, 1.0])
print(np.round(terminal_price(100.0, 0.05, 0.2, 1.0, z), 3))         # the price for z = -1, 0, 1
print(np.round(terminal_price(100.0, 0.05, 0.2, 1.0, -z), 3))        # and for the opposite draws: mirror images

Por que parear ajuda

Se z é grande e positivo, a trajetória termina alta e um call paga muito; com -z termina baixa e o call paga pouco. Seus pagamentos são negativamente correlacionados, então a média do par é bem mais estável que cada um. Para um pagamento monótono como o de um call, a média do par (f(z) + f(-z)) / 2 tem variância (1 + corr) / 2 vezes a de um pagamento, com corr negativa. O preço continua sendo a média de uma quantidade não viesada, então não há viés.

import numpy as np

def terminal_price(s0, rate, sigma, t, z):
    return s0 * np.exp((rate - 0.5 * sigma ** 2) * t + sigma * np.sqrt(t) * z)

z = np.random.default_rng(1).standard_normal(20000)
up = np.maximum(terminal_price(100.0, 0.05, 0.2, 1.0, z) - 100.0, 0.0)       # call payoffs for z
down = np.maximum(terminal_price(100.0, 0.05, 0.2, 1.0, -z) - 100.0, 0.0)    # ... and for -z
print(round(float(np.corrcoef(up, down)[0, 1]), 3))                           # negative: the payoffs offset each other

O estimador e sua barra de erro

Trate cada par como uma observação: pair = (payoff(z) + payoff(-z)) / 2. A estimativa é a média descontada dos pagamentos dos pares, e o erro padrão é o desvio-padrão descontado (ddof=1) dos pares sobre sqrt(número de pares). A contagem é de *pares*, pois os dois membros de um par não são independentes. Para julgar o ganho, compare com uma estimativa simples que use o mesmo número de pagamentos de sorteios independentes.

import numpy as np

def terminal_price(s0, rate, sigma, t, z):
    return s0 * np.exp((rate - 0.5 * sigma ** 2) * t + sigma * np.sqrt(t) * z)

def mc_estimate(payoffs, discount):
    return discount * payoffs.mean(), discount * payoffs.std(ddof=1) / np.sqrt(len(payoffs))

def call_payoff(z):
    return np.maximum(terminal_price(100.0, 0.05, 0.2, 1.0, z) - 100.0, 0.0)

rng = np.random.default_rng(3)
z, w = rng.standard_normal(5000), rng.standard_normal(5000)
discount = np.exp(-0.05)

plain = mc_estimate(call_payoff(np.concatenate([z, w])), discount)           # 10,000 independent payoffs
anti = mc_estimate((call_payoff(z) + call_payoff(-z)) / 2, discount)         # 5,000 pairs: the same 10,000 evaluations
print([round(float(v), 4) for v in plain], [round(float(v), 4) for v in anti])
print(round(float(plain[1] / anti[1]), 2), "times smaller error")

# A call profits from pairing; a straddle (symmetric in z) does not
import numpy as np

def terminal_price(z):
    return 100.0 * np.exp((0.05 - 0.5 * 0.2 ** 2) + 0.2 * z)

rng = np.random.default_rng(4)
z = rng.standard_normal(20000)
for name, payoff in (("call", lambda p: np.maximum(p - 100.0, 0.0)), ("straddle", lambda p: np.abs(p - 100.0))):
    pairs = (payoff(terminal_price(z)) + payoff(terminal_price(-z))) / 2
    plain = payoff(terminal_price(np.concatenate([z, rng.standard_normal(20000)])))
    ratio = (plain.std(ddof=1) / np.sqrt(len(plain))) / (pairs.std(ddof=1) / np.sqrt(len(pairs)))
    print(name, round(float(ratio), 2))                       # the call gains (1.41); the straddle gets worse (0.78)

• Use cada sorteio z junto com -z; faça a média dos dois pagamentos em um pagamento de par, bem menos ruidoso para um pagamento monótono.
• O erro padrão é calculado sobre o número de pares, com ddof=1: discount * std(pairs) / sqrt(pairs).
• Compare com uma estimativa simples do mesmo custo; não ajuda em pagamentos simétricos em z, e combina com variáveis de controle.

A seguir: as variáveis de controle, usar uma quantidade de valor exato conhecido para cancelar o ruído.

Lab: o que é avaliado

• Transformar sorteios em preços com terminal_price25 pts
• Reportar a estimativa e seu erro com mc_estimate25 pts
• Parear cada sorteio com seu oposto com antithetic_call_price25 pts
• Medir o ganho com variance_ratio25 pts