*
* BONDSPLINE.RPF
* Estimates a yield curve using a cubic spline approximation to the
* discount function, as described in McCulloch(1971), "Measuring the
* Term Structure of Interest Rates", Journal of Business, vol 44, pp
* 19-31.
*
compute nbonds=20
all nbonds
open data bonds.xls
data(format=xls,org=cols) / coupon value maturity
*
* The data set consists of US Treasury bonds, which have semi-annual
* coupons, but with the coupon stated in annual terms. The maturity
* series also provides annual information. The coupon is divided by two
* and the maturity multiplied by 2 to give values for the actual coupon
* period.
*
set coupon   = coupon / 2
set maturity = maturity * 2
*
* This models P(n)=discount factor for maturity n as a cubic spline.
*
* These are for the unrestricted linear, quadratic and cubic terms. (The
* constant is 1).
*
dec vect[series] baseregs(3)
*
* These are for the knots
*
compute nknots=4
dec vect knots(nknots)
compute knots=||2.0,4.0,6.0,8.0||
dec vect[series] knotregs(nknots)
clear(zeros) knotregs baseregs
*
* This is for the value of the bond adjusted for the constant terms in
* P(n) and for the partial initial coupon.
*
set valueadj = value
*
* Loop over the bonds to generate the required set of regressors.
*
do i=1,nbonds
   compute mdate=maturity(i)
   *
   * Walk backwards through the payments. The last one will include the
   * face value as well as the final coupon.
   *
   compute cdate=mdate
   while (cdate > 0.0) {
      if cdate==mdate
         compute receipt=coupon(i)+100.0
      else
         compute receipt=coupon(i)
      compute valueadj(i)=valueadj(i)-receipt
      do j=1,3
         compute baseregs(j)(i)=baseregs(j)(i)+receipt*cdate^j
      end do j
      do k=1,nknots
         compute knotregs(k)(i)=knotregs(k)(i)+receipt*%max(cdate-knots(k),0.0)^3
      end do k
      compute cdate=cdate-1
   }
   *
   * Adjust for simple interest payable by the purchaser for the initial
   * coupon. cdate will be -(fraction of period).
   *
   compute valueadj(i)=valueadj(i)-coupon(i)*cdate
end do i
*
linreg valueadj
# baseregs knotregs
*
* Graph the estimated yield curve from maturities of 0->10 periods (5
* years). After the present value function is evaluated, the maturity is
* switched to years and the annual yield computed.
*
function cspline m
type real cspline m
*
compute cspline=1.0+%beta(1)*m+%beta(2)*m^2+%beta(3)*m^3
do k=1,nknots
   compute cspline=cspline+%beta(3+k)*%max(m-knots(k),0.0)^3
end do k
end spline
*
set testm 1 50  = .20*t
set pvalue 1 50 = cspline(testm)
*
* Convert the discount factor into an annualized yield curve
*
set testm 1 50  = testm/2
set yield 1 50  = -100.0*log(pvalue)/testm
scatter(style=line,window="Yield Curve",hlabel="Years",vlabel="Yield")
# testm yield