Problem 1

Let x1,x2,x3,x4,x5,x6,x7,x8,x9,x10 represent the initial investment in each of the ten bonds. The objective is to minimize the initial cost and matching the liability cash flows for the coming 8 years.

\(min \ 102x_1 + 99x_2 + 101x_3 + 98x_4 + 98x_5 + 104x_6 + 100x_7 + 101x_8 + 102x_9 + 94x_{10}\)

\(s.t.\)

\(Year\ 1:\ (100+5)x_1 + 3.5x_2 + 5x_3 + 3.5x_4 + 4x_5 + 9x_6 + 6x_7 + 8x_8 + 9x_9 + 7x_{10}\ =\ 12,000\)

\(Year\ 2:\ (100+3.5)x_2 + (100+5)x_3 + 3.5x_4 + 4x_5 + 9x_6 + 6x_7 + 8x_8 + 9x_9 + 7x_{10}\ =\ 18,000\)

\(Year\ 3:\ (100+3.5)x_4 + 4x_5 + 9x_6 + 6x_7 + 8x_8 + 9x_9 + 7x_{10}\ =\ 20,000\)

\(Year\ 4:\ (100+4)x_5 + 9x_6 + 6x_7 + 8x_8 + 9x_9 + 7x_{10}\ =\ 20,000\)

\(Year\ 5:\ (100+9)x_6 + (100+6)x_7 + 8x_8 + 9x_9 + 7x_{10}\ =\ 16,000\)

\(Year\ 6:\ (100+8)x_8 + 9x_9 + 7x_{10}\ =\ 15,000\)

\(Year\ 7:\ (100+9)x_9 + 7x_{10}\ =\ 12,000\)

\(Year\ 8:\ (100+7)x_{10}\ =\ 10,000\)

Problem 2

library(lpSolveAPI)

my.lp <- make.lp(8,10) # 8 Constraints, 10 decision variables
print (my.lp)
## Model name: 
##   a linear program with 10 decision variables and 8 constraints
Colnames=c('x1','x2','x3','x4','x5','x6','x7','x8','x9','x10')
Rownames=c('Y1','Y2','Y3','Y4','Y5','Y6','Y7','Y8')
dimnames(my.lp)<-list(Rownames,Colnames)

add.constraint(my.lp,c(105,3.5,5,3.5,4,9,6,8,9,7), "=",12000)  #Y1
add.constraint(my.lp,c(0,103.5,105,3.5,4,9,6,8,9,7), "=",18000) #Y2
add.constraint(my.lp,c(0,0,0,103.5,4,9,6,8,9,7), "=",20000) #Y3
add.constraint(my.lp,c(0,0,0,0,104,9,6,8,9,7), "=",20000) #Y4
add.constraint(my.lp,c(0,0,0,0,0,109,106,8,9,7), "=",16000) #Y5
add.constraint(my.lp,c(0,0,0,0,0,0,0,108,9,7), "=",15000) #Y6
add.constraint(my.lp,c(0,0,0,0,0,0,0,0,109,7), "=",12000) #Y7
add.constraint(my.lp,c(0,0,0,0,0,0,0,0,0,107), "=",10000) #Y8

set.objfn(my.lp, c(102,99,101,98,98,104,100,101,102,94))
lp.control(my.lp,sense='min')
## $anti.degen
## [1] "none"
## 
## $basis.crash
## [1] "none"
## 
## $bb.depthlimit
## [1] -50
## 
## $bb.floorfirst
## [1] "automatic"
## 
## $bb.rule
## [1] "pseudononint" "greedy"       "dynamic"      "rcostfixing" 
## 
## $break.at.first
## [1] FALSE
## 
## $break.at.value
## [1] -1e+30
## 
## $epsilon
##       epsb       epsd      epsel     epsint epsperturb   epspivot 
##      1e-10      1e-09      1e-12      1e-07      1e-05      2e-07 
## 
## $improve
## [1] "dualfeas" "thetagap"
## 
## $infinite
## [1] 1e+30
## 
## $maxpivot
## [1] 250
## 
## $mip.gap
## absolute relative 
##    1e-11    1e-11 
## 
## $negrange
## [1] -1e+06
## 
## $obj.in.basis
## [1] TRUE
## 
## $pivoting
## [1] "devex"    "adaptive"
## 
## $presolve
## [1] "none"
## 
## $scalelimit
## [1] 5
## 
## $scaling
## [1] "geometric"   "equilibrate" "integers"   
## 
## $sense
## [1] "minimize"
## 
## $simplextype
## [1] "dual"   "primal"
## 
## $timeout
## [1] 0
## 
## $verbose
## [1] "neutral"
write.lp(my.lp,'Project1.lp',type='lp')
solve(my.lp)
## [1] 0
print (get.objective(my.lp))
## [1] 93944.5
print (get.variables(my.lp))
##  [1]  62.13613   0.00000 125.24293 151.50508 156.80776 123.08007   0.00000
##  [8] 124.15727 104.08986  93.45794

Problem 3

Testing of our function dedicate_g3(P,C,M,L):

P = c(102,99,101,98,98,104,100,101,102,94)
C = c(5,3.5,5,3.5,4,9,6,8,9,7)
M = c(1,2,2,3,4,5,5,6,7,8)
L = c(12000,18000,20000,20000,16000,15000,12000,10000)

Portfolio<-dedicate_g3(P,C,M,L)

print (get.objective(Portfolio))
## [1] 93944.5
print (get.variables(Portfolio))
##  [1]  62.13613   0.00000 125.24293 151.50508 156.80776 123.08007   0.00000
##  [8] 124.15727 104.08986  93.45794

Problem 4

Get US Treasury Bonds quotes from http://www.wsj.com/mdc/public/page/2_3020-treasury.html

* STARTING DATE: pulled on Feb 8, 2016, which should be treated as the STARTING DATE.

* MATURITY DATES: chosen to match the date of liability stream. For some of the maturity dates there are multiple bonds available, we just randomly pick one bond for each maturity date.

* PRICES: calculated the average of bid and ask price (\(price\ =\ (bid\ +\ ask)\ /\ 2\)) as the price of a bond, which is a widely used approach in academics.

raw_quotes = read.csv("US Treasury Bonds.csv")
print (raw_quotes)
##      Maturity Coupon      Bid    Asked    Chg Asked.yield
## 1   6/30/2016  0.500 100.0625 100.0781 0.0234       0.299
## 2   6/30/2016  1.500 100.4609 100.4766 0.0391       0.277
## 3   6/30/2016  3.250 101.1328 101.1484 0.0078       0.302
## 4  12/31/2016  0.625 100.0000 100.0391  unch.       0.605
## 5  12/31/2016  0.875 100.2813 100.2969 0.0625       0.540
## 6  12/31/2016  3.250 102.3672 102.3828 0.0469       0.562
## 7   6/30/2017  0.625 100.0078 100.0234 0.0938       0.608
## 8   6/30/2017  0.750 100.1719 100.1875 0.0859       0.614
## 9   6/30/2017  2.500 102.5859 102.6016 0.0938       0.617
## 10 12/31/2017  0.750 100.1328 100.1484 0.1094       0.671
## 11 12/31/2017  1.000 100.6172 100.6328 0.1016       0.663
## 12 12/31/2017  2.750 103.8750 103.8906 0.1406       0.675
## 13  6/30/2018  1.375 101.5156 101.5313 0.1641       0.728
## 14  6/30/2018  2.375 103.9375 103.9531 0.1797       0.704
## 15 12/31/2018  1.375 101.5781 101.5938 0.2188       0.816
## 16 12/31/2018  1.500 101.9453 101.9609 0.2422       0.812
## 17  6/30/2019  1.000 100.2656 100.2813 0.2656       0.916
## 18  6/30/2019  1.625 102.3359 102.3516 0.2734       0.919
## 19 12/31/2019  1.125 100.3750 100.3906 0.3359       1.022
## 20 12/31/2019  1.625 102.2656 102.2813 0.3047       1.025
## 21  6/30/2020  1.625 102.2109 102.2266 0.3828       1.104
## 22  6/30/2020  1.875 103.3125 103.3281 0.4063       1.096
## 23 12/31/2020  1.750 102.8125 102.8281 0.4453       1.153
## 24 12/31/2020  2.375 105.6328 105.6484 0.4609       1.183
## 25  6/30/2021  2.125 104.4688 104.4844 0.4688       1.262
## 26 12/31/2021  2.125 104.4297 104.4453 0.5938       1.338
P4 = (raw_quotes$Bid + raw_quotes$Asked) / 2
C4 = raw_quotes$Coupon
M4 = c(1,1,1,
       2,2,2,
       3,3,3,
       4,4,4,
       5,5,
       6,6,
       7,7,
       8,8,
       9,9,
       10,10,
       11,
       12)
L4 = c(9,9,10,10,6,6,9,9,10,10,5,3) * 1000000

Portfolio_4<-dedicate_g3(P4,C4,M4,L4)

print (get.objective(Portfolio_4))
## [1] 87582418
print (get.variables(Portfolio_4))
##  [1]     0.00     0.00 70495.75     0.00     0.00 72786.86     0.00
##  [8]     0.00 85152.43     0.00     0.00 87281.24     0.00 49681.48
## [15]     0.00 50861.41     0.00 81624.33     0.00 82950.73     0.00
## [22] 94298.68     0.00 96066.78 48348.36 29375.76
## $duals
##  [1] 0.9795700 0.9606915 0.9535912 0.9335738 0.9265462 0.9342083 0.9161176
##  [8] 0.9007770 0.8760568 0.8374643 0.8312077 0.8135292 1.6235174 1.0423975
## [15] 0.0000000 2.7377364 2.5222210 0.0000000 2.8478237 2.6501922 0.0000000
## [22] 3.9126529 3.4401963 0.0000000 2.3321226 0.0000000 0.3438726 0.0000000
## [29] 2.0573866 0.0000000 1.8618878 0.0000000 0.9937331 0.0000000 2.9413229
## [36] 0.0000000 0.0000000 0.0000000
## 
## $dualsfrom
##  [1]  1.721314e+06  1.484757e+06  1.271876e+06  1.031852e+06  9.138590e+05
##  [6]  8.375669e+05  7.049273e+05  5.701324e+05  3.933224e+05  1.651638e+05
## [11]  6.242350e+04  0.000000e+00 -1.000000e+30 -1.000000e+30 -1.000000e+30
## [16] -2.862950e+06 -3.164313e+06 -1.000000e+30 -4.108333e+06 -4.401785e+06
## [21] -1.000000e+30 -3.957480e+06 -4.522834e+06 -1.000000e+30 -8.102940e+06
## [26] -1.000000e+30 -4.129947e+07 -1.000000e+30 -8.394116e+06 -1.000000e+30
## [31] -1.066315e+07 -1.000000e+30 -2.172617e+07 -1.000000e+30 -8.896867e+06
## [36] -1.000000e+30 -1.000000e+30 -1.000000e+30
## 
## $dualstill
##  [1] 1.000000e+30 2.402382e+08 3.181250e+08 2.978167e+08 3.471764e+08
##  [6] 3.501622e+08 3.318506e+08 3.370969e+08 2.996823e+08 2.441281e+08
## [11] 2.722331e+08 2.402941e+08 7.242473e+04 7.171119e+04 1.000000e+30
## [16] 7.468565e+04 7.450055e+04 1.000000e+30 8.673912e+04 8.663150e+04
## [21] 1.000000e+30 8.901387e+04 8.879354e+04 1.000000e+30 5.017155e+04
## [26] 1.000000e+30 5.092412e+04 1.000000e+30 8.212943e+04 1.000000e+30
## [31] 8.336087e+04 1.000000e+30 9.453065e+04 1.000000e+30 9.665687e+04
## [36] 1.000000e+30 1.000000e+30 1.000000e+30

The first 12 duals correspond to the 12 constraints we have, and the rest 26 correspond to the non-negative constraints.

The first 12 duals can be interpreted as follows: if liability in period i is reduced by $1, the intial investment in period 0 will reduce by \(dual amount\). For example, for the first dual, if liability in year 1 reduces by $1, our initial payout would reduce by $0.97.