Feynman-Kac, a ponte entre equações diferenciais e probabilidade

Existem duas maneiras de precificar uma opção, resolver uma equação diferencial parcial ou calcular uma média sobre cenários simulados. Elas parecem não ter nada a ver uma com a outra e dão o mesmo número. Feynman-Kac é o teorema que explica por quê.

Leitura de cerca de 15 minutos · exige o lema de Itô

Esta página fecha uma série Os três textos anteriores construíram as peças. O lema de Itô deu a ferramenta de cálculo. Black-Scholes como equação do calor percorreu a rota da equação diferencial. Girsanov e a medida neutra ao risco percorreu a rota da esperança. Aqui as duas rotas se encontram.

Imagine dois profissionais recebendo a mesma tarefa, precificar um derivativo. O primeiro monta uma equação diferencial parcial, discretiza um domínio, roda um esquema numérico e entrega um número. O segundo simula alguns milhões de trajetórias de preço, calcula o payoff em cada uma, tira a média e desconta. Entrega outro número.

Os dois números são iguais, até a precisão dos respectivos métodos. Existe um teorema por trás disso, e ele diz algo mais forte que a igualdade de dois cálculos. Diz que uma classe inteira de equações diferenciais é, por baixo, um valor esperado.

O caminho

O enunciado, e o que cada peça faz
A demonstração, em três passos
Por que tudo se resume a matar o drift
A conta que confere, três rotas e um número
O caso particular que você já conhece
O que Feynman-Kac não faz

1. O enunciado

Comece com um processo qualquer com arrasto e ruído, do tipo que aparece em todo modelo de preço.

$$ dX_s = \mu(X_s,s)\,ds + \sigma(X_s,s)\,dW_s $$

Agora considere a seguinte equação diferencial parcial para uma função $ u(x,t) $, com uma condição imposta no fim do intervalo.

$$ \frac{\partial u}{\partial t} + \mu\,\frac{\partial u}{\partial x} + \frac{1}{2}\sigma^{2}\,\frac{\partial^{2} u}{\partial x^{2}} - r\,u = 0, \qquad u(x,T) = \psi(x) $$

Feynman-Kac afirma que a solução dessa equação tem uma fórmula pronta, e que ela é uma esperança.

$$ \boxed{\;u(x,t) = \mathbb{E}\!\left[\,e^{-r(T-t)}\,\psi(X_T)\;\middle|\;X_t = x\,\right]\;} $$

A correspondência é peça a peça.

O termo $ \mu\,u_x $ da equação é o arrasto do processo. Ele diz para onde o cenário tende a caminhar.
O termo $ \tfrac{1}{2}\sigma^{2}u_{xx} $ é a difusão. Ele espalha, e é exatamente o termo que apareceu na equação do calor.
O termo $ -r\,u $ é o desconto, e reaparece do outro lado como o fator $ e^{-r(T-t)} $.
A condição terminal $ \psi $ é o payoff, e reaparece como a função aplicada em $ X_T $.

Uma convenção que vai importar na seção 5. A esperança é tomada sob a medida em que o processo $ X $ foi escrito, seja ela qual for. O teorema não escolhe medida, ele usa a que você entregou.

2. A demonstração

A prova cabe em três passos. O motor é o lema de Itô, mais duas hipóteses técnicas que aparecem no caminho, e eu aviso as duas onde elas entram. O truque inteiro está em escolher o objeto certo para olhar.

Passo 1, construa o candidato

Suponha que $ u $ resolve a equação, e que é duas vezes derivável em $ x $ e uma vez em $ t $, que é o mínimo para o lema de Itô poder ser aplicado. Essa é a primeira hipótese técnica. Defina, ao longo de uma trajetória do processo, a quantidade

$$ Y_s = e^{-r(s-t)}\,u(X_s,\,s), \qquad s \in [t,\,T] $$

Em palavras, pegue a função solução, avalie ela no ponto onde o processo está naquele instante, e desconte até hoje. É o valor presente da solução vista ao longo do caminho.

Passo 2, aplique o lema de Itô

Primeiro a função sozinha. O lema de Itô aplicado a $ u(X_s,s) $ dá

$$ du = u_s\,ds + u_x\,dX_s + \tfrac{1}{2}u_{xx}\,(dX_s)^{2} $$

e aqui entra a regra que a página do lema de Itô estabeleceu, $ dW\,dW = dt $, de onde $ (dX_s)^{2} = \sigma^{2}ds $. Substituindo também $ dX_s = \mu\,ds + \sigma\,dW_s $,

$$ du = \Big(u_s + \mu\,u_x + \tfrac{1}{2}\sigma^{2}u_{xx}\Big)ds + \sigma\,u_x\,dW_s $$

Agora o desconto. Como $ e^{-r(s-t)} $ é determinístico, ele varia de forma previsível, com $ d\big(e^{-r(s-t)}\big) = -r\,e^{-r(s-t)}ds $, e não gera termo cruzado de segunda ordem. Logo a regra do produto vale na forma usual.

$$ dY_s = -r\,e^{-r(s-t)}u\,ds + e^{-r(s-t)}\,du $$

Juntando os dois,

$$ dY_s = e^{-r(s-t)}\Big[\underbrace{u_s + \mu\,u_x + \tfrac{1}{2}\sigma^{2}u_{xx} - r\,u}_{\text{o lado esquerdo da EDP}}\Big]ds \;+\; e^{-r(s-t)}\,\sigma\,u_x\,dW_s $$

Atenção a um detalhe de notação que derruba muita gente. O $ u_s $ ali é a derivada parcial de $ u $ na casa do tempo, com $ x $ congelado. Ele não é a variação total de $ u $ ao longo do caminho, porque a variação total é justamente o que o lado esquerdo inteiro está medindo.

Olhe agora para o colchete. Ele é o lado esquerdo da equação do enunciado, avaliado no ponto onde o processo está naquele instante. E a hipótese é que $ u $ resolve aquela equação em todo ponto do domínio, então ela vale em $ (X_s,s) $ seja qual for a trajetória sorteada. Aquele colchete vale zero caminho a caminho.

$$ dY_s = e^{-r(s-t)}\,\sigma\,u_x\,dW_s $$

Passo 3, colha o resultado

O que sobrou não tem termo em $ ds $. Sobrou só ruído, e integral de Itô tem média zero, porque ela soma incrementos $ dW $ de média zero multiplicados por pesos que já estavam decididos antes de o incremento acontecer. Cada parcela tem média zero e a soma também. Logo a esperança de $ Y $ não muda com o tempo.

$$ \mathbb{E}\big[\,Y_T \mid X_t = x\,\big] = Y_t $$

A ressalva vem aqui, não depois

Drift zero entrega, a rigor, um martingale local, e martingale local não garante $ \mathbb{E}[Y_T] = Y_t $. Para fechar o passo é preciso uma hipótese de integrabilidade sobre $ \sigma\,u_x $, que controla o quanto o termo de ruído pode crescer. Nos casos desta página ela vale com folga, e em finanças costuma valer, mas seria desonesto usar a conclusão por quatro seções e só avisar no fim. O terceiro item da seção 6 volta ao assunto.

Falta apenas traduzir as duas pontas. No instante inicial o desconto é 1 e o processo está em $ x $, então $ Y_t = u(x,t) $. No instante final a condição terminal entra em ação, $ u(X_T,T) = \psi(X_T) $, então $ Y_T = e^{-r(T-t)}\psi(X_T) $. Substituindo, a fórmula aparece.

O que realmente aconteceu aqui

Repare que em nenhum momento resolvemos a equação. O que fizemos foi construir um objeto cujo drift é exatamente o lado esquerdo da equação, e então usar a hipótese para zerá-lo. A EDP não foi resolvida, foi usada como certificado de que uma certa quantidade não tem tendência.

3. Por que tudo se resume a matar o drift

Essa estrutura já apareceu antes nesta série, e vale comparar as duas aparições com cuidado, porque elas são parecidas o bastante para serem confundidas.

Em Girsanov, o arrasto $ \mu $ do ativo foi trocado por $ r $, mudando a probabilidade. Aqui, quem fica sem arrasto é outro objeto, o preço descontado $ Y $, e a ferramenta não é trocar medida, é impor uma equação. Os dois atacam alvos diferentes e por isso não competem. Eles se encaixam, e a seção 5 mostra exatamente o encaixe, Girsanov entrega o processo já sob $ \mathbb{Q} $ e Feynman-Kac converte esse processo em equação.

O que os une é o destino. Nos dois casos o produto final é um martingale, e martingale é o conceito central de precificação porque preço justo descontado não pode ter tendência previsível sob a medida de precificação. Repare na qualificação, porque ela é onde mora o erro comum. Sob a medida real $ \mathbb{P} $ o preço tem tendência, e essa tendência é o prêmio de risco, remuneração legítima de quem carrega a posição. Não há dinheiro na mesa ali. O que não pode existir é tendência previsível no preço descontado sob $ \mathbb{Q} $, e é esse objeto que a demonstração da seção 2 constrói.

Conferido numericamente, o martingale

Simulei 400 mil trajetórias e acompanhei $ \mathbb{E}[Y_s] $ ao longo do tempo, com $ \mu = 0{,}35 $, $ \sigma = 0{,}60 $, $ r = 0{,}05 $, partindo de $ x = 1 $ com dois anos até o vencimento. Se a demonstração estiver certa, esse número não pode ter tendência nenhuma.

instante $ s $	$ \mathbb{E}[Y_s] $ medido	erro padrão
0,10	3,26659	0,00093
0,50	3,26534	0,00209
1,00	3,26657	0,00301
1,50	3,26405	0,00373
2,00	3,26274	0,00437

O valor inicial $ Y_t $ vale $ 3{,}26646 $. A série oscila em torno dele dentro do erro de medição e não anda para lado nenhum, que é exatamente o que significa não ter drift. Repare também no erro padrão crescendo com o tempo, porque as trajetórias se espalham. A média fica parada, a dispersão não.

4. A conta que confere, três rotas e um número

O teorema afirma que uma equação diferencial e uma esperança devolvem o mesmo valor. Nenhuma conta numérica prova um teorema, a prova está na seção 2. O que a conta faz é procurar erro, e não achar erro em quatro pontos com três métodos é informação real sobre se o enunciado foi escrito corretamente.

Escolhi um caso construído de propósito para ter forma fechada. Vale ser honesto sobre o que essa terceira coluna é. Ela sai da rota probabilística, de $ \mathbb{E}[X_T^{2}] $, então ela não é um juiz independente das duas outras, é a rota de Monte Carlo resolvida exatamente em vez de simulada. A independência real está entre a coluna das diferenças finitas, que nunca sorteia nada, e as outras duas. O processo é um browniano com arrasto, $ dX = \mu\,dt + \sigma\,dW $, o payoff é $ \psi(x) = x^{2} $, e existe desconto a uma taxa $ r $. Nesse caso a solução sai na mão,

$$ u(x,t) = e^{-r\tau}\Big[(x + \mu\tau)^{2} + \sigma^{2}\tau\Big], \qquad \tau = T - t $$

Note que esse exemplo não é Black-Scholes, o que é proposital. Serve para mostrar que o alcance do teorema é uma família de equações, e não apenas opções.

Com $ \mu = 0{,}35 $, $ \sigma = 0{,}60 $, $ r = 0{,}05 $ e dois anos até o vencimento, as três rotas dão isto.

$ x $	forma fechada	Monte Carlo	diferenças finitas
−2,0	2,18066	2,17941	2,18064
0,0	1,09485	1,09435	1,09485
1,5	5,03090	5,03107	5,03089
4,0	20,63934	20,64210	20,63934

A coluna do Monte Carlo usou oito milhões de trajetórias por célula e não olhou para equação nenhuma, só simulou o processo e tirou média descontada. A coluna das diferenças finitas resolveu a EDP para trás no tempo em 36 mil passos, numa malha explícita truncada em $ x = \pm 14 $, e nunca simulou uma trajetória sequer. A fronteira merece nota, porque o payoff $ x^{2} $ cresce sem limite e não dá para fixar valor nas pontas. Usei extrapolação linear da primeira derivada nas duas bordas, equivalente a impor segunda derivada constante ali. Os dois métodos não compartilham nenhum passo de cálculo.

O maior afastamento do Monte Carlo em relação à forma fechada foi de 1,61 erro padrão, que é ruído estatístico comum. O maior erro das diferenças finitas foi de 1,5 × 10⁻⁵ em valor absoluto, que é erro de discretização esperado. O que sustenta o teste não é o tamanho dos números, é a assinatura deles. O erro do Monte Carlo encolhe com a raiz do número de trajetórias e o das diferenças finitas encolhe com o refino da malha, cada um obedecendo à sua própria taxa. Se o enunciado estivesse errado, sobraria uma diferença fixa que nenhum dos dois refinos removeria.

5. O caso particular que você já conhece

Agora troque as peças do enunciado geral por peças específicas. Tome o processo sob a medida neutra ao risco, aquela construída na página de Girsanov, em que a ação cresce à taxa livre de risco.

$$ dS = r S\,dt + \sigma S\,d\tilde{W}, \qquad \mu = rS, \quad \sigma_{\text{proc}} = \sigma S $$

O til em $ \tilde{W} $ é o mesmo da página de Girsanov, e marca que esse browniano é browniano sob $ \mathbb{Q} $, não sob $ \mathbb{P} $. Toda a seção 2 vale sem alterar uma letra, porque o teorema nunca perguntou qual medida você trouxe. Ele só exige coerência, o processo e a esperança precisam morar na mesma.

Substituindo na equação geral do enunciado, com $ \psi $ sendo o payoff da opção,

$$ \frac{\partial V}{\partial t} + rS\,\frac{\partial V}{\partial S} + \frac{1}{2}\sigma^{2}S^{2}\,\frac{\partial^{2} V}{\partial S^{2}} - rV = 0 $$

que é, sem tirar nem pôr, a equação de Black-Scholes. E o lado direito do teorema vira

$$ V(S,t) = e^{-r(T-t)}\,\mathbb{E}_{\mathbb{Q}}\big[\,\text{payoff}(S_T)\,\big] $$

que é, sem tirar nem pôr, a fórmula de precificação neutra ao risco.

A série fecha aqui

A primeira página desta trilha mostrou que a equação de Black-Scholes é a equação do calor disfarçada, e a resolveu por mudança de variáveis. A terceira mostrou como chegar ao preço por mudança de medida, sem tocar em equação diferencial. As duas entregavam o mesmo número e a explicação ficou prometida. É esta. As duas rotas são a mesma identidade matemática lida em duas direções, e Feynman-Kac é o enunciado dessa identidade.

6. O que Feynman-Kac não faz

Não vale para equação não linear

A demonstração depende de o colchete do passo 2 ser exatamente o lado esquerdo da equação. Isso só acontece porque a equação é linear em $ u $ e nas suas derivadas. Modelos com custo de transação, com restrição de liquidez ou com volatilidade dependente da própria carteira geram equações não lineares, e ali a ponte quebra. Existem generalizações, mas elas vêm em camadas. As equações diferenciais estocásticas retrógradas cobrem a não linearidade que aparece na função e na primeira derivada, o caso chamado semilinear. Os três exemplos que acabei de citar são piores do que isso, porque a não linearidade neles ataca o termo de segunda ordem, e aí só a versão de segunda ordem dá conta. Em qualquer dos casos, são outro assunto e bem mais pesadas.

Não diz qual rota é melhor

Ele afirma que as duas dão o mesmo número, e fica em silêncio sobre custo. Na prática a escolha é ditada pela dimensão. Com uma ou duas variáveis de estado, resolver a EDP costuma ser mais rápido e mais preciso, e na tabela da seção 4 as diferenças finitas chegaram cerca de duas ordens de grandeza mais perto da forma fechada. Essa comparação não foi normalizada por tempo de máquina, então leia como precisão obtida em cada configuração, não como eficiência medida. À medida que as variáveis de estado se acumulam, a malha cresce exponencialmente e o Monte Carlo vira a escolha padrão, porque a taxa de convergência dele, a raiz do número de simulações, não depende da dimensão. A taxa, não o custo. A constante que multiplica essa taxa pode piorar muito com a dimensão, e cada trajetória fica mais cara de gerar. O que não acontece é a explosão exponencial da malha, e é essa assimetria que decide. É aqui que o teorema deixa de ser curiosidade e vira decisão de engenharia.

O passo do martingale pede cuidado

Já avisado na seção 2, e vale fechar aqui. A ressalva de integrabilidade sobre $ \sigma\,u_x $ é da mesma família do cuidado que aparece na condição de Novikov, e existe pelo mesmo motivo, porque um objeto que parece bem-comportado pontualmente pode não ser no agregado. Uma condição suficiente e comum é exigir que a esperança do quadrado de $ \sigma\,u_x $, integrada no tempo, seja finita. Quando $ \sigma $ e $ u_x $ têm crescimento no máximo linear, como no exemplo da seção 4 e como em Black-Scholes, isso vale sem esforço.

Veja a difusão que está dos dois lados da ponte

O termo de segunda ordem da equação e o espalhamento das trajetórias simuladas são o mesmo fenômeno visto de dois ângulos. Construí dois simuladores que mostram cada lado, a EDP resolvida ao vivo e as trajetórias se abrindo com a raiz do tempo. Gratuitos e sem cadastro.

Abrir o laboratório da EDP

Para fixar

Feche a aba e responda no papel.

Por que a demonstração precisa descontar, ou seja, por que $ Y_s $ leva o fator $ e^{-r(s-t)} $ em vez de ser apenas $ u(X_s,s) $?
Em que passo exato a hipótese de que $ u $ resolve a equação é usada?
Se o payoff fosse $ \psi(x) = x $ em vez de $ x^{2} $, qual seria a forma fechada?
Você precisa precificar um derivativo que depende de sete ativos. Qual rota escolhe?

Respostas da terceira e da quarta. Com $ \psi(x) = x $ sai $ u(x,t) = e^{-r\tau}(x + \mu\tau) $, porque a esperança de $ X_T $ é $ x + \mu\tau $. O termo $ \sigma^{2}\tau $ some porque payoff linear passa a variância batido, ao contrário do payoff quadrático, em que a variância entra via $ \mathbb{E}[X_T^{2}] = (\mathbb{E}X_T)^{2} + \text{Var}\,X_T $. Na quarta, Monte Carlo, porque uma malha em sete dimensões é inviável e o erro da simulação não piora com a dimensão.