Z Skrypty dla studentów Ekonofizyki UPGOW

Stochastyczne równania różniczkowe

Google

W tym rozdziale zostaną opisane metody numeryczne, które służa do rozwiązywania stochastycznych równań różniczkowych typu:

\(\frac{dX(t)}{dt} = F(X(t), t) + G(X(t), t)\Gamma(t)\)

gdzie F i G to dowolne funkcje, a \(\Gamma(t)\) jest procesem losowym. Najczęstszym przypadek to taki w którym \(\Gamma(t)\) to biały szum Gaussowski. Tak zapisane równanie nie jest precyzyjnie określone ze względu na dylemat Stratonowicza-Ito. Dlatego poprawne jest zapisanie równanie Ito w postaci:

\(dX(t)= F(X(t), t)dt + G(X(t), t) dW(t)\;\)

Nie zmienia to ogólności, gdyż jak wiadomo każde równanie zapisane w interpretacji Stratonowicza ma swój odpowiednik Ito. Dla potrzeb metod numerycznych będziemy rozpatrywać zawsze równania Ito, a jeśli pojawią się równania Stratonowicza to będziemy je transpormować do postaci Ito.

Schemat Eulera dla równań stochastycznych

Najprostszą metodą aproksymacji numerycznej równania stochastycznego jest podobnie jak w przypadku równań różniczkowych zwyczajnych jest schemat Eulera. Część deterministyczną równania stochastycznego traktujemy w taki sam sposób jak w schemacie Eulera dla równań różniczkowych zwyczajnych. Niech h oznacza krok całkowania i oś czasowa będzie zdyskretyzowana na przedzialy \(t_{i-1},t_{i},t_{i+1}\) oraz \(h=t_{i}-t_{i-1}\). Wtedy część deterministyczna równania stochastycznego przyjmuje postać:

\(X(t_i) = X(t_{i-1}) + \int_{t_{i-1}}^{t_{i}} F(X(t), t) dt \simeq X(t_{i-1}) + F(X(t_{i-1}, t_{i-1}) h \)

Aby całkować część stochastyczną potrzebujemy formuły na przyrost skończony procesu Wienera:

\( \int_{t_{i-1}}^{t_{i}} \Gamma(t) dt =\;\int_{t_{i-1}}^{t_{i}} dW(t) = W(t_{i})-W(t_{i-1})\)

Wiemy, że proces Wienera jest procesem o przyrostach niezależnych, które są gaussowską zmienna losową o zerowej wartości średniej i wariancji \(2(t_{i} − t_{i-1})=2 h\).

Tak więc widać, że w schemacie Eulera całkę typu \( \int_{t_{i-1}}^{t_{i}} \Gamma(t) dt \) należy zastąpić w każdym kroku całkowania gaussowską zmienną losową o wariancji proporcjonalnej do kroku całkowania h. Ponieważ z reguły dysponujemy gaussowskich generatorem liczb losowych o jednostkowej wariancji N(0,1), można zapisać:

\( \int_{t_{i-1}}^{t_{i}} \Gamma(t) dt = \sqrt{2 h} N(0,1) \)

Ten zapis pokazuje też ważną cechę przy obliczaniu aproksymacji rozwiązań stochastycznych - najniższy rząd w h jest nie O(h) ale \(O(h^{1/2})\).

Ponadto z takiego sformułowania widać też, że zmiany procesu Wienera w stosunku do przyrostów czasu są rozbieżne w granicy \(h\to 0\).

Korzystając w powyższych faktów, możemy zapisać pełny schemat Eulera dla równania stochastycznego (Ito):

\(X(t_i) = X(t_{i-1}) + F(X(t_{i-1}, t_{i-1}) h  + \sqrt{h}G(X(t_{t-1}), t_{t-1})  N(0,1)\).

Gaussowskie zmienne losowe możemy otrzymać np. korzystając z algorytmu Box-a Mullera.

Schemat Eulera dla układu równań stochastycznych

Schemat Eulera można uogólnić na układy równań stochastycznych. Niech \(\mathbf X(t)\) będzie wektorem o składowych \[\mathbf X(t)=( X^1(t),X^2(t),...,X^n(t))\]. Układ równań stochastycznych (Ito) można zapisać w ogólnej postaci:

\[d X^i(t)= F^i(\mathbf X(t), t)dt + \sum_{j=1}^{n} G^{i,j}(\mathbf X(t), t) dW^j(t)\;\]\; j=1,2,...,n,

gdzie \(W^i(t),\;W^j(t)\) są niezależnymi procesami Wienera dla \(i\neq j\), \(F^i\) oznacza wektor drytfu a \(G^{i,j}\) jest macierzą \(n \times n\) funkcji.

Wtedy schemat Eulera ma postać:

\( X^j(t_i) = X^j(t_{i-1}) + F^j (\mathbf X(t_{i-1}), t_{i-1}) h + \sqrt{h} \sum_{k=1}^{n} G^{j,k}(\mathbf X(t), t) N^k(0,1)\; \)

Proces Wienera

Proces Wienera jest rozwiązaniem następującego stochastycznego równania różniczkowego:

\(dX(t)= 2 D dW(t)\;\).

Jego realizacja jest funkcją ciągłą, ale jednocześnie nigdzie nie jest różniczkowalna.

Proces Wienera

Stosując schemat Eulera można wygenerować pojedyńczą trajektorię takiego procesu:

h=0.01;
N=400;
x(1)=0;
D=1;
for i=2:N
  x(i)=x(i-1)+ sqrt(2*D*h)* normrnd (0,1,1,1);
endfor
 
plot((1:N)*h,x,'-')