"""ir.py: Script that compares various Interest Rate term structure models. A common model used in the financial industry for modelling the short rate (think overnight rate, but actually an infinitesimally short amount of time) is the Vasicek model. # Although it is unlikely to perfectly fit the yield curve, it has some nice properties that make it a good model to work with. The dynamics of the Vasicek model are describe below. In this model, the parameters are constants, and the random motion is generated by the Q measure Brownian motion . An important property of the Vasicek model is that the interest rate is mean reverting to , and the tendency to revert is controlled by . Also, this process is a diffusion process, hence Markovian, which will lead to some nice closed form formulas. Finally, the future value of the interest rate is normally distributed with the distribution . """ __author__ = "Alexander Baker" __copyright__ = "Copyright 2012, The Hatch Project" __credits__ = ["Alexander Baker"] __license__ = "GPL" __version__ = "1.0.1" __maintainer__ = "Alexander Baker" __email__ = "baker.alexander@gmail.com" __status__ = "Production" import numpy as np from numpy.random import standard_normal from numpy import array, zeros, sqrt, shape from pylab import plot, title, xlabel, ylabel, show r0 = 0.05 # current UK funding rate theta = 0.10 # 1 % long term interest rate k = 0.3 beta = 0.03 ## simulate short rate paths n = 10000 # MC simulation trials T = 24. # total time m = 100. # subintervals dt = T/m # difference in time each subinterval r = np.zeros(shape=(n, m), dtype=float) # matrix to hold short rate paths for j in np.arange(0,n): # number of MC simulations r[j,0] = r0 for i in np.arange(1,m): #trials per simulation r[j,i] = r[j,i-1] + k*(theta-r[j,i-1])*dt + beta*sqrt(dt)*standard_normal(); #print j,i, r[j,i] plot(np.arange(0, T, dt), r[j]) ## plot paths t = np.arange(0, T, dt) rT_expected = theta + (r0-theta)*pow(np.e,-k*t) rT_stdev = sqrt( pow(beta,2)/(2*k)*(1-pow(np.e,-2*k*t))) print 'expected', rT_expected, 'std', rT_stdev plot(t, rT_expected, '-+r') plot(t, rT_expected+2*rT_stdev, '-b') plot(t, rT_expected-2*rT_stdev, '-b') plot(t, rT_expected+4*rT_stdev, '-g') plot(t, rT_expected-4*rT_stdev, '-g') print shape(t), shape(r) title('Simulations %d Steps %d r0 %.2f alpha %.2f beta %.2f sigma %.2f' % (int(n), int(m), r0, k, theta, beta)) xlabel('steps') ylabel('short term interest rate') show() #rT.stdev <- sqrt( beta^2/(2*k)*(1-exp(-2*k*t))) #matplot(t, r[,1:10], type="l", lty=1, main="Short Rate Paths", ylab="rt") #abline(h=theta, col="red", lty=2) #lines(t, rT.expected, lty=2) #lines(t, rT.expected + 2*rT.stdev, lty=2) #lines(t, rT.expected - 2*rT.stdev, lty=2) #points(0,r0)