%% Simulate Future Term Structures Using a CIR Model
% Using the short rate, simulate the rate dynamics and term structures in the future using a CIR model. The CIR model is expressed as
% dr(t)=\alpha(b-r(t))dt+\sigma\sqrt{r(t)}dW(t)
% The exponential affine form of the bond price is
% B(t,T)=e^{-A(t,T)r(t)+C(t,T)}
% where
% A(t,T) = \frac{2(e^{\gamma (T-t)}-1)}{(\gamma+\alpha)(e^{\gamma(T-t)}-1)+2\gamma}
% B(t,T) = \frac{2\alpha b}{\sigma^2} \log\left(\frac{2\gamma e^{(\alpha + \gamma) (T-t)/2}}{(\gamma+\alpha)(e^{\gamma(T-t)}-1)+2\gamma}\right)
% and
% \gamma=\sqrt{\alpha^2+2\sigma^2}
% Define the parameters for the cir object.

alpha = .1;
b = .05;
sigma = .05;
r0 = .04;

% Define the function for bond prices.

gamma = sqrt(alpha^2 + 2*sigma^2);
A_func = @(t, T) ...
    2*(exp(gamma*(T-t))-1)/((alpha+gamma)*(exp(gamma*(T-t))-1)+2*gamma);
C_func = @(t, T) ...
    (2*alpha*b/sigma^2)*log(2*gamma*exp((alpha+gamma)*(T-t)/2)/((alpha+gamma)*(exp(gamma*(T-t))-1)+2*gamma));
P_func = @(t,T,r_t) exp(-A_func(t,T)*r_t+C_func(t,T));

% Create a cir object.

obj = cir(alpha,b,sigma,'StartState',r0)

% Define the simulation parameters.

nTrials = 100;
nPeriods = 5;   % Simulate future short over the next five years
nSteps = 12;    % Set intermediate steps to improve the accuracy

% Simulate the short rates. The returning path is a (NPeriods + 1)-by-NVars-by-NTrials three-dimensional time-series array. For this example, the size of the output is 6-by-1-by-100.

rng('default');    % Reproduce the same result
rPaths = simByTransition(obj,nPeriods,'nTrials',nTrials,'nSteps',nSteps);
size(rPaths)
rPathsExp = mean(rPaths,3);

% Determine the term structure over the next 30 years.

maturity = 30;
T = 1:maturity;
futuresTimes = 1:nPeriods+1;

% Preallocate simTermStruc
simTermStructure = zeros(nPeriods+1,30);
for i = futuresTimes
    for t = T
        bondPrice = P_func(i,i+t,rPathsExp(i));
        simTermStructure(i,t) = -log(bondPrice)/t;
    end
end
plot(simTermStructure')
legend('Current','1-year','2-year','3-year','4-year','5-year')
title('Projected Term Structure for Next 5 Years')
ylabel('Long Rate Maturity R(t,T)')
xlabel('Time')

% Copyright 2018 The MathWorks, Inc.