Derivativos de Ações e Commodities¶

Precificação com o Método de Diferenças Finitas¶

_{Uirá Caiado. 29 de Junho, 2016}

Resumo¶

Neste projeto vou implementar o método de diferenças finitas explícito para precificar diferentes derivativos. Como o preço de opções geralmente é descrito por equações diferenciais de difusão ou parabólicas, o método de diferenças finitas se mostra adequado, uma vez que é utilizado justamente para encontrar soluções numéricas para equações diferencias. Depois de implementar o modelo, vou compará-lo com os resultados obtidos pela solução analítica de cada instrumento.

1. Introdução¶

Nesta seção dou uma breve descrição sobre o método de diferenças finitas e declaro o problema que será abordado.

1.1. O Método Utilizado¶

Como colocado por Wilmott, encontrar soluções fechadas para precificação de muitos derivativos pode ser difícil, ou mesmo inviável. Como o preço destes instrumentos frequentemente são descritos por equações diferencias parciais, frequentemente utilizam-se métodos numéricos para encontrar uma solução, como árvores binomiais, simulação de monte carlo e diferenças finitas.

O método de diferenças finitas encontra o valor de um derivativo calculando-o em todo domínio (factivel) de preços do instrumento base, incluíndo a passagem de tempo até seu vencimento. Assim, são similares às árvores binomiais. Porém, ao invéz e discretizar os preços do ativo e a passagem do tempo em uma estrutura de árvore, discretiza em um grid.

Diferenças finitas também são muitos boas para lidar com problemas com poucas dimensões e equações diferenciais não lineares (o preço e tempo já são duas, e o Wilmott sugere até quatro). Para muitas dimensões, a implementação começa a ficar complicada e pouco eficiente, coisa que o Monte Carlo lida melhor. Porém, se há exercício antecipado e não linearidade, o método de diferenças finitas acaba se mostrando como solução mais viável.

1.2. O Problema¶

Considerando um ativo-objeto cuja dinâmica do preço satisfaz a seguinte EDE:

$$\frac{\mathrm{d} S_t}{S_t}=\mu\cdot \mathrm{d}t + \sigma\cdot \mathrm{d}W_{t}$$

Deseja-se calcular o preço justo de um derivativo um derivativo com característica europeia cujo payoff é descrito por uma função qualquer $V_T=V(T, S_T)$, onde $T$ é o vencimento do derivativo. $S_T$ é o preço do ativo-objeto em $T$ e é possível negociar qualquer quantidade dele em qualquer instante. Não há custo de transação (corretagem, emolumento, bid-ask spread, etc) e posições vendidas a descoberto no subjacente são permitidas, não havendo custos associados.

Pede-se que se implemente o algoritmo de diferenças finitas explícito e se calcule o preço justo, o Delta e o Gamma de cada instrumento abaixo. Os valores encontrados devem ser comparados com o resultado de suas expressões analítica. A simulação deve ser feita para os payoffs abaixo. $K$ é o Strike da opção.

• $V(T, S_T)=ln(S_T)$
• $V(T, S_T)=(ln(S_T))^2$
• $V(T, S_T)=(S_T-K)^2$
• $V(T, S_T)=\mathbf{1}_{S_T > K}$
• $V(T, S_T)=max(S_T-K, 0)$

2. Diferenças Finitas¶

Abaixo, apresentarei os elementos necessários para implementar o método diferenças finitas, sendo eles: a diferenciação das derivadas parciais necessárias usando o Grid; a discretização da condição final para cada derivativo; e suas respectivas condições de contorno (nos limites do grid).

2.1. Diferenciação no Grid¶

De acordo com Wilmott, o grid utilizado pelo método de diferenças finitas tem passos de tempo e de preço (ou log do preço) geralmente homogênios. Porém, não há restrição para a forma do Grid, desde que os ajustes necessários sejam feitos. Como exposto por

Considerando que podemos discretizar $t$ como $t=T-k\delta t$ e do $S$ como $S = i\delta S$, onde $i$ e $k$ são seus respectivos passos no grid, podemos escrever o valor da opção em cada ponto do grig como sendo:

$$V^{k}_{i} = ( i \delta S, T-k \delta t )$$

Seguindo notas de aula, aplicando Taylor ao preço $V(T, S_T)$ de um derivativo genérico, posso escrever a seguinte equação diferencial parabólica:

$$\frac{\partial V}{\partial t} + a(S, t) \cdot \frac{\partial V^2}{\partial S^2} + b(S, t) \cdot \frac{\partial V}{\partial S} + c(S, t)\cdot V = 0$$

Sendo que $\frac{\partial V}{\partial t}$ também é chamado de theta ($\theta$), $\frac{\partial V}{\partial S}$ de delta ($\Delta$) e $\frac{\partial V^2}{\partial S^2}$ de gamma ($\Gamma$). Como demostrado por Wilmott, como a definição da primeira derivada $V$ em relação a $t$ é dado por

$$\frac{\partial V}{\partial t} = \underset{h \, \rightarrow \, 0}{\lim}\, \frac{V(S, t + h) - V(S, t)}{h}$$

Então podemos aproximar $\theta$ como sendo $\frac{\partial V}{\partial t} (S, t) \approx \frac{V_{i}^{k} - V_{i}^{k+1} }{\delta t}$. Note que $i$, que é relacionado ao passo do ativo, ficou fixo. Também estou ignorando os erros de aproximação, como ignorarei em todas as outras diferenciações.

O mesmo raciocínio pode ser utilizado para aproximar o $\Delta$. Porém, Wilmott ainda sugere que se utilize a diferença centrada, que é dada por $\frac{\partial V}{\partial S} (S, t) \approx \frac{V_{i+1}^{k} - V_{i-1}^{k} }{2 \delta S}$. A discretização anterior também poderia ser utilizada, mas esta oferece um erro de aproximação menor. O único problema desta abordagem é que preciso saber os valores $S + \delta S$ e $S - \delta S$ e nas fronteiras do Grid não terei esta informação. Por tanto, nestes casos utilizarei a discretrização utilizando apenas um dos lados.

Por fim, para achar o $\gamma$, Wilmott subtrai o delta forward do delta backward e divide o resultado por $\delta S$, chegando na aproximação $\frac{\partial V^2}{\partial S^2} (S, t) \approx \frac{V_{i+1}^{k} - 2 V_{i}^{k} + V_{i-1}^{k} }{\delta S^2}$. Como demonstrado em notas de aula, todos estes resultados podem ser checados aplicando a expansão de Taylor.

2.2. Condição Terminal e Payoffs¶

O valor da opção no vencimento é simplemente seu payoff. Assim, não é necessário resolver nada para $T$, apenas para $S$. Usando a notação de diferenças finitas, os payoffs desejados ficam sendo da forma:

• $V_{i}^{0}=ln(i\delta S)$
• $V_{i}^{0}=(ln(i\delta S))^2$
• $V_{i}^{0}=(i\delta S-K)^2$
• $V_{i}^{0}=\mathbf{1}_{i\delta S > K}$
• $V_{i}^{0}=max(i\delta S-K, 0)$

Como o Wilmott explica, o método de diferenças finitas começa de trás para frente, será destes valores que a iteração começará, como se estivéssemos calculando o preço de um derivativo por árvores binomiais, também de trás para frente.

2.3. Condições de Contorno¶

Quando estivermos percorrendo o Grid, necessitaremos definir o preço do derivativo em seus extremos, quando $S=0$ e $S=I \delta S$, onde $I$ é o ponto mais alto do Grid. Esta condição depende do instrumento que está sendo precificado. Para todos dos contratos, utilizaremos a seguinte condição para quando $S=0$:

$$V^{k}_{0} = (1 - r \delta t)V^{k-1}_{0}$$

Para a condição de contorno superior, utilizaremos para a maioria dos contratos:

$$V^{k}_{I} = 2V^{k}_{I-1} - V^{k}_{I-2}$$

Esta condição é adequada pois, a medida que $S \rightarrow \infty$, o $\Gamma$ da maioria dos contratos tende a zero. Porém, isso não é verdade para o contrato cujo payoff é descrito por $V_{i}^{0}=(i\delta S-K)^2$. Diferenciando novamente o delta deste contrato, chegamos que o gamma será dado por $2 \cdot e^{(r + \sigma^2)(T-t)}$. Neste caso, partindo da equação discretizada do $\Gamma$ e resolvendo para $V^{k}_{I}$, tenho que:

$$\begin{equation} \begin{aligned} \Gamma &= \frac{V_{I+1}^{k} - 2 V_{I}^{k} + V_{I-1}^{k} }{\delta S^2}\\ 2 \cdot e^{(r + \sigma^2)(\delta t)} &= \frac{V_{I+1}^{k} - 2 V_{I}^{k} + V_{I-1}^{k} }{\delta S^2}\\ V_{I}^{k} &= 2 \delta S^2 e^{(r + \sigma^2)(\delta t)} + 2V_{I-1}^{k} - V_{I-2}^{k}\\ V_{I}^{k} &\approx 2 \delta S^2 + 2V_{I-1}^{k} - V_{I-2}^{k}\\ \end{aligned} \end{equation}$$

3. Implementando o Modelo¶

Nesta seção vou decsrever o método de diferenças finitas explícito e implementar os códigos necessários.

3.1. Método Explícito¶

Seguindo notas de aula, substituindo as diferenciações encontradas na equação parabólica mencionada anteriormente e reescrevendo os outros termos com a notação de diferenças finitas, temos que:

$$\begin{equation} \begin{aligned} \frac{\partial V}{\partial t} + a(S, t) \cdot \frac{\partial V^2}{\partial S^2} + b(S, t) \cdot \frac{\partial V}{\partial S} + c(S, t)\cdot V &= 0 \\ \frac{V_{i}^{k} - V_{i}^{k+1} }{\delta t} + a_{i}^{k} \cdot \frac{V_{i+1}^{k} - 2 V_{i}^{k} + V_{i-1}^{k} }{\delta S^2} + \\b_{i}^{k} \cdot \frac{V_{i+1}^{k} - V_{i-1}^{k} }{2 \delta S} + c_{i}^{k}\cdot V_{i}^{k} &= 0 \end{aligned} \end{equation}$$

Rearranjando equação acima para isolar $V_{i}^{k+1}$ e renomeando alguns termos, ficamos com:

$$V_{i}^{k+1} = A_{i}^{k}V_{i-1}^{k} + \left( 1 + B_{i}^{k} \right)V_{i}^{k} + C_{i}^{k}V_{i+1}^{k}$$

Onde:

$$\begin{equation} \begin{aligned} A_{i}^{k} &= \nu_1 a_{i}^{k} - 0.5\nu_2 b_{i}^{k} \\ B_{i}^{k} &= -2\nu_1 a_{i}^{k} - \delta t c_{i}^{k} \\ C_{i}^{k} &= \nu_1 a_{i}^{k} + 0.5\nu_2 b_{i}^{k} \end{aligned} \end{equation}$$

Sendo que $\nu_1=\frac{\delta t}{\delta S ^2}$ e $\nu_2=\frac{\delta t}{\delta S}$. A equação acima está definida apenas entre $i=1, ..., I-1$. OS pontos restantes necessários para discretização vem das condições de contorno. Como nós conhecemos o valor terminal em $V_{i}^{0}$, podemos calcular o valor de $V_{i}^{1}$ e assim por diante. Como o valor do instrumento em $k+1$ só depende dos valores dele em $k$, chamamos este método de método explícito.

3.2. Outros Detalhes¶

a. Estabilidade¶

Para que a solução deste método seja estável, é necessário que satisfaça $\delta t \, \leqslant \frac{\delta S^2}{2\sigma^2 S^2}$ e $\delta S \leqslant \frac{2a}{|b|}$. Vou utilizar a primeira delas explicitamente na implementação, porém vou fazer a seguinte transformação:

$$\delta t \, \leqslant \frac{\delta S^2}{2\sigma^2 S^2} = \frac{1}{\sigma^2} \left(\frac{\delta S}{S}\right)^2 = \frac{1}{\sigma^2 I^2}$$

Onde $I = \delta S \times S$.

b. Interpolação Bilinear¶

Devido a discretização do modelo, pode ser que seja necessário obter um valor que esteja entre os nós do Grid. Nestes casos, adotarei uma interpolação de duas dimensões (Preço e tempo) chamada interpolação bilinear. Ela consiste em realizar uma interpolação linear primeiro em uma direção e depois na outra. Considerando um quadrante específico do Grid, em que há quatro preços de opção, basicamente se divide o quadrante em quatro (onde os pontos de divisão são os valores desejados) e se utiliza a área de cada um deles como peso para os valores das opções. Assim, terminamos com a seguinte fórmula:

$$\frac{\sum_{i=1}^{4} A_i V_i}{\sum_{i=1}^{4} A_i}$$

Onde $A$ é a área dos retângulos e $V$ o preço da opção em cada canto do quadrante principal. Cada $A_i$ é oposto ao $V_i$ correspondente (o $A_2$ é oposto ao $V_2$, por exemplo).

3.3. Criando o Grid¶

Primeiro, vou definir as classes para a criação do Grid. Como quero acessar os valores posteriormente, vou manter todos em uma estrutura, ainda que este provavelmente não seja o método mais adequado pensando em velocidade.

Como comentei, não é uma estrtutura eficiente. Demorou quase 5 segundos para criar uma estrtutura de 10.000 nós. Porém acredito que isso me ajudará adiante.

3.3. Testando Método¶

Seguindo implementação do livro, também há outra opção para se encontrar $V_{i}^{k+1}$. Substituindo a notação da fórmula de precificação e assumindo que possuimos os valores para $\Delta$, $\Gamma$ e $V_{i}^{k}$, podemos isolar o $\theta$, ficando com

$$\begin{equation} \begin{aligned} \theta + a_{i}^{k} \cdot \Gamma + b_{i}^{k} \cdot \Delta + c_{i}^{k}\cdot V_{i}^{k} &= 0\\ -a_{i}^{k} \cdot \Gamma - b_{i}^{k} \cdot \Delta - c_{i}^{k}\cdot V_{i}^{k} &= \theta \end{aligned} \end{equation}$$

Assim, para uma opção que não paga dividendos, que o ativo base segue um movimento brawniano geométrico, e esta sendo avaliado no mundo neutro a risco, temos que o parâmetro $c=-r$, $b=rS$ e $a=0.5 \cdot \sigma^2 S^2$. Assim, fico com o theta sendo da forma:

$$\theta = r V_{i}^{k} - rS_{i}^{k}\Delta - 0.5\sigma^2 S^2\Gamma$$

Então, o preço do derivativo em $k+1$ pode ser dado por:

$$V_{i}^{k+1} = V_{i}^{k} - \delta t \theta$$

Para comparar os resultados obtidos, também serão calculados os preços e as gregas de cada opção analiticamente, seguindo as equações de preficição demostradas em trabalho anterior. Abaixo vou imprimir os resultados obtidos para uma Call Européia para comparação com exemplos do Wilmott.

O título das colunas se refere ao tempo, começando no vencimento (tempo zero) até $t$ (o dia em que a opção está sendo avaliada). As linhas se referem ao preço do ativo objeto. Note que começa a aparecer uma diferença em relação ao livro quando o ativo está valendo $R\$ 110.00$ ou menos. Como o método utiliza os três nós subjacentes anteriores para calcular o valor do nó posterior, o nó $(\,90.0 \,; 0.111)$ não deveria apresentar valor algum de fato. Porém, no exemplo do Wilmott, apresentou.

4. Precificando Instrumentos¶

Nesta seção vou precificar alguns instrumentos e comparar com suas soluções analíticas.

4.1. Call Européia¶

Partindo de implementação anterior, apenas resta imeplementar o gamma analítico da equação, que é dado abaixo:

$$\frac{\partial V^2}{\partial S^2} = \frac{\Delta}{\partial S} = \frac{N'(d1)}{S \sigma \sqrt{t}}$$

Abaixo vou comparar os resultados analíticos com o que obtive do método de diferenças finitas.

4.2. Contrato Log¶

Derivando gamma analítico da equação a paritr do Delta $\Delta = \frac{e^{-r\cdot t}}{S}$, tenho que:

$$\begin{equation} \begin{aligned} \frac{\partial V^2}{\partial S^2} &= \frac{\Delta}{\partial S} = -\frac{e^{-r\cdot t}}{S^2} \end{aligned} \end{equation}$$

Abaixo vou comparar os resultados analíticos com o que obtive do método de diferenças finitas. A única diferença que é necessária tomar nesta implementação é quando o valor do ativo no final é $0$. Nestes casos, fiz o valor da opção também ser $0$.

4.3. Contrato Log Quadrático¶

Derivando gamma analítico da equação e aplicando a regra de quociente, tenho que, se

$$\Delta = \frac{2e^{-r\cdot t}}{S}\left[ln(S)+ \left(r - \frac{\sigma^2}{2} \right) \cdot t \right]$$

Então

$$\begin{equation} \begin{aligned} \frac{\partial V^2}{\partial S^2} &= \frac{\Delta}{\partial S}= 2e^{-r\cdot t} \frac{\partial}{\partial S} \left [ \frac{ln(S)+ t \left(r - \frac{\sigma^2}{2} \right)}{S} \right ]\\ &= 2e^{-r\cdot t} \frac{S \frac{\partial}{\partial S} \left ( ln(S)+ t \left(r - \frac{\sigma^2}{2} \right) \right) - \left( ln(S)+ t \left(r - \frac{\sigma^2}{2} \right) \right) \frac{\partial}{\partial S} (S)}{S^2}\\ &= 2e^{-r\cdot t} \frac{S \cdot \frac{1}{S} - ln(S)- tr + \frac{t \sigma^2}{2}}{S^2}\\ &= \frac{e^{-r t}}{S^2}\left[2 + \sigma^2 t - 2ln(S) -2rt\right] \end{aligned} \end{equation}$$

Abaixo vou comparar os resultados analíticos com o que obtive do método de diferenças finitas.

4.4. Exótico Quadrático¶

Derivando gamma analítico da equação a partir de $\Delta = 2S e^{(r + \sigma^2)} - 2K$, tenho que:

$$\begin{equation} \begin{aligned} \frac{\partial V^2}{\partial S^2} &= 2 e^{(r + \sigma^2)} \end{aligned} \end{equation}$$

Neste caso, também preciso mudar as condições de contorno. Abaixo vou comparar os resultados analíticos com o que obtive do método de diferenças finitas.

4.5. Opção Digital¶

Por último, vou compara o resultado do método explícito com o resultado analítico de uma opção digital. Seguindo a derivação do gamma analítico a daqui, tenho que:

$$\begin{equation} \begin{aligned} \frac{\partial V^2}{\partial S^2} &= - \frac{e^{-rt} d_1 N(d_2) }{S^2 \sigma^2 t } \end{aligned} \end{equation}$$

Abaixo vou plotar os resultados para comparação

A solução não aderiu muito bem às soluções analíticas. Vou tentar melhorar a quantidade de discretizaçòes do ativo subjacente para tentar melhorar este resultado

Curiosamente, a melhorar foi marginal. A implementação que fiz também se mostrou inadequada, uma vez que quando criei um grid com $100$ passos no ativo, demorou 4 minutos e meio para realizar todos os cálculos.

5. Conclusão¶

Neste trabalho o método de diferenças finitas explícito foi comparado à precificação de instrumentos com solução analítica. De maneira geral, o método conseguiu aproximar bem o preço das opções e se mostrou bastante flexível, precisando de poucos ajustes para precificar cada instrumento. Apesar de não testado aqui, provavelmente a inclusão de não linearidades nos contratos também não exigiria muitos ajustes, como a utilização de parâmetros incertos, ou mesmo exercício antecipado.

6. Últimas Considerações¶

A decisão de usar classes no lugar de um array dificultou a implementação. A aplicação do método criando um objeto para cada nó do Grid também se mostrou muito lenta. Porém, foi interessante notar que o simples aumento da discretização do grid não melhorou tanto a aderência aos valores calculados pelas soluções analíticas. Talvez a mudança do formato do grid mostrasse melhores resultados.

Style notebook and change matplotlib defaults