View source for Runge–Kutta method (SDE)

In [[mathematics]] of stochastic systems, the '''Runge–Kutta method''' is a technique for the approximate [[numerical analysis|numerical solution]] of a [[stochastic differential equation]]. It is a generalisation of the [[Runge–Kutta methods|Runge–Kutta method]] for [[ordinary differential equation]]s to stochastic differential equations (SDEs).  Importantly, the method does not involve knowing derivatives of the coefficient functions in the SDEs.

==Most basic scheme==

Consider the [[Itō diffusion]] <math>X</math> satisfying the following Itō stochastic differential equation

:<math>{{d} X_{t}} = a(X_{t}) \, {d} t + b(X_{t}) \, {d} W_{t},</math>

with [[initial condition]] <math>X_0=x_0</math>, where <math>W_t</math> stands for the [[Wiener process]], and suppose that we wish to solve this SDE on some interval of time <math>[0,T]</math>. Then the basic '''Runge–Kutta approximation''' to the true solution <math>X</math> is the [[Markov chain]] <math>Y</math> defined as follows:<ref>P. E. Kloeden and E. Platen. ''Numerical solution of stochastic differential equations'', volume 23 of Applications of Mathematics. Springer--Verlag, 1992.</ref>

* partition the interval <math>[0,T]</math> into <math>N</math>  subintervals of width <math>\delta=T/N>0</math>:

:<math>0 = \tau_{0} < \tau_{1} < \dots < \tau_{N} = T;</math>

* set <math>Y_0:=x_0</math>;
* recursively compute <math>Y_n</math> for <math>1\leq n\leq N</math> by

:<math>Y_{n + 1} := Y_{n} + a(Y_{n}) \delta + b(Y_{n}) \Delta W_{n} + \frac{1}{2} \left( b(\hat{\Upsilon}_{n}) - b(Y_{n}) \right) \left( (\Delta W_{n})^{2} - \delta \right) \delta^{-1/2},</math>

where <math>\Delta W_{n} = W_{\tau_{n + 1}} - W_{\tau_{n}}</math>
and <math>\hat{\Upsilon}_{n} = Y_{n} + a(Y_n) \delta + b(Y_{n}) \delta^{1/2}.</math>
The [[random variables]] <math>\Delta W_{n}</math> are [[independent and identically distributed]] [[normal distribution|normal random variables]] with [[expected value]] zero and [[variance]] <math>\delta</math>.

This scheme has strong order 1, meaning that the approximation error of the actual solution at a fixed time scales with the time step <math>\delta</math>. It has also weak order 1, meaning that the error on the statistics of the solution scales with the time step <math>\delta</math>. See the references for complete and exact statements.

The functions <math>a</math> and <math>b</math> can be time-varying without any complication. The method can be generalized to the case of several coupled equations; the principle is the same but the equations become longer.

==Variation of the Improved Euler is flexible==

A newer Runge—Kutta scheme also of strong order 1 straightforwardly reduces to the Improved Euler scheme for deterministic ODEs.
<ref>A. J. Roberts. Modify the improved Euler scheme to integrate stochastic differential equations. [https://arxiv.org/abs/1210.0933], Oct 2012.</ref>  
Consider the vector stochastic process <math>\vec X(t)\in \mathbb R^n</math> that satisfies the general Ito SDE
:<math>
d\vec X=\vec a(t,\vec X)\,dt+\vec b(t,\vec X)\,dW,
</math>
where drift <math>\vec a</math> and volatility <math>\vec b</math> are sufficiently smooth functions of their arguments.
Given time step <math>h</math>, and given the value <math>\vec X(t_k)=\vec X_k</math>, estimate <math>\vec X(t_{k+1})</math> by <math>\vec X_{k+1}</math> for time <math>t_{k+1}=t_k+h</math> via
:<math>
\begin{array}{rl}
&\vec K_1=h\vec a(t_k,\vec X_k)+(\Delta W_k-S_k\sqrt h)\vec b(t_k,\vec X_k),
\\&\vec K_2=h\vec a(t_{k+1},\vec X_k+\vec K_1)+(\Delta W_k+S_k\sqrt h)\vec b(t_{k+1},\vec X_k+\vec K_1),
\\&\vec X_{k+1}=\vec X_k+\frac12(\vec K_1+\vec K_2),
\end{array}
</math>
* where <math>\Delta W_k=\sqrt hZ_k</math> for normal random <math>Z_k\sim N(0,1)</math>; 
* and where <math>S_k=\pm1</math>, each alternative chosen with probability <math>1/2</math>.

The above describes only one time step.
Repeat this time step <math>(t_m-t_0)/h</math> times in order to integrate the SDE from time <math>t=t_0</math> to <math>t=t_m</math>.

The scheme integrates Stratonovich SDEs to <math>O(h)</math> provided one sets <math>S_k=0</math> throughout (instead of choosing <math>\pm 1</math>).

==Higher order Runge-Kutta schemes==

Higher-order schemes also exist, but become increasingly complex.
Rößler developed many  schemes for Ito SDEs,<ref>{{Cite journal | doi = 10.1137/060673308| title = Second Order Runge–Kutta Methods for Itô Stochastic Differential Equations| journal = SIAM Journal on Numerical Analysis| volume = 47| issue = 3| pages = 1713–1738| year = 2009| last1 = Rößler | first1 = A. }}</ref><ref>{{Cite journal | doi = 10.1137/09076636X| title = Runge–Kutta Methods for the Strong Approximation of Solutions of Stochastic Differential Equations| journal = SIAM Journal on Numerical Analysis| volume = 48| issue = 3| pages = 922–952| year = 2010| last1 = Rößler | first1 = A. }}</ref>
whereas Komori developed schemes for Stratonovich SDEs.<ref>{{Cite journal | doi = 10.1016/j.apnum.2006.02.002| title = Multi-colored rooted tree analysis of the weak order conditions of a stochastic Runge–Kutta family| journal = Applied Numerical Mathematics| volume = 57| issue = 2| pages = 147–165| year = 2007| last1 = Komori | first1 = Y. | url = https://kyutech.repo.nii.ac.jp/?action=repository_action_common_download&item_id=1627&item_no=1&attribute_id=17&file_no=1}}</ref><ref>{{Cite journal | doi = 10.1016/j.cam.2006.03.010| title = Weak order stochastic Runge–Kutta methods for commutative stochastic differential equations| journal = Journal of Computational and Applied Mathematics| volume = 203| pages = 57–79| year = 2007| last1 = Komori | first1 = Y. }}</ref><ref>{{Cite journal | doi = 10.1016/j.cam.2006.06.006| title = Weak second-order stochastic Runge–Kutta methods for non-commutative stochastic differential equations| journal = Journal of Computational and Applied Mathematics| volume = 206| pages = 158–173| year = 2007| last1 = Komori | first1 = Y. }}</ref> Rackauckas extended these schemes to allow for adaptive-time stepping via Rejection Sampling with Memory (RSwM), resulting in orders of magnitude efficiency increases in practical biological models<ref>{{Cite journal | doi = 10.3934/dcdsb.2017133| title = ADAPTIVE METHODS FOR STOCHASTIC DIFFERENTIAL EQUATIONS VIA NATURAL EMBEDDINGS AND REJECTION SAMPLING WITH MEMORY| journal = Discrete and Continuous Dynamical Systems - Series B| volume = 22| issue = 7| pages = 2731–2761| year = 2017| last1 = Rackauckas| first1 = Christopher | last2 = Nie | first2 = Qing }}</ref>, along with coefficient optimization for improved stability<ref>{{cite arXiv |last1 = Rackauckas| first1 = Christopher | last2 = Nie | first2 = Qing |date=2018|title=Stability-optimized high order methods and stiffness detection for pathwise stiff stochastic differential equations |eprint=1804.04344}}</ref>.

==References==
{{reflist}}

{{DEFAULTSORT:Runge-Kutta Method (Sde)}}
[[Category:Numerical differential equations]]
[[Category:Stochastic differential equations]]