402-FIZeroCoupon.Rmd

---
editor_options:
  markdown:
    wrap: 72
output: pdf_document
---

# The Term Structure of Interest Rates

```{r, echo=FALSE, results='hide', message=FALSE}
library(tufte)
library(kableExtra)
source("formatters.R", local=knitr::knit_global())
```

```{r, init-ZeroCoupon,echo=TRUE,results='hide', message=FALSE}
  library(fBasics)
  library(empfin)
  library(YieldCurve)
  library(RQuantLib)
  data(FedYieldCurve, package='YieldCurve')
  data(ECBYieldCurve, package='YieldCurve')
  data(LiborRates, package='empfin')
```

`r newthought("The relationship")` between yield and maturity is called the "term
structure of interest rate". Since there are many ways to measure yield,
there are also many types of term structures. The term structure of
interest rate is important because, everything else being equal,
maturity is the main determinant of yield. In addition, investors and
economists believe that the shape of this curve may reflect the market's
expectation of future rates.

## Term Structure Shapes

The term structure of interest rates can assume a variety of shapes, as
illustrated by the following 2 examples.

Figure \@ref(fig:plot-ECB-1) displays a sample of yield curves, published by
the European Central Bank (ECB). These curves represent the yield to
maturity for AAA rated euro-zone government bonds. We show quarterly
data from December 2006 to June 2009. Three patterns are visible:

-   an almost flat curve (August 2007)
-   an inverted curve at the short end, and upward sloping at the long
    end (January 2008)
-   an evenly upward-sloping curve (December 2006, June 2009)

The data is extracted from the package, data set :

```{r, label=ECB-1, echo=T}
  tau <- c(3/12, 6/12, 1:30)
  d1 <- as.Date('2006-12-29')
  d2 <- as.Date('2009-07-24')
  nbObs <- dim(ECBYieldCurve)[1]
  dtObs <- seq(d1, d2, length.out=nbObs)
  indx <- seq(1, nbObs, 90)
  YC <- ECBYieldCurve[indx,]
  dtObs <- dtObs[indx]
```
 
```{r plot-ECB-1, echo=FALSE, eval=TRUE, fig.cap="Euro Sovereign Zero Yield Curves"}
   tau <- c(3/12, 6/12, 1:30)
  nbSample <- length(indx)
  mat = dtObs[1] + (tau *365) 
  plot(mat, YC[1,],type='l', ylim=c(0,5), ylab='Euro Sovereign Bond Yield', col=1, lwd=2)

  for(i in seq(2,nbSample)) {
    mat = dtObs[i] + (tau *365) 
    lines(mat, YC[i,],type='l', col=i, lwd=2)
  }

  legend("bottomright",legend=as.character(dtObs), col=seq(1,nbSample),lty=1, lwd=2)
```


Figure \@ref(fig:fig-plot-fed) presents a similar plot for the US Treasury yield,
over the period 1982-2009. This illustrates the wide range of slope that
a yield curve can assume over a long period of time. The plot uses the
data set , and is obtained as follows:

```{r plot-fed, eval=FALSE}
  tau <- c(3/12, 6/12, 1, 2, 3, 5, 7, 10)
  dtObs <- as.Date(timeSequence(from='1982-01-01', to='2012-11-30', by='month'))
  nbObs <- dim(FedYieldCurve)[1]
  indx <- seq(1, nbObs, 12)
  # monthly sample
  YC <- FedYieldCurve[indx,]
  dtObs <- dtObs[indx]
  nbSample <- length(indx)
  mat = dtObs[1] + (tau*365) 
  plot(mat, YC[1,],type='l', ylim=c(0,15), xlim=c(dtObs[1], dtObs[nbSample]+10*365),
       col=1, lwd=2, ylab='Treasury Yield', xlab='Maturity')
  for(i in seq(2,nbSample)) {
    mat = dtObs[i] + (tau *365) 
    lines(mat, YC[i,],type='l', col=i, lwd=2)
  }
```

```{r fig-plot-fed,fig=T, echo=F}
<<plot-fed>>
```

## Building a Zero-Coupon Bond Curve

In a stylized manner, the calculation of a bond zero-coupon yield curve
is straight-forward. Assume a set of bonds with cash flows $F_{i,j}$
occurring at regular intervals, where $i$ is the index of the bond and
$j$ the index of the cash flow date. We look for an interest rate
function $z(t)$ that prices the bonds exactly:

$$P_i = \sum_{j=1}^{n_i} F_{i,j} e^{-z(t_j)t_j}$$

The function $z(t)$ is called the continuously compounded zero-coupon
yield curve, or "continuous zero curve" in short.

How to perform the calculation is best explained with an example. Assume
that we observe 4 bonds, as summarized in the following table:

```{r, eval=TRUE, echo=FALSE}
  Maturity <- seq(4)
  BondYield <- c(5.5, 6.0, 7.5, 8.5)
  Coupon <- c(9,7,5,4)/100
 df <- data.frame(Maturity, Coupon, BondYield)
 knitr::kable(df, align="l") %>%
   kable_styling(full_width=FALSE)
```

Let $y_i$ be the bond yield for bond $i$ maturing in $i$ years, and
$P_i$ the corresponding present value.

```{r, zc-1, echo=FALSE}
  y <- c(5.5, 6.0, 7.5, 8.5)/100
  c <- c(9,7,5,4)/100
  P <- sapply(seq(4), function(i) SimpleBondPrice(c[i], i, yield=y[i]))
  z <- rep(0,4)
  r <- rep(0,4)
```

The zero-coupon rates $z_i$ for each maturity are computed recursively,
one maturity at a time. The one year zero-rate is obtained by converting
the bond yield into a continuously compounded rate:

$$e^{-z_1} = \frac{1}{(1+y_1)}$$

or,

$$z_1 = \ln(1+y_1)$$

```{r, zc-2, echo=F}
  z[1] <-  log(1+y[1])
```

which gives $z_1 = `r round(z[1]*100,2) `\%$.

The zero coupon rate for year $i$, with the rates up to year $(i-1)$
assumed known, is obtained by solving for $z_i$:

$$P_i = 100 \left( \sum_{j=1}^{i-1} c e^{-j z_j} + (1+c)e^{-i z_i} \right)$$

which yields:

$$z_i = -\frac{1}{i} \log \left( \frac{1}{1+c} \left( \frac{P_i}{100} - c \sum_{j=1}^{i-1} e^{-j z_j} \right) \right)$$

This recursive calculation is performed by the function :

```{r, zc-3, echo=TRUE}
  ZCRate <- function(c, z, i) { 
    # compute the ZC rate for maturity i, given the zero-coupon rates 
    # for prior maturities
    zc <- -(1/i)*log( 1/(1+c[i]) * ( P[i]/100 - c[i] *
          sum(exp(-seq(i-1)*z[1:(i-1)]))))
    zc
  }

  for(i in seq(2,4)) {
    z[i] <- ZCRate(c, z, i)
    }
```

A third term structure of interest is the par yield curve. For a given
maturity, the par yield is the coupon rate of a bond priced at par (i.e.
priced at $100\%$). Formally, given a zero-coupon curve, the par yield
$r_i$ for maturity $i$ is such that:

$$1 = r_i \sum_{j=1}^{i} e^{-j z_j} + e^{-i z_i}$$

Back to our example, the par yields are computed by the following
function, noting that the one year par rate is by definition the
corresponding bond yield:

```{r, zc-4, echo=TRUE}
  ParRate <- function(i) {
    # compute the par rate for maturity i,
    # given the zero-coupon rates for all maturities 
    r <- (1-exp(-i*z[i])) / sum(exp(-seq(i)*z[1:i]))
    r
  }

  r[1] <- y[1]
  for(i in seq(2,4)) {
    r[i] <- ParRate(i)
    }
```

Figure \@ref(fig:zc-11) displays the bond yield curve and the corresponding
zero-coupon curve and par curves.

```{r zc-11, echo=F, fig.cap="Three curves derived from the same bond prices: Yield to maturity, zero coupon yield  and par yield"}
  plot(seq(4), y, xlab='Maturity', ylab='Rate', type='l', col=1)
  lines(seq(4), z, col=2)
  lines(seq(4), r, col=3)
  legend("bottomright",legend=c('YTM', 'ZC', 'Par'), col=c(1,2,3),lty=1)
```

In practice, this simple method is not realistic for several reasons:

-   In most markets, we will not find a set a bonds with identical
    anniversary dates,

-   the method outlined above only provides values of $z(t)$ at maturity
    dates corresponding to the cash flow dates. An interpolation method
    is needed to obtain $z(t)$ at arbitrary maturities, and finally

-   it does not account for the fact that some bonds may trade at a
    premium or discount: for example, a bond with a small outstanding
    amount will often trade at a discount, while the current 10 year
    benchmark bond will trade at a premium.

If we insist on a zero-coupon curve that prices exactly every bond in
the sample, we will get a unnatural shape for the zero curve, and even
more so for the forward rates.

The alternative is to specify a functional form (sufficiently flexible)
for the zero-coupon curve, and to estimate its parameters by minimizing
the pricing error over a set of bonds.

Since the term structure of interest rate can equivalently be
represented in terms of spot rates, forward rates or discount factors,
one has three choices for specifying the functional form of the term
structure.

## Functional forms for the zero-coupon curve

A popular method is to model the instantaneous forward rate, $f(t)$.
Integrating this function gives the zero-coupon rate. @Nelson1987
introduced the following model for the instantaneous forward rate:

$$f(t) = \beta_0 + \beta_1 e^{-\frac{t}{\tau_1}}
+ \beta_2 \frac{t}{\tau_1} e^{-\frac{t}{\tau_1}}$$

Integrating $f(t)$, the corresponding zero-coupon rates are given by:

$$z(t) = \beta_0 + \beta_1 \frac{1-e^{-\frac{t}{\tau_1}}}{\frac{t}{\tau_1}} \\
  + \beta_2 \left( \frac{1-e^{-\frac{t}{\tau_1}}}{\frac{t}{\tau_1}} - e^{-\frac{t}{\tau_1}} \right)$$

@Svensson1994 extended this model with 2 additional parameters, defining
the instantaneous forward rate as:

$$f(t) = \beta_0 + \beta_1 e^{-\frac{t}{\tau_1}}
+ \beta_2 \frac{t}{\tau_1} e^{-\frac{t}{\tau_1}}
+ \beta_3 \frac{t}{\tau_2} e^{-\frac{t}{\tau_2}}$$

The following script, taken from the package, illustrates the fitting of
the Nelson-Siegel and Svensson models to US Treasury yields.

```{r, label=estimate-zc, echo=T}
  # maturities in months of FedYieldCurve data set
  tau <- c(3, 6, 12, 24, 36, 60, 84, 120)

  # parameter estimation
  Parameters.NS <- Nelson.Siegel( rate=FedYieldCurve[5,], maturity=tau)
  Parameters.S <- Svensson(rate=FedYieldCurve[5,], maturity=tau)

  # sample the fitted curves 
  tau.sim <- seq(1,120,2)
  y.NS <- NSrates(Parameters.NS,tau.sim)
  y.S <- Srates(Parameters.S, tau.sim, whichRate='Spot')
```

Figure \@ref(fig:plot-nss) shows a comparison of the two fitted spot curves.

```{r plot-nss, echo=FALSE,eval=TRUE, fig.cap="Fitted Treasury zero-coupon spot curves, with model of Svensson and Nelson-Siegel"}

plot(tau,as.double(FedYieldCurve[5,]), type="p", xlab='Maturity (months)', 
     ylab='Spot rate (%)')
lines(tau.sim,as.double(y.NS), col=2)
lines(tau.sim,as.double(y.S), col=3)
legend("bottomright",legend=c("actual","N-S", 'Svensson'), col=c(1,2,3),lty=1)
```
The corresponding instantaneous fitted forward curves are computed as follows:

```{r zc-nss-fwd, echo=FALSE, eval=TRUE, fig.cap=""}
f.S <- Srates(Parameters.S, tau.sim, whichRate='Forward')
forward.NS <- function(params, tau) {
MoverTau <- as.double(params$lambda)*tau
eMoT <- exp(-MoverTau)
beta = as.double(params[1,1:3])
beta[1]*rep(1, length(tau)) + beta[2]*eMoT + beta[3]*MoverTau*eMoT
}
f.NS <- forward.NS(Parameters.NS[1,],tau.sim)
```

```{r nss-fwd, echo=F,eval=TRUE, fig.cap="Fitted Treasury Zero-coupon forward curves, with Svensson model and Nelson-Siegel model"}
plot(tau.sim, f.NS, ylim=c(min(c(f.NS, f.S)), max(c(f.NS,f.S))), 
     xlab='Maturity', ylab='Forward rate', type='l', col=1)
lines(tau.sim, f.S, col=2)
legend("bottomright",legend=c("N-S", 'Svensson'), col=c(1,2),lty=1)
```

The Svensson or Nelson-Siegel curves can be easily decomposed into 3 elementary shapes that capture the main features of a yield curve: level, slope and curvature. Figure \@ref(fig:svensson-3) shows the decomposition of the Svensson fitted spot curve into these three components.
 
```{r svensson-3, echo=FALSE,eval=TRUE, fig.cap="Components of Fitted Treasury zero-coupon spot curves, with Nelson-Siegel model"}
# NS parameters
p.NS <- Parameters.NS
p.1 <- p.NS
p.1$beta_1 <- 0
p.1$beta_2 <- 0
p.2 <- p.NS
p.2$beta_0 <- 0
p.2$beta_2 <- 0
p.3 <- p.NS
p.3$beta_0 <- 0
p.3$beta_1 <- 0

y.1 <- NSrates(p.1, tau.sim)
y.2 <- NSrates(p.2, tau.sim)
y.3 <- NSrates(p.3, tau.sim)

plot(tau.sim, y.1, xlab='Maturity', ylab='Spot rate', type='l', col=1)
par(new=T)
plot(tau.sim, y.2, ylim=c(min(c(y.2, y.3)), max(c(y.2, y.3))), axes=F,xlab="",ylab="", col=2)
lines(tau.sim, y.3, col=3)
axis(side=4)
legend("bottomright",legend=c("level", 'slope', 'convexity'), col=c(1,2,3),lty=1)
```

`r newthought("Empirical justification")`
  
The functional forms chosen by Svensson or Nelson and Siegel are comforted by empirical observations. Indeed, a principal components analysis of daily changes in zero-coupon rates shows that actual changes are accounted for by three factors that have shapes consistent with the decomposition illustrated in Figure\@ref{fig:svensson-3}.
 
The data set consists in AAA euro-zone government rates, observed from December 29, 2006 to July 24, 2009.
 
```{r  pca-analysis, echo=TRUE, message=FALSE}
n <- nrow(ECBYieldCurve)
daily.zc.change <- as.matrix(ECBYieldCurve[2:n,]) - 
as.matrix(ECBYieldCurve[1:(n-1),])
pca <- princomp(daily.zc.change, center=TRUE, scale=TRUE)
```


Figure \@ref{fig:pca-var} shows the contribution of the first 6 factors to the total variance, and shows that of the bulk of the variance (about 90\%) is captured by the first 3 factors.


```{r pca-var, echo=F, fig.cap="Contribution (in %) of the first 6 principal components to total variance (norm of the covariance matrix of changes in zero-rates"}
barplot(pca$sdev[1:6]^2/sum(pca$sdev^2))
```

The factor loadings show the sensitivity of zero-coupon of various maturities to factor movements. The first three factors are displayed in Figure\@ref{fig:pca-2}. 

```{r pca-2, echo=FALSE,warning=FALSE,fig.cap="Factor loadings for first three principal components."}
mat <- c(0.25, 0.5, seq(30))
plot(mat, pca$loadings[,1], ylim=c(min(pca$loadings[,1:3]), max(pca$loadings[,1:3])), xlab='Maturity (Yr)',
ylab='', col='red', lty=1)
lines(mat, pca$loadings[,2], col='blue', lty=1)
lines(mat, pca$loadings[,3], col='green', lty=1)
```

Focusing on the shapes for maturities over 2 years, 
the first factor (red) can be interpreted as a level factor, since it affects all zero-coupons in a similar way (except for the short maturities). 
The second factor (blue) can be interpreted as a slope or rotation factor. 
Finally, the third factor (green) can be interpreted as a curvature factor. We find a good qualitative agreement with the decomposition of Figure \@ref{fig:svensson-3}. 

## Bootstrapping a Zero-Coupon Libor Curve

Given a set of LIBOR deposit rates and swap rates, we would like to compute
 the zero-coupon curve $z_i$. Again, this is best explained by an example. We use the `r to_index("LiborRates", "data")`, which contains a time series of deposit rates and swap rates for various maturities. The construction of this data set was described in Chapter \@ref{ch:data-sources}. The data for February 23, 2010, is  displayed in Table \@ref{tab:libor-1}.
  
```{r ZC-123, echo=FALSE, eval=TRUE, fig.cap="One row of the LiborRates data set."}
dtTrade <- as.timeDate("2010-02-23")
l <- ts_libor[dtTrade]/100
print(l)
```

`r newthought("Bootstrap from the deposit rates")`

The first part of the calculation is simple, and amounts to converting the 
LIBOR deposit rates ($\frac{\mbox{ACT}}{360}$) into zero-coupon continuous rates ($\frac{\mbox{ACT}}{365}$),
using the relation:

\[
(1 + \frac{d}{360} l) = e^{z \frac{d}{365}}
\]

where:

: $d$
Actual number of days from settlement date to maturity

: $l$
Deposit rate

: $z$
Zero-coupon rate

The 3 and 6 months zero-coupon rates are:
 
```{r zcRate, echo=TRUE}
dep.to.zc <- function(nb.days, deposit.rate) {
   log(1 + deposit.rate*nb.days/360)*365/nb.days}

dtTrade <- as.timeDate("2010-02-23")
dtSettlement <- addDays(dtTrade,2)
dtPayment <- timeSequence(from=dtSettlement, 
                           by='3 months',
                           length.out=3)[2:3]
 
dd = as.numeric(difftime(dtPayment, dtSettlement, units="days"))
r <- ts_libor[dtTrade, c("Libor3M", "Libor6M")]/100
zc.deposits = dep.to.zc(dd, r)
print(zc.deposits)
```

 
`r newthought("Bootstrap from the swap rates")`
 
At one year and beyond, the calculation proceeds recursively, one year at a time. The swap rate is a par yield,
meaning that the fixed leg of a swap has a PV of 100.
We can use the same logic as in the previous section,
and compute recursively the zero-coupon rate for each maturity.
 
The calculation is performed by the following script, where we linearly interpolate the swap rates for the missing maturities.
 
The linear interpolation is performed with:
 
```{r libor-interp, echo=TRUE}
 libor.mat <- c(1,2,3,4,5,7,10)
 libor.swap.rates <- ts_libor[dtTrade, 4:10]/100
 mat <- seq(10)
 libor.swap.rates <- approx(libor.mat, libor.swap.rates, xout=mat)$y
```
 
The recursive computation of the zero-coupon rates is
finally done with:
 
```{r libor-zc, echo=TRUE}
 zc.libor <- rep(0, 10)
 zc.libor[1] <- log(1+libor.swap.rates[1])
 P <- rep(100,length(mat))
 for(i in seq(2,10)) {
   zc.libor[i] <- ZCRate(libor.swap.rates, zc.libor, i) }
``` 
 
The resulting zero-curve is the 
concatenation of the zero-coupon deposit rates and the zero-coupon swap rates.
 
```{r zero-libor, echo=FALSE, fig.cap="Zero-coupon curve  bootstrapped from Libor deposits and swaps"}
 
 mat <- c(0.25, 0.50, seq(10))
 z <- c(zc.deposits, zc.libor)
 
 plot(mat, z*100, type='l',col='red', lwd=2, xlab='Maturity (Yr)', ylab='zero-coupon libor rate (%)')
```
````{comment}
`r newthought("Using QuantLib")`

 A more accurate calculation method is provided by QuantLib, 
 a pricing library written in C++, and interfaced to \RR. The following script performs the same calculation as above, using QuantLib. The main difference is the ability to define in interpolation function between nodes. Here,
we choose a log-linear interpolation scheme on the discount function.
 
```{r ZC, echo=T} 
dt = .25
params <- list(tradeDate=dtTrade,
settleDate=dtSettlement,
dt=dt,
interpWhat='discount',
interpHow='loglinear')
l <- ts_libor[dtTrade]/100 
tsQuotes <- list(d1m = l$Libor1M,
d3m = l$Libor3M,
d6m = l$Libor6M,
s2y = l$Swap2Y,
s3y = l$Swap3Y,
s5y = l$Swap5Y,
s10y = l$Swap10Y,
s30y = l$Swap30Y)
 
horizon <- 10
times <- seq(0,horizon,.1)
curve <- DiscountCurve(params, tsQuotes, times)
```
 
```{r, eval=TRUE, echo=FALSE, fig.cap="Libor Zero-coupon curves by QuantLib"}
plot(curve$times, curve$zerorates*100, type='l', col=2, xlab='Maturity (Yr)', ylab='zero-coupon libor rate (%)')
```
````