Simulating Heston Spot and Variance Paths in QuantLib

# Simulation of Heston process Quantlib-Python


# Simulation of Heston process Quantlib-Python












I am wondering weather there exists some method such that one can simulate sample paths for the Heston model in Quantlib-Python. I am currently working on a project that require simulations with the Quadratic Exponental scheme with martingale correction, as calibrated parameters severely violates the Feller condition. I have tried to implement the QE scheme by coding it my self, however, with limited success. I have also tried installation of `PyQL`, but have not been able to make it work.

I have seen that one can utilize the`ql.GaussianMultiPathGenerator` in order to simulate other multi dimensional scenarios. Is it possible one can utilize this method to simulate a CIR and a OU process and correlate them somehow? Is there some other way?

## Answer by StackG (score 6, accepted)

https://quant.stackexchange.com/a/55926

The snippets below will generate spot and vol paths from QuantLib's `HestonProcess`, and generate the plots shown.

Notice that in the vol histogram, we see a peak appearing in the 0 bucket - due to Feller not being well satisfied, we're seeing many vols landing in 0 and staying for a long amount of time

Snippet to generate the paths:

```
import QuantLib as ql
import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# Utility function to pull out spot and vol paths as Pandas dataframes
def generate_multi_paths_df(sequence, num_paths):
    spot_paths = []
    vol_paths = []

    for i in range(num_paths):
        sample_path = seq.next()
        values = sample_path.value()

        spot, vol = values

        spot_paths.append([x for x in spot])
        vol_paths.append([x for x in vol])

    df_spot = pd.DataFrame(spot_paths, columns=[spot.time(x) for x in range(len(spot))])
    df_vol = pd.DataFrame(vol_paths, columns=[spot.time(x) for x in range(len(spot))])

    return df_spot, df_vol

today = ql.Date(1, 7, 2020)
v0 = 0.01; kappa = 1.0; theta = 0.04; rho = -0.3; sigma = 0.4; spot = 100; rate = 0.0

# Set up the flat risk-free curves
riskFreeCurve = ql.FlatForward(today, rate, ql.Actual365Fixed())
flat_ts = ql.YieldTermStructureHandle(riskFreeCurve)
dividend_ts = ql.YieldTermStructureHandle(riskFreeCurve)

heston_process = ql.HestonProcess(flat_ts, dividend_ts, ql.QuoteHandle(ql.SimpleQuote(spot)), v0, kappa, theta, sigma, rho)

timestep = 8
length = 2
times = ql.TimeGrid(length, timestep)
dimension = heston_process.factors()

rng = ql.GaussianRandomSequenceGenerator(ql.UniformRandomSequenceGenerator(dimension * timestep, ql.UniformRandomGenerator()))
seq = ql.GaussianMultiPathGenerator(heston_process, list(times), rng, False)

df_spot, df_vol = generate_multi_paths_df(seq, 10000)
df_spot.head()
```

Snippet to generate the plots:

```
# Plot the first ten paths for spot and vol, and the distribution of the final path step across all paths
plt.figure(figsize=(20, 10))

plt.subplot(2, 2, 1)
plt.plot(df_spot.iloc[0:10].transpose())
plt.title("Sample Spot Paths")

plt.subplot(2, 2, 2)
plt.hist(df_spot[2.0], bins=np.linspace(0, 250, 51))
plt.title("Spot, t=2Y")

plt.subplot(2, 2, 3)
plt.plot(np.sqrt(df_vol.iloc[0:10]).transpose())
plt.title("Sample Vol Paths")

plt.subplot(2, 2, 4)
plt.hist(np.sqrt(df_vol[2.0]), bins=np.linspace(0, 0.8, 17))
plt.title("Instantaneous Vol, t=2Y")
```