# Heston model with jumps in both variance and underlying dynamic

How can I build on Matlab a Heston model using characteristic function adding jumps in both variance and underlying dynamic ? Suppose that the number of jumps is Poisson-distributed but the jump size is log-normal for the underlying jumps and exponential for the variance. This is the current code in Matlab, so how can I add jumps?

function y = call_heston_cf(s0, v0, vbar, a, vvol, r, rho, t, k)

%   Heston call value using characteristic functions.
%   y = call_heston_cf(s0, v0, vbar, a, vvol, r, rho, t, k)

%   Inputs:
%   s0: stock price
%   v0: initial volatility (v0^2 initial variance)
%   vbar: long-term variance mean
%   a: variance mean-reversion speed
%   vvol: volatility of the variance process
%   r: risk-free rate
%   rho: correlation between the Weiner processes of the stock price and its variance
%   t: time to maturity
%   k: option strike
%   chfun_heston: Heston characteristic function

%  1st step: calculate pi1 and pi2
%  Inner integral 1
int1 = @(w, s0, v0, vbar, a, vvol, r, rho, t, k) real(exp(-i.*w*log(k)).*chfun_heston(s0, v0, vbar, a, vvol, r, rho, t, w-i)./(i*w.*chfun_heston(s0, v0, vbar, a, vvol, r, rho, t, -i)));

% inner integral1
int1 = integral(@(w)int1(w,s0, v0, vbar, a, vvol, r, rho, t, k),0,100);

% numerical integration
pi1 = int1/pi+0.5; %final pi1

%  Inner integral 2:
int2 = @(w, s0, v0, vbar, a, vvol, r, rho, t, k) real(exp(-i.*w*log(k)).*chfun_heston(s0, v0, vbar, a, vvol, r, rho, t, w)./(i*w));

int2 = integral(@(w)int2(w,s0, v0, vbar, a, vvol, r, rho, t, k),0,100);
int2 = real(int2);

pi2 = int2/pi+0.5; % final pi2

%  2rd step: calculate call value
y = s0*pi1-exp(-r*t)*k*pi2;

end


You typically assume that the jump component and the diffusion (sto vol) component are independent. This means that the characteristic function of the log stock price $$\ln(S_t)$$ is given by \begin{align*} \varphi_{\ln(S_t)}(u) = \varphi_{StoVol}(u) \cdot \varphi_{Jumps}(u). \end{align*}
You can simply look up the characteristic functions of the Merton and the Kou model. For example, \begin{align*} \varphi_t^\text{Kou}(u) &= \exp\left( \ln\left(S_0e^{(r-q+\omega)t}\right)iu\right) \cdot \varphi_{\sigma W_t}(u) \cdot \varphi_{N_t}(u) \\ &= \exp\left( \ln\left(S_0e^{(r-q+\omega)t}\right)iu-\frac{1}{2}\sigma^2u^2t+t\lambda\left(\frac{p\zeta}{\zeta-iu}+\frac{q\eta}{\eta+iu}-1\right)\right), \end{align*} where \begin{align*} \omega &= -\frac{1}{2}\sigma^2-\lambda\left(\frac{p\zeta}{\zeta-1}+\frac{q\eta}{\eta+1}-1\right). \end{align*}
The formula for integrating the characteristic function $$\Pi_1$$ and $$\Pi_2$$ remain unchanged. And you can then get the option price in a 'Black-Scholes-style'. Fourier methods can be sped up using the Carr Madan (1999) approach, by controlling for auxiliary models and making use of the cosine representation of the density function.
• @AttilioMeden Precisely as your current code. Just where you use chfun_heston'', you need to implement a different characteristic function. This will obviously depend on a few more model parameters. But the calculation of $\Pi_1$, $\Pi_2$ and the final option price remains the same. – KeSchn Nov 12 '19 at 21:40