Calibration of Theta, A(t) and B(t) of Hull White 1Factor model

Question

I’m simulating interest rates via the HullWhite One factor model. To simulate the short rate I’m using the code from the Quantlib Python Cookbook, chapter 15 and beyond (By Goutham Balaraman and Luigi Ballabio). The code can be found here:

The code generates the short rate. Unfortunately the code ends with having only the simulated short rate while my goal is to create a whole new term structure based on the simulated curve.

To create a term structure I need to calibrate theta based on the initial term structure. Besides that I also need to calculate A(t) and B(t) such that I’m able to create the term structure.

Questions: Does somebody have an algorithm to calculate theta, A(t) and B(t)?

Answer 1

You can check out here a blog post on simulating the yield term structure for the HullWhite model.

The basic idea is that once you have the paths for the short rate, you can simply integrate (approximately) the short rate throughout each path to obtain the discount factors.

The average of the simulations should match the initial term structure.

Here is an example based on that excellent blog post, with an improvement on execution time by vectorizing the integration.

import QuantLib as ql
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
import copy 

nPaths = 500
years = 30

startDate = ql.Date(3, 12, 2018)
endDate = startDate + ql.Period(years, ql.Years)
tenor = ql.Period(1, ql.Days)
schedule = ql.MakeSchedule(startDate, endDate, tenor)

dates = [dt for dt in schedule]

times = [ql.Actual360().yearFraction(startDate, dt) for dt in dates]

curve = ql.YieldTermStructureHandle(ql.FlatForward(startDate, 0.04875825, ql.Actual365Fixed()))
reversionSpeed = 0.01
rateVolatility = 0.001
process = ql.HullWhiteProcess(curve, reversionSpeed, rateVolatility)

periods = (endDate - startDate) + 1
sequenceGenerator = ql.UniformRandomSequenceGenerator(periods, ql.UniformRandomGenerator())
gaussianSequenceGenerator = ql.GaussianRandomSequenceGenerator(sequenceGenerator)
pathGenerator = ql.GaussianPathGenerator(process, years, periods, gaussianSequenceGenerator, False)

paths = np.zeros(shape = (nPaths, periods))
for i in range(nPaths):
    path = pathGenerator.next().value()
    paths[i, :] = np.array([path[j] for j in range(periods)])
    
dfs = np.zeros(shape=paths.shape)
dt = years / len(schedule)
integral = 0
for j in range(periods):
    if j == 0:
        dfs[:, 0] = 1
        integral = copy.deepcopy(paths[:, 0])
    else:
        integral += paths[:, j]        
        dfs[:, j] += np.exp(-integral * dt)
    
simulatedCurve = ql.DiscountCurve(dates, dfs.mean(0), ql.Actual365Fixed(), ql.NullCalendar())
simulatedDfs = np.array([simulatedCurve.discount(dt) for dt in dates])
dfs_curve = np.array([curve.discount(dt) for dt in dates])

fig, ax = plt.subplots(1,2, figsize=(10,3))
plt.tight_layout()

ax[0].set_title('Discount Factors')
ax[0].plot(times, simulatedDfs, linestyle = 'dashed', label = 'simulated curve')
ax[0].plot(times, dfs_curve, linestyle = 'solid', label = 'original curve')
ax[0].legend()

ax[1].set_title('Paths')
ax[1].plot(paths.T, linewidth=0.5)

The output would be: