Welcome to THETAWIKI. If you like to create or edit a page please make sure to login or register an account. All registered users please make sure to provide a valid email address.

# Jump Diffusion

## The Merton Jump-Diffusion Model

Merton (1976)  jump-diffusion model is an extension to the Geometric Brownian Motion model, with the underlying asset exhibit jumps in addition to have continuous diffusion paths. The asset price evolves as

$d S_t = \mu S_t d t \ + \ \sigma S_t \ d W_t \ + \ (\eta -1)d q$

where $\mu$ is the drift rate, $\sigma$ the volatility of $S$ and $d W_t$ the increment of a Wiener process. The independent Poisson process $d q$ has a value of $1$ with probability $\lambda d t$ and 0 otherwise.

The following is a ThetaML implementation of the Merton Jump-Diffusion Model.

## Merton jump diffusion Process simulation

### Implementation as ThetaML

model JumpDiffusionProcess
%Jump diffusion model for stock price movement
%Example input parameters:
%S0 = 100, r = 0.03, sigma = 0.2, lambda = 1
%nu = -0.1, delta = 0.05, t = 1, ht = 1/12
import S0      "Initial stock price"
import r       "Risk-free interest rate"
import sigma   "Stock price volatility"
import lambda  "Jump intensity"
import nu      "Mean of jumps"
import delta   "Std. of jumps"
import t       "Simulation time horizon"
import ht      "Discretization time step"
export S       "Simulated stock prices"

S = S0         %initialize stock prices

loop t/ht      %loop t/ht times
theta ht   %theta advances time by ht units
P = poissrnd(lambda*ht,size(rand()))
U = exp(P*nu + sqrt(P)*delta*randn())
S = S * exp((r-lambda*(exp(nu+0.5*delta^2)-1)-0.5*sigma^2)*ht + sigma*sqrt(ht)*randn())*U
end
end

Note that this implementation requires the Matlab statistics toolbox. The model JumpDiffusionProcess uses Matlab function poissrnd to sample a poission distributed random variable.

Here is an implementation as a Theta.m function without using the Matlab statistics toolbox:

### Implementation as Theta.m

function state = Theta(dt, state)

if nargin == 0
state.S.value = 100;
state.EUR.value = 1;
else
r = 0.05;

mu1 = -0.2;
lambda = 0.1;
sigma = 0.2;
gamma = 0.45;

mu = r -(exp(mu1 + 0.5*gamma^2) - 1) * lambda;

max_dt = lambda;
n = max(dt/max_dt, 1);

sub_dt = dt/n;
dq = zeros(size(state.S));
for i = 1 : n
jumpsize = mu1 + gamma*randn(size(state.S));%
u = rand(size(state.S));
dq(u<lambda*sub_dt) = dq(u<lambda*sub_dt) + jumpsize(u<lambda*sub_dt);
end

state.S = state.S .* exp((mu-0.5*sigma^2)*dt + sqrt(dt)*sigma*randn(size(state.S))+dq);
state.EUR = state.EUR* exp(-r*dt);
end
end

## Analytic pricing of European calls

European call and put options with the underlying following a Merton jump-diffusion process have analytical solutions. Below is a Matlab implementation of the analytical solution for European calls.

function c = CPMertonJD
% function c = CPMertonJD calculates the analytical price
% of a European call under jump diffusion.

So = 100;     % Initial Stock price
X = 100;      % Strike
T = 1;        % Time to maturity
delta = 0.4;  % std of lognormal jump process
nu = 0;       % mean of lognormal jump process
lambda = 5;   % intensity of jumps
vola = 0.4;   % vola of stock price
r = 0.05;     % interest rate

% The price is gven as a series of terms,
% compute the first N terms for the result
N = 50;
K = exp(nu + 0.5*delta^2) - 1;
c = 0;
for n = 0 : N
sigma_n = sqrt(vola^2 + n*delta^2/T);
r_n = r - lambda*K + n*log(1+K)/T;
d1 = (log(So/X) + (r_n + 0.5*sigma_n^2)*T)/(sigma_n*sqrt(T));
d2 = d1 - sigma_n*sqrt(T);
f_n = So * normcdf(d1,0,1) - X * exp(-r_n*T) * normcdf(d2,0,1);
c = c + exp(-lambda*(1+K)*T) * (lambda*(1+K)*T)^n * f_n/factorial(n);
end
end

Some observations about Monte Carlo convergence:

• The error decreases with increasing number of paths.
• The error is less if the variance of the jump-size is less.
• If the jump intensity (lambda) is high, then the error decreases at a smaller rate with increasing number of paths.