Physics and Mathematics of Disordered Systems

Fokker-Planck equation for a jump diffusion process

with 7 comments

One of the simplest stochastic processes is a Brownian motion with drift, or a diffusion process:

\displaystyle   \begin{array}{rl}  \dot{x}(t) =& \mu + \sigma\,\xi(t), \\  x(t) =& \mu t + \sigma \,W(t).  \end{array}   (1)

Here, \xi(t) is Gaussian white noise with mean zero and variance \left\langle \xi(t_1)\xi(t_2)\right\rangle = \delta(t_1-t_2). Its integral W(t) = \int_0^t \mathrm{d} t'\,\xi(t') is a Brownian motion.
Continuous ItΓ΄ stochastic processes such as eq. (1) are insufficient for applications where the random variable x may jump suddenly (such as in avalanches). A natural extension of (1) for modelling this is a so-called jump-diffusion process. Let us suppose that our jump sizes s are positive, independent and identically distributed with density Q(s). Then, the jump diffusion process y(t) is

\displaystyle   y(t) = x(t) + \sum_{i=0}^{N(t)} s_i = \mu t + \sigma\, W(t) + \sum_{i=0}^{N(t)} s_i,   (2)

where s_i, i\in \mathbb{N} are i.i.d. jump sizes as above, and N(t) is the number of jumps encountered up to time t. For simplicitly, let us assume that jumps occur independently with rate \lambda, i.e. that the probability to have a jump in a time interval \mathrm{d}t is \lambda\, \mathrm{d}t. Then, N(t) is a Poisson process with rate \lambda.

It is well-known that the diffusion process x(t) in (1) is equivalently described by a partial differential equation for the distribution P(x) of x, the Fokker-Planck equation (FPE)

\displaystyle \partial_t P_t(x) = \frac{\sigma^2}{2}\partial_x^2 P_t(x) - \mu \partial_x P_t(x).   (3)

This representation is useful e.g. for first-passage problems: they correspond to various boundaries introduced in the PDE (3). So how does one generalise the Fokker-Planck (3) to the case of the jump-diffusion process y(t) in (2)? I will explain in the following that the answer is

\displaystyle \partial_t P_t(y) = \frac{\sigma^2}{2}\partial_y^2 P_t(y) - \mu \partial_y P_t(y)- \lambda P_t(y) + \lambda\,\int_0^\infty \mathrm{d}s \,Q(s)P_t(y - s) ,   (4)

and then discuss a specific example.

1. Deriving the jump-diffusion FPE

Let us consider a time step from t to t + \mathrm{d}t. The probability for a jump to occur during this interval is \lambda\,\mathrm{d}t, so

\displaystyle P_{t+\mathrm{d}t}(y) = (1-\lambda \mathrm{d}t)\left\langle P_t\left(y-\mu \,\mathrm{d}t - \left[W(t+\mathrm{d}t)-W(t)\right] \right) \right\rangle_W + \lambda \mathrm{d}t\,\left\langle P_t(y-s) \right\rangle_s,   (5)

where \left\langle \cdot \right\rangle_W denotes averaging over all realizations of the Brownian motion W, and \left\langle \cdot \right\rangle_s denotes averaging over the distribution Q(s) of the jump size s. Since the jump term is already multiplied by the jump probability \propto \mathrm{d}t, the drift and noise contributions there are of higher order in \mathrm{d}t and were dropped.

The averaging over the jump size in (5) yields a convolution with the jump size distribution Q(s):

\displaystyle \left\langle P_t(y-s) \right\rangle_s = \int_0^\infty \mathrm{d}s \,Q(s)P_t(y - s).

The average over the noise W in (5) is the same as for standard diffusion. During the interval \mathrm{d}t, the increment of the noise term in (2) is

\displaystyle W(t+\mathrm{d}t)-W(t) = \int_t^{t+\mathrm{d}t}\mathrm{d}t'\,\xi(t') =: V\sqrt{\mathrm{d}t},

where the last equality is a definition for V.
Since W is a Brownian motion, V is normally distributed:

\displaystyle \mathcal{N}(V) = \frac{1}{\sqrt{2\pi}}\exp\left(-\frac{V^2}{2}\right).

Thus, the average over W in (5) is

\displaystyle  \begin{array}{rl}  \displaystyle \left\langle P_t\left(y-\mu \,\mathrm{d}t - \left[W(t+\mathrm{d}t)-W(t)\right] \right) \right\rangle_W = & \displaystyle \int_{-\infty}^\infty \mathrm{d}V\,\mathcal{N}(V)\, P_t\left(y-\mu \,\mathrm{d}t - V \sqrt{\mathrm{d}t} \right) \\   = & \displaystyle P_t(y) - \mu\, \mathrm{d}t\, \partial_y P_t(y)- \sqrt{\mathrm{d}t} \partial_y P_t(y) \int_{-\infty}^\infty \mathrm{d}V\, \mathcal{N}(V)\, V  \\  & \displaystyle + \frac{1}{2}\mathrm{d}t \partial_y^2 P_t(y)\, \int_{-\infty}^\infty \mathrm{d}V\, \mathcal{N}(V)\, V^2 + \mathcal{O}(\mathrm{d}t)^2 \\  =& \displaystyle P_t(y) - \mu\, \mathrm{d}t\, \partial_y P_t(y) + \frac{1}{2}\mathrm{d}t\,\partial_y^2 P_t(y)  \end{array} .

In the last line we used the fact that the normal distribution has mean zero and variance 1, and dropped all terms of higher order than \mathrm{d}t. Inserting all this into (5), we obtain the jump-diffusion FPE (4).

2. Example: Exponential jumps

Now let us consider a simple example, where the jump sizes are distributed exponentially:

\displaystyle  Q(s) = \alpha\,\exp(-\alpha s),\quad\quad \alpha > 0  .

We will compute the distribution of y(t), starting from y(0)=0.
For this, we solve the jump-diffusion FPE (4) with the initial condition P_{t=0}(y) = \delta(y). Let us take a Fourier transform of (4) from y to k:

\begin{array}{rl}  \displaystyle  \tilde{P}_t(k) := & \displaystyle \int_{-\infty}^\infty \mathrm{d}y\,e^{i k y}\,P_t(y), \\  \displaystyle  \partial_t \tilde{P}_t(k) = & \displaystyle -\frac{\sigma^2}{2}k^2\,\tilde{P}_t(k) + i\mu k\, \tilde{P}_t(k) - \lambda \tilde{P}_t(k) + \lambda\tilde{Q}(k)\tilde{P}_t(k)  \end{array}  .

The different spatial modes decouple, and instead of a partial differential equation we obtain an ordinary differential equation for each value of the Fourier variable k.
For the above exponential form of the jump size distribution Q,

\displaystyle  \tilde{Q}(k) = \frac{\alpha}{\alpha - i k}  .

Furthermore, the initial condition P_{t=0}(y) = \delta(y) gives \tilde{P}_{t=0}(k) = 1. Hence, the solution $\tilde{P}$ reads

\displaystyle  \tilde{P}_t(k) = \exp\left[ \left( -\frac{\sigma^2}{2}k^2+  i\mu k + i\frac{\lambda k}{\alpha - i k} \right) t\right]  .

While it does not seem feasible to invert this Fourier transform to obtain a closed expression for P_t(y) (but if anyone has an idea, let me know!), \tilde{P}_t(y) is already enough to determine the moments of y(t). Taking derivatives, we obtain for example

\displaystyle  \begin{array}{rl}  \left\langle y(t) \right\rangle = & -i\partial_k \big|_{k=0} \tilde{P}_t(k) = \left(\frac{\lambda}{\alpha}+\mu\right)t \\  \left\langle y(t)^2 \right\rangle^c = & -\partial_k^2 \big|_{k=0} \tilde{P}_t(k) + \left[\partial_k \big|_{k=0} \tilde{P}_t(k)\right]^2  \\  = & \left(\frac{2\lambda}{\alpha^2}+\sigma^2\right)t, \\  ... &  \end{array}  .

Similarly, solving the FPE (4) with an absorbing boundary allows computing first-passage times (or at least their moments) for our jump-diffusion process.

Have fun, and do let me know if you have any comments!


Written by inordinatum

November 20, 2013 at 8:04 pm

7 Responses

Subscribe to comments with RSS.

  1. Hi there,
    Any chance you could give some references to literature where you studied this?
    Thanks πŸ™‚


    February 21, 2015 at 11:31 am

    • Hi Alex,

      Glad to see that you’re interested πŸ™‚

      I got the idea for this post when reading the paper by Kou and Wang, First Passage Times of a Jump Diffusion Process. They consider a slightly more complicated case than here, namely a two-sided exponential distribution, and then go beyond the FPE to consider first-passage times. But the first few pages of their paper are also useful as an introduction.

      Another useful reference which uses similar ideas geared towards finance applications is the paper by Kou, A Jump-Diffusion Model for Option Pricing.

      Have fun, and let me know if you find other useful references!



      February 25, 2015 at 9:08 pm

  2. Good work. Is there any published or working paper of yours on this topic I can quote?


    October 30, 2015 at 8:21 pm

  3. Hi there
    I stumbled across your page and found this very useful for a piece of work I’m currently doing. Thank you very much. I do have one question: In the equation right above Section 2, it seems you dropped the d_t_P_t(y) term i.e. the first derivative of P(y) with respect to time. Can I ask why you have done this? Or have you included it among the terms you say are of the order (dt^2)? If so how come? The coefficient of the time derivative should be dt not (dt^2). Hope you see this message. Thanks again


    February 20, 2017 at 1:51 pm

    • Hi Sean,

      glad that you read the post and like it πŸ™‚ In the equation you mean, I’m looking just at the right-hand side of eq. (5) – more precisely the average over the Brownian motion (the part in angular brackets). The time derivative appears when you expand the left-hand side of eq. (5) in terms of \mathrm{d}t (which then leads to the left-hand side of eq. (4) when you collect the terms of order \mathrm{d}t).

      Does that help?


      February 26, 2017 at 3:46 pm

  4. Thanks for the response. yes it does help


    February 28, 2017 at 11:39 pm

    • I had a question. How can I estimate Pt ? I mean is there any algorithm to estimate Pt with jump term given the fooker planck equation?


      April 22, 2017 at 7:57 am

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: