aboutsummaryrefslogtreecommitdiffstats
path: root/R/cds_functions.R
diff options
context:
space:
mode:
Diffstat (limited to 'R/cds_functions.R')
-rw-r--r--R/cds_functions.R343
1 files changed, 343 insertions, 0 deletions
diff --git a/R/cds_functions.R b/R/cds_functions.R
new file mode 100644
index 00000000..7404328d
--- /dev/null
+++ b/R/cds_functions.R
@@ -0,0 +1,343 @@
+source("cds_utils.R")
+
+setClass("defaultcurve", representation(dates="Date", hazardrates="numeric"))
+setClass("defaultprepaycurve", representation(prepayrates="numeric"), contains="defaultcurve")
+setClass("creditcurve", representation(issuer="character", startdate="Date",
+ recovery="numeric", curve="defaultprepaycurve"))
+
+survivalProbability1 <- function(startdate, date, survival.curve) {
+ #based on a flat hazard rate curve
+ T <- yearFrac(startdate, survival.curve$dates)
+ Tmat <- yearFrac(startdate, date)
+ for ( i in 1:length(survival.curve$dates) ){
+ if ( date >= survival.curve$dates[i] ) {
+ next
+ }else{
+ if( i > 1 ) {
+ w <- ( Tmat - T[i-1] ) / (T[i] - T[i-1])
+ logprob <- - (1-w) * survival.curve$hazardrates[i-1] * T[i-1] -
+ w * survival.curve$hazardrates[i] * T[i]
+ }else{
+ logprob <- - Tmat * survival.curve$hazardrates[1]
+ return( exp(as.numeric(logprob)) )
+ }
+ }
+ }
+ ## if date is greater than last survival.curve date, keep the hazard rate flate
+ logprob <- - yearFrac(startdate, date) * survival.curve$hazardrates[i]
+ return( exp(as.numeric(logprob)) )
+}
+
+survivalProbability.exact <- function(credit.curve, date) {
+ #based on a forward hazard rate curve
+ curve <- credit.curve@curve
+ T <- c(0, yearFrac(credit.curve@startdate, curve@dates))
+ dT <- diff(T)
+ Tmat <- yearFrac(credit.curve@startdate, date)
+ logprob <- 0
+ for ( i in 1:length(dT) ){
+ if ( date > curve@dates[i] ) {
+ logprob <- logprob - curve@hazardrates[i] * dT[i]
+ }else{
+ if( i > 1 ){
+ logprob <- logprob - curve@hazardrates[i] * (Tmat - T[i])
+ }else{
+ logprob <- logprob - curve@hazardrates[1] * Tmat
+ }
+ break
+ }
+ }
+ return( exp(as.numeric(logprob)) )
+}
+
+PD <- function(sc){
+ ## computes the default probability associated with the survival curve
+ if(class(sc)=="defaultprepaycurve"){
+ T <- c(0, yearFrac(today(), sc@dates))
+ dT <- diff(T)
+ return( 1-cumprod(exp(-sc@hazardrates * dT)) )
+ }else{
+ stop("not a default curve")
+ }
+}
+
+couponleg <- function(cs, sc, accruedondefault=TRUE){
+ ## computes the pv of the risky coupon leg based on a given coupon schedule
+ ## and a survival curve. Also called premium leg or fixed leg.
+ dT <- diff(c(0, yearFrac(today(), cs$dates)))
+ PD <- cumprod(exp(-sc@hazardrates[1:length(dT)] * dT))
+ if(accruedondefault){
+ PDadj <- 0.5 * (c(1, PD[-length(PD)]) + PD)
+ }else{
+ PDadj <- PD
+ }
+ return( crossprod(PDadj, cs$coupons * cs$df) )
+}
+
+couponleg.flat <- function(cs, h, accruedondefault=TRUE){
+ T <- yearFrac(today(), cs$dates)
+ PD <- exp(- h * T )
+ if(accruedondefault){
+ PDadj <- 0.5 * (c(1, PD[-length(PD)]) + PD)
+ }else{
+ PDadj <- PD
+ }
+ return( crossprod(PDadj, cs$coupons * cs$df) )
+}
+
+couponleg.prepay <- function(cs, dpc, accruedondefault=TRUE){
+ ## computes the pv of the risky coupon leg from the coupon schedule,
+ ## a hazard rate curve, and a prepay curve. We assume the poisson process driving
+ ## default and prepayment are independent, so the intensity of either event
+ ## happenning is the sum of the two.
+
+ dT <- diff(c(0, yearFrac(today(), cs$dates)))
+ SP <- cumprod(exp( - (dpc@hazardrates[1:length(dT)] + dpc@prepayrates[1:length(dT)] ) * dT))
+ if(accruedondefault){
+ SPadj <- 0.5 * (c(1, SP[-length(SP)]) + SP)
+ }else{
+ SPadj <- SP
+ }
+ return( crossprod(SPadj, cs$coupons * cs$df) )
+}
+
+couponleg.prepay.flat <- function(cs, h){
+ dpc <- new("defaultprepaycurve", dates = cs$dates, hazardrates=rep(h, length(cs$dates)),
+ prepayrates=rep(k(h), length(cs$dates)))
+ couponleg.prepay(cs, dpc)
+}
+
+dcouponleg <- function(cs, sc, index, accruedondefault=TRUE){
+ ## derivative of couponleg with respect to hazardrate
+ dT <- diff(c(0, yearFrac(today(), cs$dates)))
+ PD <- cumprod(exp(-sc@hazardrates[1:length(dT)] * dT))
+ dPD <- cumsum(index * dT) * PD
+ if(accruedondefault){
+ dPDadj <- 0.5 * (c(0, dPD[-length(PD)]) + dPD)
+ }else{
+ dPDadj <- dPD
+ }
+ return( - crossprod( dPDadj, cs$coupons * cs$df) )
+}
+
+k <- function(h, alpha=0.25, beta=15){
+ ## prepay rate as a function of the hazardrate
+ ## this is a decreasing function
+ ## alpha is the maximum prepay rate
+ return ( alpha * exp(- beta * h) )
+}
+
+dcouponleg.prepay <- function(cs, dpc, index, beta, accruedondefault=TRUE){
+ ## derivative of couponleg.prepay with respect to hazardrate
+ ## If k is the prepay rate, it assumes dk/dh = - beta * k,
+ ## which is the case if k(h) = alpha * exp(-beta *h)
+ dT <- diff(c(0, yearFrac(today(), cs$dates)))
+ SP <- cumprod(exp( - (dpc@hazardrates[1:length(dT)] + dpc@prepayrates[1:length(dT)] ) * dT))
+ dSP <- -cumsum(index * dT * (1 - beta * dpc@prepayrates[1:length(dT)])) * SP
+ if(accruedondefault){
+ dSPadj <- 0.5 * (c(0, dSP[-length(SP)]) + dSP)
+ }else{
+ SPadj <- dSP
+ }
+ return( crossprod(dSPadj, cs$coupons * cs$df) )
+}
+
+dcouponleg.prepay.flat <- function(cs, h, index, beta, accruedondefault=TRUE){
+ hvec <- rep(h, length(cs$dates))
+ kvec <- k(h)
+ dpc <- new("defaultprepaycurve", dates=cs$dates, hazardrates=hvec
+ prepayrates=kvec )
+ return( dcouponleg.prepay(cs, dpc, index, beta, accruedondefault))
+}
+
+cdsduration <- function(sc, maturity, accruedondefault=TRUE){
+ # computes the risky PV01, also called risky annuity of a cds
+ cs <- couponSchedule(nextIMMDate(today()), maturity, "Q", "FIXED", 1)
+ couponleg(cs, sc, accruedondefault)
+}
+
+defaultleg.flat <- function(cs, h, recovery){
+ ## Computes the pv of the default leg of a cds based on a given
+ ## coupon schedule, flat hazard rate, and recovery.
+ T <- yearFrac(today(), cs$dates)
+ Q <- exp(-h * T) * cs$df
+ Qmid <- 1/2 * (c(1, Q[-length(Q)]) + Q)
+ r <- (1 - recovery) * crossprod(h * Qmid, dT)
+ return( r )
+}
+
+ddefaultleg.flat <- function(cs, h, recovery){
+ ## derivative of defaultleg.flat with respect to hazardrate
+ T <- yearFrac(today(), cs$dates)
+ dT <- diff(c(0, T))
+ dQ <- - T * exp(-h * T) * cs$df
+ Q <- exp(-h * T) * cs$df
+ Qmid <- 1/2 * (c(1, Q[-length(Q)]) + Q)
+ dQmid <- 1/2 * (c(0, dQ[-length(dQ)]) + dQ)
+ dr <- (1-recovery) * (crossprod(Qmid, dT) + h * crossprod(dQmid, dT))
+ return( dr )
+}
+
+defaultleg <- function(cs, sc, recovery){
+ ## Computes the pv of the default leg of a cds based on a given
+ ## coupon schedule, hazard rate curve, and recovery.
+ T <- yearFrac(today(), cs$dates)
+ dT <- diff(c(0, T))
+ Q <- cumprod(exp(-sc@hazardrates[1:length(dT)] * dT)) * cs$df
+ Qmid <- 1/2 * (c(1, Q[-length(Q)]) + Q)
+ r <- (1 - recovery) * crossprod(sc@hazardrates[1:length(dT)] * Qmid, dT)
+ return( r )
+}
+
+ddefaultleg <- function(cs, sc, recovery, index){
+ ## derivative of defaultleg with respect to hazardrate
+ T <- yearFrac(today(), cs$dates)
+ dT <- diff(c(0,T))
+ Q <- cumprod(exp(-sc@hazardrates[1:length(dT)] * dT)) * cs$df
+ dQ <- - cumsum(index * dT) * Q
+ Qmid <- 1/2 * (c(1, Q[-length(Q)]) + Q)
+ dQmid <- 1/2 *(c(0, dQ[-length(dQ)]) + dQ)
+ dr <- (1-recovery) * (crossprod(index * Qmid, dT) + crossprod(sc@hazardrates * dQmid, dT))
+ return( dr )
+}
+
+defaultleg.prepay <- function(cs, sc, recovery){
+ ## Computes the pv of the default leg of a cds based on a given
+ ## coupon schedule, hazard rates curve, prepay curves, and recovery.
+ T <- yearFrac(today(), cs$dates)
+ dT <- diff(c(0, T))
+ Q <- cumprod(exp(- (dpc@hazardrates[1:length(dT)]+dpc@prepayrates[1:length(dT)]) * dT)) * cs$df
+ Qmid <- 1/2 * (c(1, Q[-length(Q)]) + Q)
+ r <- (1 - recovery) * crossprod(dpc@hazardrates[1:length(dT)] * Qmid, dT)
+ +crossprod(dpc@prepayrates[1:length(dT)] * Qmid, dT)
+ return( r )
+}
+
+cdspv <- function(cs, sc, R){
+ return ( couponleg(cs, sc) - defaultleg(cs, sc, R) )
+}
+
+cdsspread <- function(sc, maturity, recovery){
+ ## computes exact cds running spread for a cds
+ ## should be very close to hazardrate * (1-recovery)
+ cs <- couponSchedule(nextIMMDate(today()), maturity, "Q", "FIXED", 1)
+ defaultleg(cs, sc, recovery)/couponleg(cs, sc)
+}
+
+dcdspv <- function(cs, sc, recovery, index){
+ ## derivative of cdspv with respect to hazardrate
+ return ( dcouponleg(cs, sc, index) - ddefaultleg(cs, sc, recovery, index) )
+}
+
+cdshazardrate.flat <- function(upfront, running, maturity, R=0.4){
+ ## computes the implied hazard rate of the cds based on the upfront
+ ## and running quotes, as well as maturity and recovery
+ cs <- couponSchedule(nextIMMDate(today()), maturity, "Q", "FIXED", running)
+ lambda <- 0.05
+ eps <- 1e-8
+ while(abs(cdspv(lambda, cs, R) + upfront) > eps){
+ lambda <- lambda - (upfront + cdspv(lambda, cs, R))/dcdspv(lambda, cs, R)
+ }
+ return(lambda)
+}
+
+cdshazardrate <- function(quotes, R=0.4){
+ ## bootstrap the implied hazard rate curve of the cds based on the upfront
+ ## and running quotes, as well as maturity and recovery
+ previous.maturity <- today()
+ hvec <- c()
+ previous.hvec <- c()
+ eps <- 1e-8
+ previous.cs <- data.frame()
+ for(i in 1:nrow(quotes)){
+ maturity <- quotes$maturity[i]
+ cs <- couponSchedule(nextIMMDate(today()), maturity, "Q", "FIXED", quotes$running[i])
+ new.h <- 0.05
+ flength <- nrow(cs) - nrow(previous.cs)
+ hvec <- c(previous.hvec, rep(new.h, flength))
+ sc <- new("defaultprepaycurve", dates=cs$dates, hazardrates=hvec)
+ index <- c(rep(0, length(previous.hvec)), rep(1, flength))
+
+ while(abs(cdspv(cs, sc, R) + quotes$upfront[i]) > eps){
+ new.h <- new.h - (quotes$upfront[i] + cdspv(cs, sc, R))/dcdspv(cs, sc, R, index)
+ hvec <- c(previous.hvec, rep(new.h, flength))
+ sc@hazardrates <- hvec
+ }
+ previous.hvec <- hvec
+ previous.maturity <- maturity
+ previous.cs <- cs
+ }
+ return(sc)
+}
+
+
+bonddp <- function(collateral, R=0.7){
+ ## bootstrap the implied hazard rate curve of a bond based on the upfront
+ ## and running quotes, as well as maturity and recovery
+ previous.maturity <- today()
+ hvec <- c()
+ previous.hvec <- c()
+ eps <- 1e-8
+ previous.cs <- data.frame()
+ for(i in 1:nrow(quotes)){
+ maturity <- quotes$maturity[i]
+ cs <- couponSchedule(collateral$nextpaydate, collateral$maturity,
+ collateral$frequency, collateral$fixedorfloat,
+ collateral$grosscoupon*0.01, collateral$spread*0.01)
+ new.h <- 0.05
+ flength <- nrow(cs)-nrow(previous.cs)
+ hvec <- c(previous.hvec, rep(new.h, flength))
+ sc <- new("survivalcurve", dates=cs$dates, hazardrates=hvec)
+ index <- c(rep(0, length(previous.hvec)), rep(1, flength))
+
+ while(abs(cdspv(cs, sc, R) + quotes$upfront[i]) > eps){
+ new.h <- new.h - (quotes$upfront[i] + cdspv(cs, sc, R))/dcdspv(cs, sc, R, index)
+ hvec <- c(previous.hvec, rep(new.h, flength))
+ sc@hazardrates <- hvec
+ }
+ previous.hvec <- hvec
+ previous.maturity <- maturity
+ previous.cs <- cs
+ }
+ return(sc)
+}
+
+indexpv <- function(portfolio, index, epsilon=0){
+ ## computes the intrinsic index pv of a portfolio of cds
+ r <- rep(0, length(portfolio))
+ cs <- couponSchedule(nextIMMDate(today()), index$maturity, "Q", "FIXED", index$coupon)
+ for(i in 1:length(portfolio)){
+ if(epsilon!=0){
+ tweakedcurve <- portfolio[[i]]@curve
+ tweakedcurve@hazardrates <- tweakedcurve@hazardrates * (1 + epsilon)
+ r[i] <- cdspv(cs, tweakedcurve, index$recovery)
+ }else{
+ r[i] <- cdspv(cs, portfolio[[i]]@curve, index$recovery)
+ }
+ }
+ return( 1+mean(r) )
+}
+
+tweakcurves <- function(portfolio, index){
+ ## computes the tweaking factor
+ epsilon <- 0
+ f <- function(epsilon, ...){
+ abs(indexpv(portfolio, index, epsilon)-index$indexref)
+ }
+ epsilon <- optimize(f, c(-0.5, 0.5), portfolio, index, tol=1e-6)$minimum
+ portfolio.new <- portfolio
+ for(i in 1:length(portfolio)){
+ portfolio.new[[i]]@curve@hazardrates <- portfolio[[i]]@curve@hazardrates * (1 + epsilon)
+ }
+ return( portfolio.new )
+}
+
+SPmatrix <- function(portfolio, index){
+ cs <- couponSchedule(nextIMMDate(today()), index$maturity, "Q", "FIXED", index$coupon)
+ SP <- matrix(0, length(portfolio), length(cs$dates))
+ for(i in 1:length(portfolio)){
+ SP[i,] <- 1 - PD(portfolio[[i]]@curve)[1:length(cs$dates)]
+ }
+ return( SP )
+}
+