From 8966b9753f876250182ef75a8c4294d523f2126a Mon Sep 17 00:00:00 2001
From: Guillaume Horel <guillaume.horel@serenitascapital.com>
Date: Thu, 8 May 2014 13:40:58 -0400
Subject: implement skew mapping by doing tranche loss proportion mapping

---
 R/tranche_functions.R | 49 +++++++++++++++++++++++++++++++++++++++++++------
 1 file changed, 43 insertions(+), 6 deletions(-)

(limited to 'R/tranche_functions.R')

diff --git a/R/tranche_functions.R b/R/tranche_functions.R
index b575c33..a5b99d8 100644
--- a/R/tranche_functions.R
+++ b/R/tranche_functions.R
@@ -651,7 +651,7 @@ BClossdistC <- function(defaultprob, issuerweights, recov, rho, Z, w,
     L <- matrix(0, N, dim(defaultprob)[2])
     R <- matrix(0, N, dim(defaultprob)[2])
     rho <- rep(rho, length(issuerweights))
-    r <- .C("BClossdist", defaultprob, dim(defaultprob)[1], dim(defaultprob)[2],
+    r <- .C("BCloss_recov_dist", defaultprob, dim(defaultprob)[1], dim(defaultprob)[2],
             as.double(issuerweights), as.double(recov), as.double(Z), as.double(w),
             as.integer(length(Z)), as.double(rho), as.integer(N), as.logical(defaultflag), L=L, R=R)
     return(list(L=r$L,R=r$R))
@@ -718,17 +718,54 @@ BCtranche.delta <- function(portfolio, index, coupon, K1, K2, rho1, rho2, Z, w,
     return( delta )
 }
 
-BCstrikes <- function(defaultprob, issuerweights, recov, cs, Kmodified, rho, Z, w, N=101) {
+EL <- function(defaultprob, issuerweights, recov, cs){
+    ## computes the expected loss of a portfolio (time discounted)
+    ## given the default curves and recovery
+    ## should be very close to the protection leg of the portfolio of cds
+    ELvec <- crossprod( issuerweights * (1-recov), defaultprob)
+    return( crossprod(cs$df, diff(c(0, as.numeric(ELvec)))) )
+}
+
+ELtrunc <- function(defaultprob, issuerweights, recov, cs, K, rho, Z, w, Ngrid){
+    ## computes the expected loss of a portfolio below strike K
+    ## could be written faster by using a truncated version of lossdistrib
+    rho <- rep(rho, length(issuerweights))
+    L <- matrix(0, Ngrid, dim(defaultprob)[2])
+    L <- .C("BCloss_dist", defaultprob, dim(defaultprob)[1], dim(defaultprob)[2],
+            as.double(issuerweights), as.double(recov), as.double(Z), as.double(w),
+            as.integer(length(Z)), as.double(rho), as.integer(Ngrid), FALSE, L=L)$L
+    support <- seq(0, 1, length=Ngrid)
+    trancheloss <- pmin(support, K)
+    ELvec <- crossprod(trancheloss, L)
+    return( crossprod(cs$df, diff(c(0, as.numeric(ELvec)))) )
+}
+
+BCstrikes <- function(defaultprob, issuerweights, recov, cs, K, rho, Z, w, N=101) {
     ## computes the strikes as a percentage of expected loss
     ## Kmodified is the current attachment points (adjusted for losses)
-    EL <- c()
+    el <- EL(defaultprob, issuerweights, recov, cs)
+    ELvec <- c()
     for(i in 2:length(K)){
-        EL <- c(EL, -BCtranche.pv(defaultprob, issuerweights, recov, cs,
-                                  K[i-1], K[i], rho[i-1], rho[i], Z, w, N)$pl)
+        ELvec <- c(ELvec, ELtrunc(defaultprob, issuerweights, recov, cs, K[i], rho[i], Z, w, N))
     }
-    return(cumsum(EL*diff(Kmodified))/sum(EL*diff(Kmodified)))
+    return( ELvec/el )
 }
 
+skewmapping <- function(defaultprob1, issuerweights1, recov1, cs1, rhofun,
+                        defaultprob2, issuerweights2, recov2, cs2, K2, Z, w, N=101) {
+
+    EL1 <- EL(defaultprob1, issuerweights1, recov1, cs1)
+    EL2 <- EL(defaultprob2, issuerweights2, recov2, cs2)
+    f <- function(x, ...){
+        return(abs(ELtrunc(defaultprob1, issuerweights1, recov1, cs1,
+                           0, x, 0, rhofun(x), Z, w, N)/EL1-
+                   ELtrunc(defaultprob2, issuerweights2, recov2, cs2,
+                           0, K2, 0, rhofun(x), Z, w, N)/EL2))
+    }
+    return( optimize(f, interval=c(0,1)) )
+}
+
+
 delta.factor <- function(K1, K2, index){
     ## compute the factor to convert from delta on current notional to delta on original notional
     ## K1 and K2 original strikes
-- 
cgit v1.2.3-70-g09d2