diff options
Diffstat (limited to 'R')
| -rw-r--r-- | R/TAGS | 321 | ||||
| -rw-r--r-- | R/tranche_functions.R | 531 |
2 files changed, 691 insertions, 161 deletions
@@ -1,4 +1,8 @@ +./backtest.R,31 +aargs8,156 +aalldates28,851 + ./backtest_tranches.R,69 aargs1,0 bbps18,543 @@ -67,15 +71,16 @@ pprep.data62,2278 ssqlstr65,2441 tticker_company75,2878 -./build_portfolios.R,167 +./build_portfolios.R,184 aargs4,38 ccode.dir5,77 -iindex16,431 -ccalibration.date35,941 -gglobal.params38,1013 -ccusipdata40,1096 -ccashspread.discount41,1122 -ccurrdealnames42,1147 +ccalibration.date36,972 +iindex38,1043 +iindex39,1072 +gglobal.params40,1121 +ccusipdata42,1204 +ccashspread.discount43,1230 +ccurrdealnames44,1255 ./build_scenarios.R,234 hhostname5,56 @@ -116,28 +121,26 @@ EELmat56,1917 rrisk.numbers85,3129 ddata93,3487 -./calibrate_tranches_BC.R,47 -aargs2,32 -uupdateflag23,806 -ttenor24,826 +./calibrate_tranches_BC.R,35 +ooption_list12,321 +aargs21,880 -./calibrate_tranches_MF.R,120 -aargs6,74 -ccode.dir14,238 -nn.int33,769 -NNgrid35,828 -iindex36,841 -iindex37,905 -iindex39,973 -ddist42,1100 +./calibrate_tranches_MF.R,106 +aargs5,69 +ccode.dir13,233 +iindex32,889 +iindex33,925 +iindex34,967 +iindex36,1035 +ddist39,1162 ./calibration.R,141 bbuildSC3,35 gget.cdsSchedule14,418 sset.singlenamesdata27,918 sset.tranchedata52,2173 -bbuild.skew100,4484 -bbuild.MFdist132,5931 +bbuild.skew105,4655 +bbuild.MFdist137,6102 ./cds_functions_generic.R,1451 sabstractcurve12,258 @@ -188,19 +191,19 @@ ccdshazardrate458,21474 bbondhazardrate.shaped492,22797 ttweakportfolio521,23947 iindexpv541,24665 -iindexduration568,25871 -iindexspread574,26100 -iindextheta587,26601 -pportfoliospread597,26962 -pportfolioduration623,28122 -ttweakcurves628,28347 -ssurvivalProbability1648,29172 -ssurvivalProbability.exact671,30090 -SSP693,30792 -SSPmatrix700,31025 -DDP2716,31589 -SSPmatrix2733,32349 -ccreditcurve.maturity750,33035 +iindexduration580,26347 +iindexspread586,26576 +iindextheta599,27077 +pportfoliospread609,27438 +pportfolioduration635,28598 +ttweakcurves640,28823 +ssurvivalProbability1660,29652 +ssurvivalProbability.exact683,30570 +SSP705,31272 +SSPmatrix712,31505 +DDP2728,32069 +SSPmatrix2745,32829 +ccreditcurve.maturity762,33515 ./cds_utils.R,179 ttoday3,22 @@ -238,9 +241,9 @@ ccreditIndex12,192 `$<-.creditIndex19,474 cc.creditIndex26,606 pprint.creditIndex30,751 -lload.index54,1666 -ccsvheaders61,1911 -ttocsv80,2730 +lload.index55,1742 +ccsvheaders68,2148 +ttocsv88,3024 ./deal_pricer.R,92 ddrv2,21 @@ -296,31 +299,32 @@ cclipr107,2749 ssclipr111,2798 iinverse115,2865 -./intex_deal_functions.R,657 +./intex_deal_functions.R,679 hhostname8,127 ggetdealdata18,309 -ggetcollateral32,1108 -llistdealnames42,1441 -ccusip.data47,1595 -rrecovery59,2045 -ddealnamefromcusip79,2831 -ccusipsfromdealnames87,3092 -ffithazardrate.fast94,3360 -vvanillabondprice106,3892 -ddvanillabondprice121,4564 -ffithazardrate129,4916 -sstackcurve141,5379 -bbuildSC.matured167,6617 -bbuildSC179,7088 -bbuildSC.portfolio247,10417 -ccdrfromscenarios275,11965 -rrecoveryfromscenarios285,12353 -rrecoveryfromscenarios.fast310,13326 -sseverityfromscenarios316,13528 -gget.reinvassets334,14220 -ggetdealschedule347,14641 -iintexportfolio.forwardprice368,15496 -ccompute.reinvprices398,17113 +ggetcollateral32,1092 +llistdealnames42,1425 +ccusip.data47,1579 +rrecovery59,2029 +ddealnamefromcusip79,2815 +ccusipsfromdealnames87,3076 +ffithazardrate.fast94,3340 +vvanillabondprice106,3872 +ddvanillabondprice121,4544 +ffithazardrate129,4896 +sstackcurve141,5359 +bbuildSC.matured167,6597 +bbuildSC179,7068 +bbuildSC.portfolio247,10397 +ccdrfromscenarios275,11943 +rrecoveryfromscenarios285,12331 +rrecoveryfromscenarios.fast310,13304 +sseverityfromscenarios316,13506 +gget.reinvassets334,14198 +ggetpayday347,14619 +ggetdealschedule366,15326 +iintexportfolio.forwardprice385,15999 +ccompute.reinvprices415,17616 ./latestprices.R,453 cconn2,15 @@ -361,25 +365,36 @@ aaddcurves21,819 mmultcurves29,1100 sspreadfrombondprice37,1382 -./load_cds_data.R,71 +./load_cds_data.R,220 ddownload.cdscurves4,39 wwrite.tranchedata21,749 bbbgConn35,1319 +hhy2239,1440 +hhy2140,1536 +hhy1941,1632 +hhy1742,1728 +iig2143,1824 +iig1944,1920 +iig945,2016 +ddrv48,2133 +ddbCon49,2163 +hhy1752,2290 -./load_cf.R,267 +./load_cf.R,281 aargs7,122 -iindex28,785 -ccalibration.date30,814 -ddm33,887 -ssanitize.column34,895 -pprocesszipfiles41,1112 -ggetconfig69,2250 -ggetdealcf79,2585 -ggetcusipcf167,6889 -ccompute.delta227,10078 -ccfdata255,11227 -ccusipdata259,11432 -ssave.dir260,11489 +iindex29,843 +iindex30,872 +ccalibration.date32,915 +ddm35,988 +ssanitize.column36,996 +pprocesszipfiles43,1213 +ggetconfig71,2351 +ggetdealcf81,2686 +ggetcusipcf163,6560 +ccompute.delta224,9807 +ccfdata252,10956 +ccusipdata256,11161 +ssave.dir257,11218 ./load_futures_data.R,180 bbbgCon11,202 @@ -430,72 +445,39 @@ lloan.data268,2110 mmaturities69,2170 ttest270,2222 -lossdistrib/R/tranche_functions.R,1711 -GGHquad13,733 -llossdistrib22,903 -llossdistrib.fft34,1170 -llossdistrib245,1554 -llossdistrib2.truncated74,2344 -rrecovdist109,3524 -llossdist.joint142,4594 -llossdist.prepay.joint181,6111 -llossdistC242,8707 -llossdistCblas248,8973 -llossdistCZ254,9248 -llossdistC.truncated262,9589 -rrecovdistC268,9851 -llossdistC.joint274,10066 -llossdistC.jointblas280,10355 -llossdistC.jointZ286,10653 -llossdistC.prepay.jointblas299,11210 -llossdistC.prepay.jointZ306,11545 -llossrecovdist322,12218 -llossrecovdist.term343,12999 -llossrecovdist.joint.term363,13720 -ddist.transform391,14724 -sshockprob413,15423 -sshockseverity432,16128 -ddshockprob441,16452 -ddqnorm445,16556 -ffit.prob449,16604 -ffit.probC474,17287 -sstochasticrecov480,17476 -sstochasticrecovC491,17926 -ppos498,18223 -ttrancheloss502,18261 -ttrancherecov506,18330 -ttranche.cl510,18407 -ttranche.cl.scenarios530,19204 -ffunded.tranche.pv548,19907 -ttranche.pl567,20585 -ttranche.pl.scenarios584,21075 -ttranche.pv600,21510 -ttranche.pv.scenarios604,21654 -aadjust.attachments609,21822 -BBClossdist615,22012 -BBClossdistC642,22921 -BBCtranche.legs653,23497 -BBCtranche.pv672,24347 -MMFtranche.pv699,25230 -aadjust.skew718,25798 -ttheta.adjust.skew767,27565 -BBCtranche.theta776,27916 -BBCtranche.delta793,28604 -MMFtranche.delta825,29978 -BBCtranche.corr01858,31354 -EEL866,31675 -BBCindex.pv881,32372 -EELtrunc899,33124 -PProbtrunc907,33544 -BBCstrikes916,33851 -ttranche.factor927,34220 -MMFupdate.prob933,34450 -MMFupdate.probC945,34854 -MMFlossrecovdist.prepay958,35362 -MMFlossdist.joint985,36519 -MMFlossdist.prepay.joint1017,37873 -MMFrecovery1044,39008 -MMFlossdist1058,39526 +lossdistrib/R/tranche_functions.R,850 +GGHquad25,1100 +llossdistrib44,1729 +llossdistrib.fft68,2578 +llossdistrib286,3252 +llossdistrib2.truncated109,3791 +rrecovdist144,4971 +llossdist.joint177,6050 +llossdist.prepay.joint216,7567 +llossdistC277,10163 +llossdistCblas283,10429 +llossdistCZ289,10704 +llossdistC.truncated297,11045 +rrecovdistC303,11307 +llossdistC.joint309,11522 +llossdistC.jointblas315,11811 +llossdistC.jointZ321,12109 +llossdistC.prepay.jointblas334,12666 +llossdistC.prepay.jointZ341,13001 +llossrecovdist357,13674 +llossrecovdist.term370,14169 +llossrecovdist.joint.term390,14890 +ddist.transform407,15573 +sshockprob429,16272 +sshockseverity448,16977 +ddshockprob457,17301 +ddqnorm461,17405 +ffit.prob465,17453 +ffit.probC490,18136 +sstochasticrecov496,18325 +sstochasticrecovC507,18775 +BBClossdist514,19072 +BBClossdistC541,19981 ./LS.R,139 rr1,0 @@ -634,26 +616,43 @@ pprotectionsize86,1551 ccl90,1614 ppl96,1766 -./tranches_RV_BC_et.R,340 -nn.int13,390 -NNgrid15,449 -aaux17,463 -aalldates22,603 -bbus.dates23,677 -oorigpv25,778 -ffirstpv26,792 -ssecondpv27,807 -tthirdpv28,823 -rrandates29,838 -iindex.name42,1085 -iindex.name143,1106 -iindextenor44,1128 -iindex1tenor45,1148 -kk50,1216 -bbigmat51,1223 -ttemp117,3637 -ttranche.names118,3699 -bbigmat119,3757 +./tranche_functions.R,805 +ppos2,21 +ttrancheloss6,59 +ttrancherecov10,128 +ttranche.cl14,205 +ttranche.cl.scenarios34,1002 +ffunded.tranche.pv52,1705 +ttranche.pl71,2383 +ttranche.pl.scenarios88,2873 +ttranche.pv104,3308 +ttranche.pv.scenarios108,3452 +aadjust.attachments113,3620 +BBCtranche.legs119,3810 +BBCtranche.pv138,4660 +MMFtranche.pv165,5543 +aadjust.skew184,6111 +ttheta.adjust.skew233,7878 +BBCtranche.theta242,8229 +BBCtranche.delta259,8922 +MMFtranche.delta291,10296 +BBCtranche.corr01324,11672 +EEL332,11993 +BBCindex.pv347,12690 +EELtrunc365,13458 +PProbtrunc373,13878 +BBCstrikes382,14185 +ttranche.factor393,14509 +MMFupdate.prob399,14739 +MMFlossrecovdist.prepay412,15205 +MMFlossdist.joint439,16362 +MMFlossdist.prepay.joint471,17716 +MMFrecovery498,18851 +MMFlossdist512,19369 + +./tranches_RV_BC.R,36 +ooption_list11,282 +aargs25,1053 ./transactions.R,98 ttransdir1,0 diff --git a/R/tranche_functions.R b/R/tranche_functions.R new file mode 100644 index 00000000..0bccdf02 --- /dev/null +++ b/R/tranche_functions.R @@ -0,0 +1,531 @@ +library(lossdistrib) +pos <- function(x){ + pmax(x, 0) +} + +trancheloss <- function(L, K1, K2){ + pos(L - K1) - pos(L - K2) +} + +trancherecov <- function(R, K1, K2){ + pos(R - 1 + K2) - pos(R - 1 +K1) +} + +tranche.cl <- function(L, R, cs, K1, K2, Ngrid=nrow(L), scaled=FALSE){ + ## computes the couponleg of a tranche + ## if scaled is TRUE, scale it by the size of the tranche (K2-K1) + ## can make use of the fact that the loss and recov distribution are + ## truncated (in that case nrow(L) != Ngrid + if(K1==K2){ + return( 0 ) + }else{ + support <- seq(0, 1, length=Ngrid)[1:nrow(L)] + size <- K2 - K1 - crossprod(trancheloss(support, K1, K2), L) - + crossprod(trancherecov(support, K1, K2), R) + sizeadj <- as.numeric(0.5 * (size + c(K2-K1, size[-length(size)]))) + if(scaled){ + return( 1/(K2-K1) * crossprod(sizeadj * cs$coupons, cs$df) ) + }else{ + return( crossprod(sizeadj * cs$coupons, cs$df) ) + } + } +} + +tranche.cl.scenarios <- function(l, r, cs, K1, K2, scaled=FALSE){ + ## computes the couponleg of a tranche for one scenario + ## if scaled is TRUE, scale it by the size of the tranche (K2-K1) + ## can make use of the fact that the loss and recov distribution are + ## truncated (in that case nrow(L) != Ngrid + if(K1==K2){ + return( 0 ) + }else{ + size <- K2 - K1 - trancheloss(l, K1, K2) - trancherecov(r, K1, K2) + sizeadj <- as.numeric(0.5 * (size + c(K2-K1, size[-length(size)]))) + if(scaled){ + return( 1/(K2-K1) * crossprod(sizeadj * cs$coupons, cs$df) ) + }else{ + return( crossprod(sizeadj * cs$coupons, cs$df) ) + } + } +} + +funded.tranche.pv <- function(L, R, cs, K1, K2, scaled = FALSE){ + if(K1==K2){ + return(0) + }else{ + size <- K2 - K1 -trancheloss(L, K1, K2) - trancherecov(R, K1, K2) + sizeadj <- as.numeric(0.5 * (size + c(K2-K1, size[-length(size)]))) + interest <- crossprod(sizeadj * cs$coupons, cs$df) + principal <- diff(c(0, trancherecov(R, K1, K2))) + principal[length(principal)] <- principal[length(principal)] + size[length(size)] + principal <- crossprod(cs$df, principal) + if(scaled){ + pv <- (interest + principal)/(K2-K1) + }else{ + pv <- (interest + principal) + } + return(pv) + } +} + +tranche.pl <- function(L, cs, K1, K2, Ngrid=nrow(L), scaled=FALSE){ + ## computes the protection leg of a tranche + ## if scaled + if(K1==K2){ + return(0) + }else{ + support <- seq(0, 1, length=Ngrid)[1:nrow(L)] + cf <- K2 - K1 - crossprod(trancheloss(support, K1, K2), L) + cf <- c(K2 - K1, cf) + if(scaled){ + return( 1/(K2-K1) * crossprod(diff(cf), cs$df)) + }else{ + return( crossprod(diff(cf), cs$df)) + } + } +} + +tranche.pl.scenarios <- function(l, cs, K1, K2, scaled=FALSE){ + ## computes the protection leg of a tranche + ## if scaled + if(K1==K2){ + return(0) + }else{ + cf <- K2 - K1 - trancheloss(l, K1, K2) + cf <- c(K2 - K1, cf) + if(scaled){ + return( 1/(K2-K1) * as.numeric(crossprod(diff(cf), cs$df))) + }else{ + return( as.numeric(crossprod(diff(cf), cs$df))) + } + } +} + +tranche.pv <- function(L, R, cs, K1, K2, Ngrid=nrow(L)){ + return( tranche.pl(L, cs, K1, K2, Ngrid) + tranche.cl(L, R, cs, K1, K2, Ngrid)) +} + +tranche.pv.scenarios <- function(l, r, cs, K1, K2){ + return( tranche.pl.scenarios(l, cs, K1, K2, TRUE) + + tranche.cl.scenarios(l, r, cs, K1, K2, TRUE)) +} + +adjust.attachments <- function(K, losstodate, factor){ + ## computes the attachments adjusted for losses + ## on current notional + return( pmin(pmax((K-losstodate)/factor, 0),1) ) +} + +BCtranche.legs <- function(index, K, rho, complement=FALSE){ + ## computes the protection leg and couponleg of a 0-K tranche + ## if complement==TRUE, computes the protection leg and coupon leg of a K-1 tranche + if((K==0 && !complement) || (K==1 && complement)){ + return(list(cl=0, pl=0)) + }else if((K==1 && !complement) || (K==0 && complement)){ + return(BCindex.pv(index)) + }else{ + dist <- BClossdistC(index$defaultprob, index$issuerweights, index$recov, rho, index$Z, index$w, index$N) + if(complement){ + return(list(cl=tranche.cl(dist$L, dist$R, index$cs, K, 1), + pl=tranche.pl(dist$L, index$cs, K, 1))) + }else{ + return(list(cl=tranche.cl(dist$L, dist$R, index$cs, 0, K), + pl=tranche.pl(dist$L, index$cs, 0, K))) + } + } +} + +BCtranche.pv <- function(index, protection=FALSE, complement=FALSE){ + ## computes the protection leg, couponleg, and bond price of a tranche + ## in the base correlation setting + ## if complement=FALSE compute the pvs starting from 0 (bottom-up skew) + ## if complement=TRUE compute the pvs starting from 1 (top-down skew) + pl <- rep(0, length(index$rho)) + cl <- rep(0, length(index$rho)) + for(i in seq_along(index$rho)){ + temp <- BCtranche.legs(index, index$K[i], index$rho[i], complement) + pl[i] <- temp$pl + cl[i] <- temp$cl + } + dK <- diff(index$K) + plvec <- diff(pl)/dK + clvec <- diff(cl)/dK*index$tranches$running + if(complement){ + plvec <- -plvec + clvec <- -clvec + } + if(protection){ + bp <- -plvec-clvec + }else{ + bp <- 1+plvec+clvec + } + return(list(pl=plvec, cl=clvec, bp=bp)) +} + +MFtranche.pv <- function(index, dist, protection=FALSE){ + n.tranches <- length(index$K)-1 + pl <- rep(0, n.tranches) + cl <- rep(0, n.tranches) + for(i in seq.int(n.tranches)){ + pl[i] <- tranche.pl(dist$L, index$cs, index$K[i], index$K[i+1]) + cl[i] <- tranche.cl(dist$L, dist$R, index$cs, index$K[i], index$K[i+1]) + } + dK <- diff(index$K) + plvec <- pl/dK + clvec <- cl/dK*index$tranches$running + if(protection){ + bp <- -plvec-clvec + }else{ + bp <- 1+plvec+clvec + } + return(list(pl=plvec, cl=clvec, bp=bp)) +} + +adjust.skew <- function(index1, index2, method="ATM"){ + #index1 is the index for which we already have computed the skew + #index2 is the index we're mapping to + # if method="ATM", do simple at the money mapping + # method="TLP", do tranche loss proportion mapping + # method="PM", do probability matching + + K1 <- index1$K[-c(1,length(index1$K))] + K2 <- index2$K[-c(1,length(index2$K))] + aux <- function(x, index1, el1, skew, index2, el2, K2){ + return(abs(ELtrunc(index1, x, skew(x))/el1- + ELtrunc(index2, K2, skew(x))/el2)) + } + aux2 <- function(x, index1, skew, index2, K2){ + return(abs(log(Probtrunc(index1, x, skew(x)))- + log(Probtrunc(index2, K2, skew(x))))) + } + if(method %in% c("ATM", "TLP")){ + el1 <- EL(index1) + el2 <- EL(index2) + } + skew <- function(x){ + #we cap the correlation at 0.99 and 0.01 + f <- splinefun(K1, index1$rho[-c(1, length(index1$rho))], "natural") + return(pmax(pmin(f(x), 0.99), 0.01)) + } + if(method=="ATM"){ + K1eq <- el1/el2 * K2 + }else if(method == "TLP"){ + K1eq <- c() + m <- max(K2) + 0.3 + for(K2val in K2){ + prog <- optimize(aux, interval=c(0,m), + index1=index1, el1=el1, skew=skew, + index2=index2, el2=el2, K2=K2val) + K1eq <- c(K1eq, prog$minimum) + } + }else if (method=="PM"){ + K1eq <- c() + m <- max(K2) + 0.25 + for(K2val in K2){ + prog <- optimize(aux2, interval=c(0, m), + index1=index1, skew=skew, index2=index2, K2=K2val) + K1eq <- c(K1eq, prog$minimum) + } + } + return(c(NA, skew(K1eq), NA)) +} + +theta.adjust.skew <- function(index, shortened=4, method="ATM"){ + #ajust the correlation skew by doing ATM mapping on the expected loss + indexshort <- index + N <- nrow(index$cs)-shortened + indexshort$defaultprob <- indexshort$defaultprob[,1:N] + indexshort$cs <- indexshort$cs[1:N,] + return(adjust.skew(index, indexshort, method)) +} + +BCtranche.theta <- function(index, shortened=4, complement=FALSE, method="ATM"){ + temp <- BCtranche.pv(index, complement=complement) + rho.adj <- theta.adjust.skew(index, shortened, method) + if(any(rho.adj[-c(1, length(rho.adj))]<=0)){ + print("probable inverted skew: no adjustment") + }else{ + index$rho <- rho.adj + } + N <- nrow(index$cs) - shortened + index$cs <- index$cs[1:N,] + index$defaultprob <- index$defaultprob[,1:N] + temp2 <- BCtranche.pv(index, complement=complement) + temp3 <- BCtranche.delta(index, complement=complement) + return(data.frame(theta=temp2$bp-temp$bp+index$tranches$running, + fw.delta=temp3$delta)) +} + +BCtranche.delta <- function(index, complement=FALSE){ + ## computes the tranche delta (on current notional) by doing a proportional + ## blip of all the curves + ## if complement is False, then computes deltas bottom-up + ## if complement is True, then computes deltas top-down + eps <- 1e-4 + index$N <- 301 ## for gamma computations we need all the precision we can get + ## we build a lit of 4 indices with various shocks + index.list <- lapply(c(0, eps, -eps, 2*eps), function(x){ + if(x==0){ + return(index) + }else{ + newindex <- index + newindex$portfolio <- tweakportfolio(newindex$portfolio, x) + newindex$defaultprob <- 1 - SPmatrix(newindex$portfolio, length(index$cs$dates)) + return(newindex) + } + }) + bp <- matrix(0, length(index$K)-1, length(index.list)) + indexbp <- rep(0, length(index.list)) + for(j in seq_along(index.list)){ + indexbp[j] <- BCindex.pv(index.list[[j]])$bp + bp[,j] <- BCtranche.pv(index.list[[j]], complement=complement)$bp + } + + deltas <- (bp[,2]-bp[,3])/(indexbp[2]-indexbp[3])*tranche.factor(index)/index$factor + deltasplus <- (bp[,4]-bp[,1])/(indexbp[4]-indexbp[1])*tranche.factor(index)/index$factor + gammas <- (deltasplus-deltas)/(indexbp[2]-indexbp[1])/100 + + return( data.frame(delta=deltas, gamma=gammas) ) +} + +MFtranche.delta <- function(index){ + ## computes the tranche delta (on current notional) by doing a proportional + ## blip of all the curves + ## if complement is False, then computes deltas bottom-up + ## if complement is True, then computes deltas top-down + eps <- 1e-4 + index$Ngrid <- 301 ## for gamma computations we need all the precision we can get + ## we build a lit of 4 indices with various shocks + index.list <- lapply(c(0, eps, -eps, 2*eps), function(x){ + if(x==0){ + return(index) + }else{ + newindex <- index + newindex$portfolio <- tweakportfolio(newindex$portfolio, x) + newindex$defaultprob <- 1 - SPmatrix(newindex$portfolio, length(index$cs$dates)) + return(newindex) + } + }) + bp <- matrix(0, length(index$K)-1, length(index.list)) + indexbp <- rep(0, length(index.list)) + for(j in seq_along(index.list)){ + indexbp[j] <- BCindex.pv(index.list[[j]])$bp + dist <- MFdist(index.list[[j]]) + bp[,j] <- BCtranche.pv(index.list[[j]], dist) + } + + deltas <- (bp[,2]-bp[,3])/(indexbp[2]-indexbp[3])*tranche.factor(index)/index$factor + deltasplus <- (bp[,4]-bp[,1])/(indexbp[4]-indexbp[1])*tranche.factor(index)/index$factor + gammas <- (deltasplus-deltas)/(indexbp[2]-indexbp[1])/100 + + return( list(deltas=deltas, gammas=gammas) ) +} + +BCtranche.corr01 <- function(index, eps=0.01, complement=FALSE){ + ##does a parallel shift of the skew and computes the change in pv + before <- BCtranche.pv(index, complement=complement) + index$rho[-1] <- index$rho[-1]+eps + after <- BCtranche.pv(index, complement=complement) + return(after$bp-before$bp) +} + +EL <- function(index, discounted=TRUE, shortened=0){ + ## computes the expected loss of a portfolio (time discounted if discounted is TRUE) + ## given the default curves and recovery + ## should be very close to the protection leg of the portfolio of cds + ## index should be a list with issuerweights, recov, defaultprob and cs parameters + ## shortened: number of quarters to shorten the maturity by + Ncol <- ncol(index$defaultprob)-shortened + ELvec <- as.numeric(crossprod(index$issuerweights * (1-index$recov), index$defaultprob[,1:Ncol])) + if(!discounted){ + return( ELvec[length(ELvec)] ) + }else{ + return( sum(index$cs$df[1:Ncol]*diff(c(0, ELvec))) ) + } +} + +BCindex.pv <- function(index, discounted=TRUE, shortened=0){ + Ncol <- ncol(index$defaultprob)-shortened + ELvec <- as.numeric(crossprod(index$issuerweights * (1-index$recov), index$defaultprob[,1:Ncol])) + size <- 1-as.numeric(crossprod(index$issuerweights, index$defaultprob[,1:Ncol])) + sizeadj <- 0.5*(c(1, size[-length(size)])+size) + if(!discounted){ + pl <- -ELvec[length(ELvec)] + cl <- as.numeric(crossprod(index$cs$coupons[1:Ncol], sizeadj)) + bp <- 1+cl+pl + }else{ + pl <- -sum(index$cs$df[1:Ncol]* diff(c(0, ELvec))) + cl <- as.numeric(crossprod(index$cs$coupons[1:Ncol], sizeadj * index$cs$df[1:Ncol] )) + bp <- 1+cl+pl + } + bp <- 1+cl*index$quotes$spread+pl + return(list(pl=pl, cl=cl, bp=bp)) +} + +ELtrunc <- function(index, K, rho){ + ## computes the expected loss of a portfolio below strike K + ## could be written faster by using a truncated version of lossdist + ## index should be a list with issuerweights, recov, defaultprob and cs parameters + dist <- BClossdistC(index$defaultprob, index$issuerweights, index$recov, rho, index$Z, index$w, index$N) + return( -tranche.pl(dist$L, index$cs, 0, K)) +} + +Probtrunc <- function(index, K, rho){ + dist <- BClossdistC(index$defaultprob, index$issuerweights, index$recov, rho, index$Z, index$w, index$N) + p <- cumsum(dist$L[,ncol(dist$L)]) + support <- seq(0, 1, length=index$N) + probfun <- splinefun(support, p, method="hyman") + return(probfun(K)) +} + + +BCstrikes <- function(index, K, rho) { + ## computes the strikes as a percentage of expected loss + ## Kmodified is the current attachment points (adjusted for losses) + el <- EL(index) + ELvec <- c() + for(i in 2:length(K)){ + ELvec <- c(ELvec, ELtrunc(index, K[i], rho[i])) + } + return( ELvec/el ) +} + +tranche.factor <- function(index){ + ## compute the factor to convert from delta on current notional to delta on original notional + ## K1 and K2 original strikes + return( diff(index$K)/diff(index$K.orig)*index$factor ) +} + +MFupdate.prob <- function(Z, w, rho, defaultprob, useC = TRUE){ + ## update the probabilites based on a non gaussian factor + ## distribution so that the pv of the cds stays the same. + p <- matrix(0, nrow(defaultprob), ncol(defaultprob)) + fit.prob <- if(useC) fit.probC else fit.prob + for(i in 1:nrow(defaultprob)){ + for(j in 1:ncol(defaultprob)){ + p[i,j] <- fit.prob(Z, w, rho[i], defaultprob[i,j]) + } + } + return( p ) +} + +MFlossrecovdist.prepay <- function(w, Z, rho, defaultprob, defaultprobmod, prepayprob, prepayprobmod, + issuerweights, recov, Ngrid=2*length(issuerweights)+1, defaultflag=FALSE){ + ## computes the loss and recovery distribution using the modified factor distribution + n.credit <- length(issuerweights) + n.int <- length(w) + Rstoch <- array(0, dim=c(n.int, n.credit, ncol(defaultprob))) + for(t in 1:ncol(defaultprob)){ + for(i in 1:n.credit){ + Rstoch[,i,t] <- stochasticrecov(recov[i], 0, Z, w, rho, defaultprob[i,t], defaultprobmod[i,t]) + } + } + parf <- function(i){ + dpshocked <- apply(defaultprobmod, 2, shockprob, rho=rho, Z=Z[i]) + ppshocked <- apply(prepayprobmod, 2, shockprob, rho=rho, Z=-Z[i]) + S <- 1 - Rstoch[i,,] + dist <- lossrecovdist.term(dpshocked, ppshocked, issuerweights, S, Ngrid, defaultflag) + } + L <- matrix(0, Ngrid, ncol(defaultprob)) + R <- matrix(0, Ngrid, ncol(defaultprob)) + for(i in 1:length(w)){ + dist <- parf(i) + L <- L + dist$L * w[i] + R <- R + dist$R * w[i] + } + return( list(L=L, R=R) ) +} + +MFlossdist.joint <- function(cl, w, Z, rho, defaultprob, defaultprobmod, issuerweights, recov, + Ngrid=2*length(issuerweights)+1, defaultflag=FALSE){ + ## rowSums(Q) is the default/loss distribution depending if + ## defaultflag is TRUE or FALSE (default setting is FALSE) + ## colSums(Q) is the recovery distribution + ## so that recovery is the y axis and L/D is the x axis + ## if we use the persp function, losses is the axis facing us, + ## and R is the axis going away from us. + n.credit <- length(issuerweights) + n.int <- lenth(w) + Rstoch <- array(0, dim=c(n.int, n.credit, ncol(defaultprob))) + for(t in 1:ncol(defaultprob)){ + for(i in 1:n.credit){ + Rstoch[,i,t] <- stochasticrecov(recov[i], 0, Z, w, rho, defaultprob[i,t], defaultprobmod[i,t]) + } + } + parf <- function(i){ + pshocked <- apply(defaultprobmod, 2, shockprob, rho=rho, Z=Z[i]) + S <- 1 - Rstoch[i,,] + dist <- lossrecovdist.joint.term(pshocked, 0, issuerweights, S, Ngrid, defaultflag) + gc() + return(dist) + } + temp <- parSapply(cl, 1:length(w), parf) + clusterCall(cl, gc) + Q <- array(0, dim=c(ncol(defaultprob), Ngrid, Ngrid)) + for(i in 1:length(w)){ + Q <- Q + w[i]*array(temp[,i], dim=c(ncol(defaultprob), Ngrid, Ngrid)) + } + return( Q ) +} + +MFlossdist.prepay.joint <- function(w, Z, rho, defaultprob, defaultprobmod, + prepayprob, prepayprobmod, issuerweights, recov, + Ngrid=2*length(issuerweights)+1, defaultflag=FALSE){ + ## rowSums is the loss distribution + ## colSums is the recovery distribution + ## so that recovery is the y axis and L is the x axis + ## if we use the persp function, losses is the axis facing us, + ## and R is the axis going away from us. + n.credit <- length(issuerweights) + n.int <- length(w) + Rstoch <- array(0, dim=c(n.credit, n.int, ncol(defaultprob))) + + for(t in 1:ncol(defaultprob)){ + for(i in 1:n.credit){ + Rstoch[i,,t] <- stochasticrecovC(recov[i], 0, Z, w, rho[i], + defaultprob[i,t], defaultprobmod[i,t]) + } + } + Q <- array(0, dim=c(ncol(defaultprob), Ngrid, Ngrid)) + for(t in 1:ncol(defaultprob)){ + S <- 1 - Rstoch[,,t] + Q[t,,] <- lossdistC.prepay.jointZ(defaultprobmod[,t], prepayprobmod[,t], + issuerweights, S, Ngrid, defaultflag, rho, Z, w) + } + return( Q ) +} + +MFrecovery <- function(index, defaultprobmod){ + n.credit <- length(index$issuerweights) + n.int <- length(index$Z) + Rstoch <- array(0, dim=c(n.credit, n.int, ncol(index$defaultprob))) + rho <- rep(0.45, n.credit) + for(t in 1:ncol(index$defaultprob)){ + for(i in 1:n.credit){ + Rstoch[i,,t] <- stochasticrecovC(index$recov[i], 0, index$Z, index$w.mod, + rho[i], index$defaultprob[i,t], defaultprobmod[i,t]) + } + } + return( Rstoch ) +} + +MFlossdist <- function(index){ + n.credit <- length(index$issuerweights) + rho <- rep(0.45, n.credit) + defaultprobmod <- MFupdate.probC(index$Z, index$w.mod, rho, index$defaultprob) + n.credit <- length(index$issuerweights) + n.int <- length(index$Z) + Rstoch <- MFrecovery(index, defaultprobmod) + Lw <- matrix(0, index$N, n.int) + Rw <- matrix(0, index$N, n.int) + L <- matrix(0, index$N, ncol(index$defaultprob)) + R <- matrix(0, index$N, ncol(index$defaultprob)) + for(t in 1:ncol(index$defaultprob)){ + S <- 1 - Rstoch[,,t] + Lw <- lossdistCZ(defaultprobmod[,t], index$issuerweights, S, index$N, 0, rho, index$Z) + Rw <- lossdistCZ(defaultprobmod[,t], index$issuerweights, 1-S, Ngrid, 0, rho, index$Z) + L[,t] <- Lw%*%index$w.mod + R[,t] <- Rw%*%index$w.mod + } + return(list(L=L, R=R)) +} |