I have simulated some data according to a Vasicek process and I am then trying to apply ordinary least squares (OLS) regression analysis to see how accurate the estimated model parameters are from the estimation. However, my estimated model parameters do not even seem close to what they should be. I have tried letting T get very large, and I have even simulated 100,000 paths and estimated the model for each path, and then averaged the estimates. Yet the results still look awry. So either I am simulating my data wrong, or the regression I am doing is somehow wrong. I have attached my python code below, so if anybody can take a look and let me know why my model parameters are off I would appreciate it very much.
The form of the model I am using is:
$$ dr_t = (a - b r_t) dt + \sigma dW.$$
To estimate my model parameters I am fitting a regression on the discrete data that gives me the following model:
$$ r_{t+1} - r_{t} = (a - b r_t) \Delta t + \sigma Z_t$$
where $ Z_t \sim N(0,\Delta t)$ (i.e. mean 0, variance $\Delta t$). From this setup I know that my estimated intercept will be $ a\Delta t$ and my estimated slope will be $-1 * b \Delta t$. Therefore, in my fit_vasicek_model function, I make the necessary transformation to recover "a" and "b".
Here is my code ( I have omitted some portions that are irrelevant to this question )
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import random
import statsmodels.formula.api as smf
import statsmodels.api as sm
import statsmodels.tools as sm_tools
class InterestRateModels:
def __init__(self, nsim, T, Vas = None):
self.N = nsim
self.T = T
self.dt = 1/T
self.Vas = Vas
def sim_vasicek(self, r0, a, b, sigma):
"""
Inputs:
-------
N = number of paths to simulate
T = number of time steps to take
r0 = initial (current) short-term rate
a, b, sigma = parameters of the model
The Model:
----------
dr_t = (a − b*r_t)*dt + σ * dBt
r_t = short rate
a = long-run average rate, a in Real Numbers
b = speed of reversion, b > 0
sigma = volatility of short rate
Bt = Gaussian (0 mean, 1 sigma) process
Returns:
--------
Creates an object attribute that is a matrix of simulated short rate paths.
"""
dt = 1.0/self.T
X = np.random.normal(0,np.sqrt(dt),size=(self.N,self.T))
self.Vas = np.zeros((self.N,self.T)); self.Vas[:,0] = r0;
for i in range(self.N):
for j in range(1,self.T):
self.Vas[i,j] = self.Vas[i,j-1] + (a-b*self.Vas[i,j-1]) * dt + sigma * X[i,j]
return
def plot_vasicek(self):
for i in range(self.N):
plt.plot(self.Vas[i,:])
return
def fit_vasicek_model(self, rate_path):
"""
Inputs:
-------
rate_path is a single interest rate path.
Returns:
--------
Parameters of the model (not the regression paraneters, as they get transformed)
Intercept = a * dt
Slope = -1 * b * dt (The slope needs to be multiplied by negative 1)
"""
n = len(rate_path)
X = rate_path[:n-1]
y = np.diff(rate_path)
X = sm_tools.tools.add_constant(X)
model = sm.OLS(y,X)
results = model.fit()
a, b1 = results.params / self.dt
b = -1 * b1
return a, b
So now when I create an instance and simulate a bunch of paths, then fit the model I parameters for a and b that are not even close to their true value. Here is my guess as to why: There needs to be a long enough of a period (i.e, T is large enough), relative to the reversion parameter, b. Meaning if b is large, T can be small. If b is small, T must be large.
# nsim=1000, T=30, a=1, b=4 #
rates_sim = InterestRateModels(nsim=1000, T=30)
rates_sim.sim_vasicek(.02, a=1, b=4, sigma=.3)
rates_sim.plot_vasicek()
b_lst = []
a_lst = []
for i in range(rates_sim.N):
rate_path = rates_sim.Vas[i,:]
a, b = rates_sim.fit_vasicek_model(rate_path)
a_lst.append(a)
b_lst.append(b)
np.mean(b_lst)
Out[273]: 7.0396842314975165
Notice that the mean of b is 7.03, not 4. So now when I let T=1000 this is what I see.
np.mean(b_lst)
Out[275]: 7.538313498938557
This is still a bit off. I did look at the plot of the rate curves and they don't totally look stationary. So letting b get really big (i.e. b=30), even for a smaller T=30 I get the following result
np.mean(b_lst)
Out[277]: 31.078385891056488
This looks like what I would expect. So, what is the general rule for reasonable convergence? I have read somewhere that we must have $\mid 1 - b \Delta t \mid <= 1$, though when $\Delta t = 1/30$ and $b=4$, this still holds true, yet the results look awful. Thanks for reading this, and I appreciate any help you filks might give.
