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 theql.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?

  • $\begingroup$ I don't understand this sentence: "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". Why aren't you using Feller condition in your calibration? $\endgroup$
    – user39119
    Jan 24, 2020 at 10:44
  • $\begingroup$ @UBM, The reason i do not enforce the Feller condition in the calibration of the model is that our we utilize the "deep calibration" algorithm developed by Horvath et al (2019), and they only limit the "training" parameter combinations to to adhere to the feller condition, not the actual calibration task of minimizing the difference between the neural network (trained) pricing function and market quotes. But surely, enforcing this criterion in terms of the calibration task would surely limit the method to find a local minimum of the loss function (squared distance)? $\endgroup$ Jan 24, 2020 at 11:02
  • $\begingroup$ @UBM, i also tried not enforcing the feller condition in the training data generation, then the calibrated parameters was more inline with the Feller condition, however still violating the feller condition. However, it is perhaps a good idea to force the condition in both training and (explicitly) calibration. Just not completely sure of its theoretical justification in terms of limiting the "search space" if you will. $\endgroup$ Jan 24, 2020 at 11:12

1 Answer 1


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

Example Heston Spot and Vol Paths

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)

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.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.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")

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.