# Simulation of Stochastic Volatility with Correlated Jumps (SVCJ) price paths

I am trying to simulate price paths for Monte Carlo option pricing of the Stochastic Volatility with correlated jumps model as presented in Dufffie et al.(2000), Eraker et. al. (2003) and Eraker (2004).

In short: I am looking for a pricing of the underlying formula under $$\mathbb{Q}$$.

I have seen that this post includes the function of the analytical solution of the call option. I am looking for a solution for the underlying.

Under the risk neutral measure $$\mathbb{Q}$$, the model is

\begin{aligned} \frac { d S _ { t } } { S _ { t } } & = \left( r - \mu ^ { * } \right) d t + \sqrt { V _ { t } } d W ( Q ) _ { t } ^ { S } + d J ( Q ) _ { t } ^ { S } \\ d V _ { t } & = \left( \kappa \left( \theta - V _ { t } \right) + \eta _ { V } V _ { t } \right) d t + \sigma _ { V } \sqrt { V _ { t } } d W ( Q ) _ { t } ^ { V } + d J ( Q ) _ { t } ^ { V } \end{aligned}

Quote for Eraker (2004):

"

$$\mu ^ { * }$$ is the jump compensator term, and where $$dW(Q)_t$$ is a standard Brownian motion under $$Q$$ defined by $$d W ( Q ) _ { t } ^ { i > } = \eta _ { i } d t + d W _ { t } ^ { i }$$ for $$i = V, J$$. Notice that the volatility is reverting at the rate $$\kappa ^ { Q } : = > \kappa - \eta ^ { V }$$ where $$\eta _ { V }$$ is a risk premium parameter associated with shocks to the volatility process. A similar risk premium, $$η_J$$, is assumed to be associated with jumps. This gives the price jump distribution under $$Q$$ as $$Z _ { t } ^ { S } | > Z _ { t } ^ { V } \sim N \left( \mu _ { y } ^ { Q } + \rho _ { J } Z _{ t } ^ { V } , \sigma _ { S } ^ { 2 } \right).$$ ... $$\eta _ { J }= \mu _ { y } ^ { Q } - \mu _ { y }$$"

Brodie and Kaya (2006) (6.2) provide a simulation scheme for the SVCJ model. Johannes and Pohlsen provide an Euler discretization scheme for the model. In an application, a Euler discretization is provided under $$\mathbb{P}$$ \begin{aligned} Y _ { t } & = \mu + \sqrt { V _ { t - 1 } } \varepsilon _ { t } ^ { y } + Z _ { t } ^ { y } J _ { t } \\ V _ { t } & = \alpha + \beta V _ { t - 1 } + \sigma _ { V } \sqrt { V _ { t - 1 } } \varepsilon _ { t } ^ { v } + Z _ { t } ^ { v } J _ { t } \end{aligned}.

where $$Y_{t+1} = log(S_{t+1}/S_{t}$$ ) is the log return, $$\alpha = \kappa \theta$$, $$\beta = 1 − \kappa$$ and $$\varepsilon _ { t } ^ { y } , \varepsilon _ { t } ^ { v } \sim N(0,1)$$ variables with correlation $$\rho$$. $$J_t$$ is a Bernoulli random variable with $$P(J_t = 1) = \lambda$$

I am trying to simulate paths of the underlying under $$\mathbb{Q}$$ for Monte Carlo option pricing. I have also tried to follow Broadie's steps, yet results were not satisfying and looked wrong.

I have simulated paths under $$\mathbb{P}$$ using the discretization about, yet I am not sure whether that is correct. #Poisson process simulation
Z = -log(runif(n))/lambda
t = seq(0, sum(Z))
X = t * 0
for (i in 1:n) {
Zc = cumsum(Z)
X[t >= Zc[i]] = i
}

svjc_price_sim    =
function(n,kappa,theta,sigmav,rho,mu_x,var_x,lambda,mu_vol,rho_j,mu,PP){

#Storage vectors stochastic volatility component
#Storage vector
Y         = c(rep(0,n))                  #log returns

#Storage vector jump component
V         = c(rep(0,n))                  #volatility
V      = mu_vol                       #set volatility mean for initialization

#Simulated volatility

for(i in 2:n){

V[i]    = kappa * theta + (1 - kappa) * V[i-1] + sigmav * sqrt(V[i - 1]) * rnorm(1) + rexp(n = 1, rate = mu_vol) * PP
a)
}

#Simulated returns

for(i in 2:n){
Y[i]    = sqrt(V[i - 1]) * rnorm(1) + rnorm(n = 1, mean = mu_x + rho_j * rexp(n = 1, rate = 1/mu_vol), sd = var_x)* PP

}

#delete first row

#Price conversion ln(St/St-1) = Rt
price    = c(rep(0,n))
price = 100           #base value

for(i in 2:n){
price[i] = (Y[i]/sqrt(252) + 1)*price[i-1]
}

return(price)

}

underlyingP = svjc_price_sim(n,kappa,theta,sigmav,rho,mu_x,var_x,lambda,mu_vol,rho_j,mu,X)

plot(underlyingP, type = "l", main = "Price path of the underlying St")


I tried the simulation under $$\mathbb{Q}$$, but the plots looked wrong.

I would be very thankful about any code/pseudocode/Euler discretization regarding the simulation of SVCJ under $$\mathbb{Q}$$. I assume that the "standard approach" is to be proceeded, that is, apply Itos lemma to $$log (S_t)$$ and then integrate-however, the model structure is quite complex.