Euler–Maruyama method

In mathematics, the Euler–Maruyama method is a method for the approximate numerical solution of a stochastic differential equation (SDE). It is a simple generalization of the Euler method for ordinary differential equations to stochastic differential equations. It is named after Leonhard Euler and Gisiro Maruyama.

Consider the stochastic differential equation (see Itō calculus)

$\mathrm{d} X_t = a(X_t) \, \mathrm{d} t + b(X_t) \, \mathrm{d} W_t,$

with initial condition X0 = x0, where Wt stands for the Wiener process, and suppose that we wish to solve this SDE on some interval of time [0, T]. Then the Euler–Maruyama approximation to the true solution X is the Markov chain Y defined as follows:

• partition the interval [0, T] into N equal subintervals of width $\Delta t>0$:
$0 = \tau_{0} < \tau_{1} < \cdots < \tau_{N} = T \mbox{ and } \Delta t = T/N;$
• set Y0 = x0;
• recursively define Yn for 1 ≤ n ≤ N by
$\, Y_{n + 1} = Y_n + a(Y_n) \Delta t + b(Y_n) \Delta W_n,$
where
$\Delta W_{n} = W_{\tau_{n + 1}} - W_{\tau_n}.$

The random variables ΔWn are independent and identically distributed normal random variables with expected value zero and variance $\Delta t$.

Example

The following Python code implements Euler–Maruyama to solve the Ornstein–Uhlenbeck process:

import numpy as np
import matplotlib.pyplot as plt
import math
import random

tBegin=0
tEnd=2
dt=.00001

t = np.arange(tBegin, tEnd, dt)
N = t.size
IC=0
theta=1
mu=1.2
sigma=0.3

sqrtdt = math.sqrt(dt)
y = np.zeros(N)
y[0] = IC
for i in xrange(1,N):
y[i] = y[i-1] + dt*(theta*(mu-y[i-1])) + sigma*sqrtdt*random.gauss(0,1)

ax = plt.subplot(111)
ax.plot(t,y)
plt.show()