MKZR:Modelowanie dynamiki instrumentów pochodnych

Z Skrypty dla studentów Ekonofizyki UPGOW

(Różnice między wersjami)
m (Model Blacka-Scholesa dla europejskiej opcji kupna)
Linia 1: Linia 1:
-
=== Geometryczny proces Wienera ===
 
-
 
-
 
-
Geometryczny proces Wienera znajduje zastosowanie w ekonomii jako model ceny akcji (patrz: [[PIZL:Przykłady zastosowań równań stochastycznych w ekonomii|Przykłady zastosowań równań stochastycznych w ekonomii]]).
 
-
 
-
 
-
 
-
==== Definicja ====
 
-
Geometryczny proces Wienera jest procesem losowym, który jest rozwiązaniem równania stochastycznego Ito:
 
-
 
-
<math>dX(t) = \mu X(t) dt + \sigma X(t) d W(t)\,</math>.
 
-
 
-
[[Plik:GeomBM.png|thumb|360px|Geometryczny proces  Wienera, dwadzieścia trajektorii. Widoczny jest eksponencjalny wzrost zmiennej losowej z pewnymi fluktuacjami.]]
 
-
Część deterministyczna część tego równania jest  liniowa, podobnie jak w przypadku  [[Ornsteina-Uhlenbecka]]. Wprzypadku <math>\sigma=0</math> rozwiązanie przedstawia exponencjalne rosnącą funkcję czasu:
 
-
 
-
<math>x(t)\simeq e^{\mu t}</math>,
 
-
 
-
Można się spodziewać, że zabużenie równania szumem będzie odwzierciedlało tą cechę.
 
-
 
-
==== Symulacja numeryczna ====
 
-
 
-
Człon stochastyczny jest proporcjonalny o wartości procesu, czyli mamy do czynienia z tzw. "szumem multiplikatywnym". Ponieważ równanie jest w interpretacji  Ito, można je bezpośrednio rozwiązać numerycznie. Interpretacja Ito dla <math>\displaystyle \sigma X(t) d W(t)</math> oznacza, że w schemacie aproksymacyjnym bierzemy wartość procesu X(t) "przed skokiem". W takim przypadku możemy zastosować  schemat numeryczny Eulera podobnie jak np. dla procesu Ornsteina-Uhlenbecka. Krok czasowy może być zapisany w postaci wektorowej jako:
 
-
 
-
<source lang='matlab'>
 
-
  x(:,i)=x(:,i-1) + mu*x(:,i-1)*h + sigma*sqrt(h)*x(:,i-1).*normrnd (0,1,M,1);
 
-
</source>
 
-
 
-
Jako warunek początkowy dla symulacji należy przyjąć wartość x>0. Łatwo zauważyć własność równania definiującego proces, że startując w <math>x(0)=0</math> rozwiązaniem jest funkcja stała  <math>\displaystyle x(t)=0</math>. Z drugiej strony równanie jest symetryczne ze względu transformację <math>x\to -x</math> i ponadto jeśli x oznacza cenę to powinno mieć wartość dodatnią. Dlatego zawsze bierzemy x>0 jako warunek początkowy.  Poniższy program generuje 20 trajektori geometrycznego procesu Wienera:
 
-
 
-
<source lang='matlab'>
 
-
clear all
 
-
close all
 
-
N=400;
 
-
M=20;
 
-
T=14;
 
-
h=T/N;
 
-
clear x
 
-
x=zeros(M,N);
 
-
x(:,1)=1*ones(M,1); # log(1)=0
 
-
sigma=.1;
 
-
mu=0.1;
 
-
for i=2:N
 
-
  x(:,i)=x(:,i-1) + mu*x(:,i-1)*h + sigma*sqrt(h)*x(:,i-1).*normrnd (0,1,M,1);
 
-
endfor
 
-
plot((1:N)*h,x,'r-')
 
-
</source >
 
-
 
-
==== Rozwiązanie analityczne: porównanie z symulacją  ====
 
-
 
-
Równanie definiujące geometryczny proces Wienera można przetransformować do równania na proces Wienera z dryfem. Transformacją jest:
 
-
 
-
<math>\displaystyle Y=\log(X)</math>.
 
-
 
-
Łatwo to zauważyć dzieląc równanie przez <math>x(t)</math> i korzystając z faktu, że
 
-
 
-
<math>\displaystyle d log(X)= dX/X</math>
 
-
 
-
i  korzystając z
 
-
[[PIZL:Stochastyczne_r%C3%B3wnania_r%C3%B3%C5%BCniczkowe#Rachunek_r.C3.B3.C5.BCniczkowy_Ito|formuły  Ito]] na dla funkcji log(x) otrzymujemy:
 
-
 
-
otrzymujemy:
 
-
 
-
:<math>d \log(X(t)) = (\mu-\frac{1}{2}\sigma^2) dt + \sigma d W(t)\,</math>.
 
-
 
-
Tak więc geometryczny proces Wienera jest równoważny procesowi Wienera z dryfem dla <math>log(X(t))</math>. Gęstośc pradwopodobieństwo  przejścia ze stanu <math>x(t_0)</math> do <math>x(t)</math> w czasie <math>t-t_0</math> dla procesu Wienera  z dryfem <math>Y(t)=log(X(t))</math> wynosi:
 
-
 
-
:<math>P_y(y,t|y_0,t_0)= \frac{1}{\sqrt{2\pi\sigma^2 t}} e^{-\displaystyle\frac{(y-(\mu-1/2\sigma^2)^2}{2\sigma^2t}}</math>
 
-
 
-
Korzystając ze reguł transformacji zmiennych losowych dla funkcji log, która jest różnowartościowa:
 
-
 
-
<math>P_y(y)=\frac{P_x(x)}{|g'(x)|}</math>
 
-
 
-
mamy:
 
-
 
-
<math>P_x(x)= \frac{P_y(\log(x))}{|x|}</math>
 
-
 
-
tak więc ostatecznie otrzymujemy formułę:
 
-
 
-
:<math>P_x(x,t|x_0,t_0)= \frac{1}{\sqrt{2\pi\sigma^2 t x^2}} e^{-\displaystyle\frac{(\log(x/x_0)-(\mu-1/2\sigma^2)^2}{2\sigma^2t}}</math>.
 
-
 
-
Zmienna losowa w geometrycznym procesie Wienera po czasie <math>t</math> ma rozkład [http://pl.wikipedia.org/wiki/Rozkład_logarytmicznie_normalny logarytmicznie normalny].
 
-
[[Plik:GeomBM_hist.png|thumb|360px|Rozkład P(x,t) dla geometrycznego procesu  Wienera w czterech następujących po sobie chwilach t =  0.7,3.5,6.2,9.1. Widać rozmywanie się piku zgodnie z odpowiednim rozkładem [http://pl.wikipedia.org/wiki/Rozkład_logarytmicznie_normalny logarytmicznie normalnym].]]
 
-
Porównajmy ten wzrór z histogramami dużej licznby realizacji w kolejnych chwilach czasu. 
 
-
W systemie GNU Octave/matlab funkcja rozkładu logarytmicznie normalnego jest dana jako '''normrnd'''. W poniższym programie startujemy z punktu <math>x(t=0)=1</math> dla którego <math>log(x)=0</math> czyli odpowiedni process Wienera z dryfem miałby warunek początkowy <math>y(t=0)=0</math>.
 
-
 
-
 
-
 
-
 
-
 
-
 
-
<source lang='matlab'>
 
-
clear all
 
-
close all
 
-
N=100;
 
-
M=20000;
 
-
T=14;
 
-
h=T/N;
 
-
clear x
 
-
x=zeros(M,N);
 
-
x(:,1)=1*ones(M,1); # log(1)=0
 
-
sigma=.1;
 
-
mu=0.1;
 
-
for i=2:N
 
-
  x(:,i)=x(:,i-1) + mu*x(:,i-1)*h + sigma*sqrt(h)*x(:,i-1).*normrnd (0,1,M,1);
 
-
endfor
 
-
 
-
hold on
 
-
for idx=1:4
 
-
  n=5+(idx-1)*20;
 
-
  t=n/N*T
 
-
  xmax=5;
 
-
  h1=.1;
 
-
  hist(x(:,n),[0:h1:xmax],1/h1)
 
-
  fplot(@(xx) lognpdf(xx,(mu-0.5*sigma^2)*t,sigma*sqrt(t)),[0.0,xmax],200,'ro-')
 
-
endfor
 
-
hold off
 
-
</source >
 
===Wycena opcji: modele z czasem dyskretnym===
===Wycena opcji: modele z czasem dyskretnym===

Wersja z 12:42, 7 cze 2010

Wycena opcji: modele z czasem dyskretnym

Binomial options pricing model


rozdział 7 P.B.

http://gillesdaniel.com/natlab/

Wycena opcji: modele z czasem ciągłym

Model Blacka-Scholesa dla europejskiej opcji kupna

Model Blacka-Scholesa zakłada że, że cena akcji \(\displaystyle S(t)\) ewoluuje zgodnie z geometrycznym ruchem Browna przy ciągłej kapitalizacji z roczną stopą procentową r:

\(dS(t) = r S(t) dt + \sigma S(t) d W(t)\,\).

W chwili \(t_0\) instument bazowy ma wartość \(S0\) a cena wykupu (lub sprzedaży) w chwili T wynosi S(T). Pytamy się o wartość instrumentu pochodnego w dowonej chwili \(t<T\). Zakładając, że chcemy utrzymać portfel bez ryzyka (delta neutralny) można pokazać, że wartość taka jest dana przez równanie Blacka-Scholesa

\( \frac{\partial V}{\partial t} + \frac{1}{2}\sigma^2 S^2\frac{\partial^2 V}{\partial S^2} + rS\frac{\partial V}{\partial S} - rV = 0. \)

gdzie r jest stopą procentową wolną od ryzyka.

Równanie to posiada analityczne rozwiązanie i jest to słynny wzór Blacka-Scholesa dla ceny opcji kupna na europejska akcję bez dywidend:

\[ C(S,t) = SN(d_1) - Ke^{-r(T - t)}N(d_2) \, \]

i ceny opcji sprzedaży na europejska akcję bez dywidend:

\[ P(S,t) = Ke^{-r(T-t)} - S + C(S,t). \ \]

gdzie:

\[ d_1 = \frac{\ln(\frac{S}{K}) + (r + \frac{\sigma^2}{2})(T - t)}{\sigma\sqrt{T - t}} \]
\[ d_2 = d_1 - \sigma\sqrt{T - t}. \]

Powyższe wzory można zaimplementować jako funkcje w matlab/GNU Octave w następujący sposób:

function [C,P] = BlackScholes(S0,K,r,T,sigma)
  d1=(log(S0./K)+(r+sigma^2/2)*T)/(sigma*sqrt(T));
  d2=d1-sigma*sqrt(T);
  C = S0.*normcdf(d1)-K*exp(-r*T)*normcdf(d2);
  P = K*exp(-r*T)*normcdf(-d2)-S0.*normcdf(-d1);
end
Wykres ceny europejskiej opcji sprzedaży dla K=51.11,r=0.1,T=5.0 oraz \(\sigma\)=0.2 od ceny aktualnej S (linia ciągła). Kropkami zaznaczono funkcje \(S-K e^{-r T}\)

Program ten zwraca dla zadanych wartości ceny aktualnej ("spot") S0, cena rozliczenia opcji K w czasie T oraz wysokość rocznej stopy procentowej wolnej od ryzyka dla terminu wygaśnięcia opcji r i volatility \(\sigma\) dwie ceny - opcji kupna C oraz sprzedaży P. Są one zwracane jako dwa argumenty, więc poprawne wywołanie takiej funkcji powinno być np.:

octave:8> [call,put]=BlackScholes(31 ,51.1 ,0.1 ,5 , 0.2)
call =  5.4876
put =  5.4813

Wzór ten ma kilka ciekawych własności. Po pierwsze zauważmy, że gdyby fluktuacje ceny akcji znikły w chwili obserwacji to wartość opcji była by równa różnicy między ceną aktualną a zdyskontowaną ceną wykupu \(S-K e^{-r T}\), gdy tylko taka różnica była by większa od zera. W innym przypadku wartość takiej opcji jest równa zeru. Na rysunku jest umieszczona zależności ceny wykupu od wartości aktualnej opcji wygenerowany następującym skryptem:

clear all
close all
T=5.
S0=31
sigma=0.2
r = 0.1
K = S0*exp(r*T)
S0_tab = 1:1:60;
plot(S0_tab,BlackScholes(S0_tab,K,r,T,sigma), S0_tab, max(S0_tab-K*exp(-r*T),0),'o' );
exact=BlackScholes(S0,K,r,T,sigma)
line([S0,S0],[0,exact])
grid on


Można przypuszczać, że wykres wzoru Blacka Scholes-a będzie dążył do kropkownej zależności dla \(\sigma \to 0\). Rzeczywiście, można się o tym przekonać rysując kilka zależności C(S) dla malejących wartości \(\sigma\).

Wykres ceny europejskiej opcji sprzedaży dla K=51.11,r=0.1,T=5.0 od ceny aktualnej S dla kilku wartości \(\sigma\)=0.01,0.1,0.2 i 0.3




clear all
close all
N=303;
M=32313;
T=13.1;
h=T/N;
x=zeros(M,1);
 
S0=11;
sigma=0.4;
K = 50;
r = 0.1;
 
x=S0*ones(M,1); 
 
% integrate in one step  
y = x .* exp(  (r-0.5*sigma^2)*T + sigma.*sqrt(T).*normrnd (0,1,M,1) ) ;
 
% integrate in N  steps mu=>r
for i=2:N
  x =x + r*x*h + sigma*sqrt(h)*x.*normrnd (0,1,M,1);
endfor
 
 
%exp(-r*T)*mean( max(K-x,0) )
 
% Call option value
Nsteps=exp(-r*T)*mean( max(x-K,0) )
onestep=exp(-r*T)*mean( max(y-K,0) )
exact=BlackScholes(S0,K,r,T,sigma)
% Black Scholes, własności
S0 = 1:1:80;
K = 55;
r = 0.08;
sigma = 0.1;
T=5.01
plot(S0,BlackScholes(S0,K,r,T,sigma), S0, S0-K*exp(-r*T),'o' );
axis([0 80 -5 35]);
grid on