A Simple Bootstrapping Algrithm For Calulating Projected Bond Yield Rates

This some mathematica code I used to implement a simple bootstrapping algorithm in an effort to compute projected bond yield rates (along with some other information related to the resulting time series) for the next 5 years. The data ragarding bond prices, coupon payments, and advertized yield rates all came from  http://www.pfin.ca/canadianfixedincome/Default.aspx, and was collected over the span of 10 consecutive business days, starting on January 13, 2016.

A note on our basic assumptions, to compute the desired yield/forward rates, we assumed that the yeild rate is constant for the first year, where the “first year” means the date the data was collected to that date in the subsequent year (e.g. January 13, 2016 to January 13, 2017), as a matter of convention we take the first year rate to be the projected rate on the date one year after the date the data was collected, (e.g. if the data was collected on January 13, 2016, then the one year rate corresponds to January 13, 2017). After the first year, we use simple linear interpolation to compute the yield rates for the following years, again taking the (n^th) year rate to always fall 365 days (or 366 for leap years) after what was taken to be the (((n-1))^th) year rate.

s = Import["Documents/1850-1.xlsx"] (*  Importing the data from excel  *)

P ={{"","Canada",4.`,"2016-Jun-01",101.37`,0.35`,101.35`,0.34`,101.35`,0.3`,101.34`,0.32`,101.33`,0.31`,101.29`,0.39`,101.24`,0.43`,101.23`,0.44`,101.22`,0.44`,101.21`,0.44`},{"","Canada",4.`,"2017-Jun-01",105.04`,0.32`,105.`,0.34`,105.04`,0.3`,104.99`,0.33`,105.02`,0.3`,104.88`,0.39`,104.78`,0.44`,104.77`,0.45`,104.8`,0.42`,104.77`,0.43`},{"","Canada",4.25`,"2018-Jun-01",109.3`,0.32`,109.24`,0.33`,109.31`,0.3`,109.22`,0.33`,109.26`,0.31`,109.02`,0.4`,108.86`,0.46`,108.82`,0.47`,108.89`,0.43`,108.85`,0.45`},{"","Canada",3.75`,"2019-Jun-01",111.3`,0.37`,111.27`,0.38`,111.4`,0.34`,111.28`,0.37`,111.3`,0.36`,111.01`,0.43`,110.8`,0.49`,110.75`,0.51`,110.87`,0.47`,110.82`,0.48`},{"","Canada",3.5`,"2020-Jun-01",113.17`,0.45`,113.14`,0.46`,113.41`,0.4`,113.27`,0.43`,113.25`,0.43`,112.99`,0.48`,112.66`,0.55`,112.42`,0.6`,112.65`,0.55`,112.62`,0.55`}}

(*  the data can be cleaned up a little ...  *)
i= 1;
While[i ≤ 5, 
P[[i]] = Delete[P[[i]],{{1},{2},{6},{8},{10},{12},{14},{16},{18},{20},{22},{24}}];
i++
]

P

(*
Here we point out that P[[i]]] is a list of data correponding to a single Government of Canada bond 
P[[i,1]] is the annual coupon payment for that bond
P[[i,2]] is the effective maturity date of the bond
P[[i, 3 ;; 12]] are the prices of the bond recorded on the 10 subsequent days
*)

(*  this is the function we will use to compute all of the desired yield/forward rates  *)
Rates[P_]:= Module[{i, M, R, r1, 
midr2, m2, b2, y2, r2,
midr3, m3, b3, y3, r3, 
midr4, m4, b4, y4, r4,
midr5, m5, b5, y5, r5},
i= 0; (*  initializing our counter  *)
M = {}; (*  This will be a list of lists, where each nested list contains the forward rates for the fives subsequent years following 2016  *)
R = {};  (*  This will be a list of lists, where each nested list contains the yield rates for the fives subsequent years following 2016  *)
While[i ≤ 9, 
(*
We start off by calculating the first year yield rate. Here we assume that the rate is constant over the first year - so to linearly interpolate for this year is trivial.
*)
(*  r1 is the (constant) first year rate  *)
r1 = ((-365)/(DateDifference[{2016, 1, 13 + i},{2016, 6, 1}]))Log[(P[[1,3 + i]]+ P[[1,1]](((DateDifference[{2015, 12, 1},{2016, 1, 13 + i}])/365)) )/(100 + ((P[[1,1]])/2))];



(*
Next we compute the second year rate. Since from the second year onward we are linearly interpolating between years. We first use our data to calculate the yield rate on June (6^th), 2017 and then linearly interpolate to get the yield rate for January (((13+i))^th), 2018.
*)
midr2 = ((-365)/(DateDifference[{2016, 1, 13 + i},{2017, 6, 1}]))Log[(( P[[2,3 + i]]+P[[2,1]](((DateDifference[{2015, 12, 1},{2016, 1, 13 + i}])/365))) - ((P[[2,1]])/2)((E^(-r1 DateDifference[{2016, 1, 13 + i},{2016, 6, 1}]/365))+ (E^(-r1 DateDifference[{2016, 1, 13 + i},{2016, 12, 1}]/365))))/(100 + ((P[[2,1]])/2))];
(*
the following two computations calculate the slope and y-intercept of the line we will use in the interpolation of the second year rate

namely, we want    y2 = m2*x + b2     knowing the rate for the dates, January (((13+i))^th), 2017 and June (6^th), 2017
*)
m2 = ((midr2 - r1)/(DateDifference[{2017, 1, 13 + i}, {2017, 6,1}]));
b2 = r1 - m2*DateDifference[{2016, 1, 13},{2017, 1, 13}];
(*
Next is the equation for the interpolation line
*)
y2[y_, m_, d_]:= m2*DateDifference[{2016, 1, 13+i}, {y, m, d}] + b2;
(*
and finally the second year rate!
*)
r2 = y2[2018, 1, 13+i];




(*
From here on we are more or less repeating a structurally similar computation to compute the rates of subsequent years using the rates calculated for the previous years. 
*)

(*  THIRD YEAR RATE  *)
midr3 = ((-365)/(DateDifference[{2016, 1, 13 + i},{2018, 6, 1}]))Log[((P[[3,3 + i]] +P[[3,1]](((DateDifference[{2015, 12, 1},{2016, 1, 13 + i}])/365))) - ((P[[3,1]])/2)((E^(-r1 DateDifference[{2016, 1, 13 + i},{2016, 6, 1}]/365))+(E^(-r1 DateDifference[{2016, 1, 13 + i}, {2016, 12, 1}]/365))+ (E^(-y2[2017, 6, 1] DateDifference[{2016, 1, 13 + i}, {2017, 6, 1}]/365))+ (E^(-y2[2017, 12, 1] DateDifference[{2016, 1, 13 + i}, {2017, 12, 1}]/365))))/(100 + ((P[[3,1]])/2))];
m3 = ((midr3 - r2)/(DateDifference[{2018, 1, 13 + i}, {2018, 6,1}]));
b3 = r2 - m3*DateDifference[{2016, 1, 13},{2018, 1, 13}];
y3[y_, m_, d_]:= m3*DateDifference[{2016, 1, 13 + i}, {y, m, d}] + b3;
r3 = y3[2019, 1, 13+i];


(*  FOURTH YEAR RATE  *)

midr4 = ((-365)/(DateDifference[{2016, 1, 13 + i}, {2019, 6, 1}]))Log[((P[[4,3 + i]] +P[[4,1]](((DateDifference[{2015, 12, 1},{2016, 1, 13 + i}])/365))) - ((P[[4,1]])/2)((E^(-r1 DateDifference[{2016, 1, 13 + i}, {2016, 6, 1}]/365))+(E^(-r1 DateDifference[{2016, 1, 13 + i}, {2016, 12, 1}]/365))+ (E^(-y2[2017, 6, 1] DateDifference[{2016, 1, 13 + i}, {2017, 6, 1}]/365))+ (E^(-y2[2017, 12, 1] DateDifference[{2016, 1, 13 + i}, {2017, 12, 1}]/365))+(E^(-y3[2018, 6, 1] DateDifference[{2016, 1, 13 + i}, {2018, 6, 1}]/365)) + (E^(-y3[2018, 12, 1] DateDifference[{2016, 1, 13 + i}, {2018, 12, 1}]/365))))/(100 + ((P[[4,1]])/2))];
m4 = ((midr4 - r3)/(DateDifference[{2019, 1, 13+i}, {2019, 6,1}]));
b4 = r3 - m4*DateDifference[{2016, 1, 13},{2019, 1, 13}];
y4[y_, m_, d_]:= m4*DateDifference[{2016, 1, 13+i}, {y, m, d}] + b4;
r4 = y4[2020, 1, 13+i];


(*  FIFTH YEAR RATE  *)


midr5 = ((-365)/(DateDifference[{2016, 1, 13 + i}, {2020, 6, 1}]))Log[((P[[5,3 + i]] +P[[5,1]](((DateDifference[{2015, 12, 1},{2016, 1, 13 + i}])/365))) - ((P[[5,1]])/2)((E^(-r1 DateDifference[{2016, 1, 13 + i}, {2016, 6, 1}]/365))+(E^(-r1 DateDifference[{2016, 1, 13 + i}, {2016, 12, 1}]/365))+ (E^(-y2[2017, 6, 1] DateDifference[{2016, 1, 13 + i}, {2017, 6, 1}]/365))+ (E^(-y2[2017, 12, 1] DateDifference[{2016, 1, 13 + i}, {2017, 12, 1}]/365))+(E^(-y3[2018, 6, 1] DateDifference[{2016, 1, 13 + i}, {2018, 6, 1}]/365)) + (E^(-y3[2018, 12, 1] DateDifference[{2016, 1, 13 + i}, {2018, 12, 1}]/365))+ (E^(-y4[2019, 6, 1] DateDifference[{2016, 1, 13 + i}, {2019, 6, 1}]/365))+(E^(-y4[2019, 12, 1] DateDifference[{2016, 1, 13 + i}, {2019, 12, 1}]/365))))/(100 + ((P[[5,1]])/2))];
m5 = ((midr5 - r4)/(DateDifference[{2020, 1, 13+i}, {2020, 6,1}]));
b5 = r4 - m5*DateDifference[{2016, 1, 13},{2020, 1, 13}];
y5[y_, m_, d_]:= m5*DateDifference[{2016, 1, 13+i}, {y, m, d}] + b5;
r5 = y5[2021, 1, 13+i];


(*
And so before we move on to computing the rates for the next days worth of data, we append the data to the appropriate lists.
*)
AppendTo[M, {((r2 - r1)/2), ((r3 - r1)/2), ((r4 - r2)/2), ((r5 - r3)/2), m5*DateDifference[{2020, 1, 13},{2021, 1, 13}]}] (*  a list of forward rates  *);
AppendTo[R, {r1, r2, r3, r4, r5}] (*  a list of yield rates  *);
i++;
];
{R, M}
]

info = Rates[P] (*  storing our computed information  *)

(*  We now calculate a time series of daily log returns - in this cell it is for the yield rates, and two cells below it is for the forward rates  *)
LRr = {};
i = 1;
While[i ≤ 9,
AppendTo[LRr, Log[Abs[((info[[1]][[i + 1]])/(info[[1]][[i]]))]]];
i++
]

LRr

(*  Daily log returns for the forward rates   *)
LFr = {};
i = 1;
While[i ≤ 9,
AppendTo[LFr, Log[Abs[((info[[2]][[i + 1]])/(info[[2]][[i]]))]]];
i++
]

LFr

(*  Finally we calculate the eigenvalues and eigenvectors corresponding to the covariance matrices of the time series daily log returns for the yield and forward rates  *)
Eigensystem[Covariance[LRr]]

Eigensystem[Covariance[LFr]]