diff options
Diffstat (limited to 'R/calibrate_tranches.R')
| -rw-r--r-- | R/calibrate_tranches.R | 274 |
1 files changed, 274 insertions, 0 deletions
diff --git a/R/calibrate_tranches.R b/R/calibrate_tranches.R new file mode 100644 index 00000000..3e2599fc --- /dev/null +++ b/R/calibrate_tranches.R @@ -0,0 +1,274 @@ +library("parallel")
+setwd("//WDSENTINEL/share/CorpCDOs/R")
+source("cds_utils.R")
+source("cds_functions_generic.R")
+source("index_definitions.R")
+source("tranche_functions.R")
+source("yieldcurve.R")
+source("optimization.R")
+
+cl <- makeCluster(6)
+
+MarkitData <- getMarkitIRData()
+L1m <- buildMarkitYC(MarkitData, dt = 1/12)
+L2m <- buildMarkitYC(MarkitData, dt = 1/6)
+L3m <- buildMarkitYC(MarkitData)
+L6m <- buildMarkitYC(MarkitData, dt = 1/2)
+setEvaluationDate(as.Date(MarkitData$effectiveasof))
+
+## calibrate HY19
+## calibrate the single names curves
+singlenames.data <- read.table(file="clipboard", sep="\t", header=T)
+nondefaulted <- singlenames.data[!singlenames.data$ticker %in% hy17$defaulted,]
+bps <- 1e-4
+
+cdsdates <- as.Date(character(0))
+for(tenor in paste0(1:5, "y")){
+ cdsdates <- c(cdsdates, cdsMaturity(tenor))
+}
+
+## clusterEvalQ(cl, {setClass("abstractcurve")
+## setClass("defaultcurve", contains="abstractcurve",
+## representation(dates="Date", hazardrates="numeric"))
+## setClass("creditcurve", representation(issuer="character", startdate="Date",
+## recovery="numeric", curve="defaultcurve"))})
+## clusterExport(cl, list("nondefaulted", "cdsdates", "cdshazardrate", "today",
+## "bps", "couponSchedule", "nextIMMDate", "DiscountCurve", "L3m"))
+## test <- parSapply(cl, 1:nrow(nondefaulted), parf)
+
+## parf <- function(i){
+## SC <- new("creditcurve",
+## recovery=nondefaulted$recovery[i]/100,
+## startdate=today(),
+## issuer=as.character(nondefaulted$ticker[i]))
+## quotes <- data.frame(maturity=cdsdates, upfront = as.numeric(nondefaulted[i,5:9])*0.01,
+## running = rep(nondefaulted$running[i]*bps,5))
+## return( cdshazardrate(quotes, nondefaulted$recovery[i]/100))
+## }
+
+hy19portfolio <- c()
+for(i in 1:nrow(nondefaulted)){
+ SC <- new("creditcurve",
+ recovery=nondefaulted$recovery[i]/100,
+ startdate=today(),
+ issuer=as.character(nondefaulted$ticker[i]))
+ quotes <- data.frame(maturity=cdsdates, upfront = as.numeric(nondefaulted[i,4:8])*0.01,
+ running=rep(nondefaulted$running[i]*bps, 5))
+ SC@curve <- cdshazardrate(quotes, nondefaulted$recovery[i]/100)
+ hy19portfolio <- c(hy19portfolio, SC)
+}
+
+issuerweights <- rep(1/length(hy19portfolio), length(hy19portfolio))
+hy19$indexref <- 0.99
+hy19portfolio.tweaked <- tweakcurves(hy19portfolio, hy19)
+
+SurvProb <- SPmatrix(hy19portfolio.tweaked, hy19)
+## load common parameters
+K <- c(0, 0.15, 0.25, 0.35, 1)
+Kmodified <- adjust.attachments(K, hy19$loss, hy19$factor)
+tranche.upf <- c(37.5625, 87.25, 101.8125, 115.0625)
+tranche.running <- c(0.05, 0.05, 0.05, 0.05)
+
+Ngrid <- 2*nrow(nondefaulted)+1
+recov <- sapply(hy19portfolio.tweaked, attr, "recovery")
+cs <- couponSchedule(nextIMMDate(today()), hy19$maturity,"Q", "FIXED", 0.05, 0)
+
+## calibrate the tranches using base correlation
+rhovec <- c()
+f <- function(rho, ...){
+ temp <- BClossdistC(SurvProb, issuerweights, recov, rho, Ngrid)
+ bp <- 100*(1+1/(Kmodified[i]-Kmodified[i-1]) *
+ (tranche.pv(temp$L, temp$R, cs, 0, Kmodified[i], Ngrid) -
+ tranche.pv(oldtemp$L, oldtemp$R, cs, 0, Kmodified[i-1], Ngrid)))
+ return( abs(tranche.upf[i-1]-bp))
+}
+
+for(i in 2:length(Kmodified)){
+ rho <- optimize(f, interval=c(0,1),
+ SurvProb, issuerweights, recov, Ngrid, tranche.upf, Kmodified, cs, oldtemp)$minimum
+ oldtemp <- BClossdistC(SurvProb, issuerweights, recov, rho, Ngrid)
+ rhovec <- c(rhovec, rho)
+}
+
+rhovec <- c(0, rhovec)
+deltas <- c()
+for(i in 2:5){
+ deltas <- c(deltas, BCtranche.delta(hy19portfolio.tweaked, hy19, 0.05, K[i-1], K[i], rhovec[i-1], rhovec[i], Ngrid))
+}
+
+##calibrate by modifying the factor distribution
+bottomup <- 1:3
+topdown <- 2:4
+n.int <- 500
+n.credit <- length(hy19portfolio)
+errvec <- c()
+quadrature <- gauss.quad.prob(n.int, "normal")
+w <- quadrature$weights
+Z <- quadrature$nodes
+w.mod <- w
+defaultprob <- 1 - SurvProb
+p <- defaultprob
+rho <- 0.45
+
+clusterExport(cl, list("shockprob", "issuerweights", "rho", "Z", "lossrecovdist.term",
+ "lossrecovdist", "lossdistC", "Ngrid",
+ "tranche.pvvec", "tranche.pv", "tranche.pl", "tranche.cl",
+ "trancheloss", "trancherecov", "pos", "Kmodified", "cs"))
+## TODO: investigate if this is the right thing w.r.t recovery
+parf <- function(i){
+ pshocked <- apply(p, 2, shockprob, rho=rho, Z=Z[i])
+ S <- 1 - Rstoch[i,,]
+ dist <- lossrecovdist.term(pshocked,, issuerweights, S, Ngrid)
+ return( tranche.pvvec(Kmodified, dist$L, dist$R, cs))
+}
+
+for(l in 1:100){
+ Rstoch <- array(0, dim=c(n.int, n.credit, ncol(SurvProb)))
+ for(t in 1:ncol(SurvProb)){
+ for(i in 1:n.credit){
+ Rstoch[,i,t] <- stochasticrecov(recov[i], 0, Z, w.mod, rho, defaultprob[i,t], p[i,t])
+ }
+ }
+
+ clusterExport(cl, list("Rstoch", "p"))
+ result <- parSapply(cl, 1:n.int, parf)
+ ## solve the optimization problem
+ program <- KLfit(100*(result[bottomup,]+1), w, tranche.upf[bottomup])
+
+
+ err <- 0
+ for(i in 1:n.credit){
+ for(j in 1:ncol(p)){
+ err <- err + abs(crossprod(shockprob(p[i,j], rho, Z), program$weight) - defaultprob[i,j])
+ }
+ }
+ errvec <- c(errvec, err)
+
+ ## update the new probabilities
+ p <- MFupdate.prob(Z, program$weight, rho, defaultprob)
+
+ errvec <- c(errvec, err)
+ w.mod <- program$weight
+ cat(err,"\n")
+}
+
+write.table(data.frame(Z=Z, w=w.mod), file=paste0("calibration-", Sys.Date(), ".csv"), col.names=T, row.names=F, sep=",")
+
+## computes MFdeltas
+newportf <- hy19portfolio.tweaked
+eps <- 1e-4
+for(i in 1:length(newportf)){
+ newportf[[i]]@curve@hazardrates <- hy19portfolio.tweaked[[i]]@curve@hazardrates * (1 + eps)
+}
+SurvProb2 <- SPmatrix(newportf, hy19)
+p2 <- MFupdate.prob(Z, w.mod, rho, 1-SurvProb2)
+dPVtranches <- MFtranche.pv(cl, cs, w.mod, rho, 1-SurvProb2, p2, issuerweights) - MFtranche.pv(cl, cs, w.mod, rho, defaultprob, p, issuerweights)
+dPVindex <- indexpv(newportf, hy19)-indexpv(hy19portfolio.tweaked, hy19)
+MFdeltas <- dPVtranches/dPVindex
+
+#global deltas
+PVtranches <- MFtranche.pv(cl, cs, w.mod, rho, 1-SurvProb, p, issuerweights, Kmodified)
+pnlindex <- 1+t(PVtranches$pv.w)%*%diff(Kmodified)-hy19$indexref
+plot(1-cumsum(w.mod),PVindex, main="pnl of going long the index", xlab="Market factor cumulative probability distribution", type="l")
+pnltranches <- sweep(1+t(PVtranches$pv.w), 2, tranche.upf/100)
+matplot(1-cumsum(w.mod), pnltranches, 2, tranche.upf/100, main="pnl of going long the tranches", xlab="Market factor cumulative probability distribution", type="l")
+global.deltas <- rep(0,4)
+for(i in 1:4){
+ global.deltas[i] <- lm(pnltranches[,i]~0+pnlindex, weights=w.mod)$coef
+}
+pnlhedged.model <- pnltranches-pnlindex%*%t(global.deltas)
+pnlhedged.market <- pnltranches-pnlindex%*%t(market.deltas)
+## generate a bunch of plots
+postscript("hedged tranches global.eps")
+matplot(1-cumsum(w.mod), pnlhedged.model, type="l", ylab="pnl of the hedged package (global deltas)", xlab="Market factor cumulative probability distribution")
+legend(0.35, y=0.75, c("0-15","15-25","25-35","35-100"), col=1:4, lty=1:4)
+dev.off()
+postscript("hedged tranches market.eps")
+matplot(1-cumsum(w.mod), pnlhedged.market, type="l", ylab="Pnl", main="pnl of the hedged package (market deltas)", xlab="Market factor cumulative probability distribution")
+legend(0.15, y=1.15, c("0-15","15-25","25-35","35-100"), col=1:4, lty=1:4)
+dev.off()
+postscript("0-15 hedged.eps")
+matplot(1-cumsum(w.mod), cbind(pnlhedged.model[,1], pnlhedged.market[,1]), type="l", ylab="Pnl", main="Pnl of the hedged package (market vs global deltas)\n 0-15 tranche", xlab="Market factor cumulative probability distribution")
+dev.off()
+postscript("15-25 hedged.eps")
+matplot(1-cumsum(w.mod), cbind(pnlhedged.model[,2], pnlhedged.market[,2]), type="l", ylab="Pnl", main="Pnl of the hedged package (market vs global deltas)\n 15-25 tranche", xlab="Market factor cumulative probability distribution")
+dev.off()
+postscript("25-35 hedged.eps")
+matplot(1-cumsum(w.mod), cbind(pnlhedged.model[,3], pnlhedged.market[,3]), type="l", ylab="Pnl", main="Pnl of the hedged package (market vs global deltas)\n 25-35 tranche", xlab="Market factor cumulative probability distribution")
+dev.off()
+postscript("35-100 hedged.eps")
+matplot(1-cumsum(w.mod), cbind(pnlhedged.model[,4], pnlhedged.market[,4]), type="l", ylab="Pnl", main="Pnl of the hedged package (market vs global deltas)\n 35-100 tranche", xlab="Market factor cumulative probability distribution")
+dev.off()
+
+## scenario based pricing
+## generate the scenarios by finding the quantiles of the loss and recovery distributions
+n.scenarios <- 100
+percentiles <- (seq(0, 1, length=n.scenarios+1)[-1]+
+ seq(0, 1, length=n.scenarios+1)[-(n.scenarios+1)])/2
+l <- matrix(0, ncol(defaultprob), n.scenarios)
+r <- matrix(0, ncol(defaultprob), n.scenarios)
+
+MFdist <- MFlossdistrib(w.mod, rho, defaultprob, p, issuerweights, Ngrid)
+MFdist.orig <- MFlossdistrib(w, rho, defaultprob, defaultprob, issuerweights, Ngrid)
+lossdist.orig <- BClossdistC(SurvProb, issuerweights, recov, rho, 101, 500)
+for(i in 1:17){
+ Lfun <- splinefun(c(0, cumsum(MFdist$L[,i])),c(0, seq(0, 1, length=Ngrid)), "monoH.FC")
+ Rfun <- splinefun(c(0, cumsum(MFdist$R[,i])),c(0, seq(0, 1, length=Ngrid)), "monoH.FC")
+ for(j in 1:99){
+ l[i, j] <- Lfun(d1[j])
+ r[i, j] <- Rfun(d1[j])
+ }
+}
+
+## joint generation of scenarios for loss and recovery distribution
+clusterExport(cl, list("lossrecovdist.joint.term", "lossdistribC.joint"))
+MFdist2 <- MFlossdistrib2(cl, w.mod, rho, defaultprob, p, issuerweights, Ngrid)
+## memory never gets released by the clusters for some reasons, needs to call gc twice
+clusterCall(cl, gc)
+clusterCall(cl, gc)
+gc()
+#L+R<=1
+test <- apply(MFdist2[21,,], 2, rev)
+test.tri <- lower.tri(test)
+sum(test[test.tri])#close to 1
+
+## D=L+R
+dist <- MFdist2[21,,]
+supportD <- outer(seq(0,1, length=Ngrid), seq(0,1, length=Ngrid), "+")
+## compute the joint density of (L/D, D)
+## u=x/(x+y)
+## v=x+y
+x <- seq(0,1, length=Ngrid) %o% seq(0,1, length=Ngrid)
+y <- sweep(-x, 2, seq(0, 1, length=Ngrid), "+")
+xgrid <- round(x/0.005)
+ygrid <- round(y/0.005)
+distchv <- matrix(0, Ngrid, Ngrid)
+for(i in 1:Ngrid){
+ for(j in 1:Ngrid){
+ distchv[i,j] <- dist[xgrid[i,j]+1, ygrid[i,j]+1]*seq(0,1, length=Ngrid)[j]
+ }
+}
+
+## compute reinvesting distribution
+beta <- 1.1
+f <- function(l, r)(1-l-r)/(1-r)^beta
+test <- outer(seq(0, 1, length=Ngrid), seq(0, 1, length=Ngrid), f)
+plot(density(x=test[-Ngrid*Ngrid], weights=dist[-Ngrid*Ngrid], from=0, to=5))
+
+## two dimensional scenarios
+n.scenarios <- 100
+percentiles <- (seq(0, 1, length=n.scenarios+1)[-1]+
+ seq(0, 1, length=n.scenarios+1)[-(n.scenarios+1)])/2
+l <- matrix(0, ncol(defaultprob), n.scenarios)
+r <- matrix(0, ncol(defaultprob), n.scenarios)
+
+MFdist <- MFlossdistrib(w.mod, rho, defaultprob, p, issuerweights, Ngrid)
+i <- 21
+Lfun <- splinefun(c(0, cumsum(MFdist$L[,i])),c(0, seq(0, 1, length=Ngrid)), "monoH.FC")
+
+Rfun <- splinefun(c(0, cumsum(MFdist$R[,i])),c(0, seq(0, 1, length=Ngrid)), "monoH.FC")
+ for(j in 1:99){
+ l[i, j] <- Lfun(d1[j])
+ r[i, j] <- Rfun(d1[j])
+ }
+}
|