From be5c7fe572975dde053ded9c84b1896424ec9ca4 Mon Sep 17 00:00:00 2001
From: Guillaume Horel <guillaume.horel@serenitascapital.com>
Date: Fri, 15 Aug 2014 15:51:00 -0400
Subject: improved function for doing the skew mapping

---
 R/tranche_functions.R | 69 +++++++++++++++++++++++++++------------------------
 1 file changed, 37 insertions(+), 32 deletions(-)

(limited to 'R')

diff --git a/R/tranche_functions.R b/R/tranche_functions.R
index 83e7199..a436de9 100644
--- a/R/tranche_functions.R
+++ b/R/tranche_functions.R
@@ -691,18 +691,41 @@ BCtranche.pv <- function(index, protection=FALSE){
     return(list(pl=plvec, cl=clvec, bp=bp))
 }
 
-adjust.skew <- function(index, el1){
-    el <- EL(index)
-    K <- index$K[-c(1,length(index$K))]
-    rhofun <- splinefun(K, index$rho[-c(1, length(index$rho))], "natural")
-    K <- el/el1*K
-    return(c(0, rhofun(K), NA))
+adjust.skew <- function(index1, index2, method="ATM"){
+    #index1 is the index of 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
+    el1 <- EL(index1)
+    el2 <- EL(index2)
+    K1 <- index2$K[-c(1,length(index1$K))]
+    K2 <- index2$K[-c(1,length(index1$K))]
+    skew <- splinefun(K1, index1$rho[-c(1, length(index1$rho))], "natural")
+    aux <- function(x, index1, el1, skew, index2, el2, K2){
+        return(abs(ELtrunc(index1, x, skew(x))/el1-
+                   ELtrunc(index2, K2, skew(x))/el2))
+    }
+    if(method=="ATM"){
+        K1eq <- el1/el2 * K2
+    }else if(method == "TLP"){
+        K1eq <- c()
+        for(K2val in K2){
+            prog <- optimize(aux, interval=c(0,1),
+                             index1=index1, el1=el1, skew=skew,
+                             index2=index2, el2=el2, K2=K2val)
+            K1eq <- c(K1eq, prog$minimum)
+        }
+    }
+    return(c(0, skew(K1eq), NA))
 }
 
-theta.adjust.skew <- function(index, shortened=4){
+theta.adjust.skew <- function(index, shortened=4, method="ATM"){
     #ajust the correlation skew by doing ATM mapping on the expected loss
-    elshort <- EL(index, shortened=shortened)
-    return(adjust.skew(index, elshort))
+    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, tradedate=Sys.Date(), shortened=4){
@@ -779,19 +802,12 @@ EL <- function(index, discounted=TRUE, shortened=0){
    }
 }
 
-ELtrunc <- function(index, K, rho, Z, w, Ngrid){
+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 lossdistrib
+    ## could be written faster by using a truncated version of lossdist
     ## index should be a list with issuerweights, recov, defaultprob and cs parameters
-    rho <- rep(rho, length(index$issuerweights))
-    L <- matrix(0, Ngrid, dim(index$defaultprob)[2])
-    L <- .C("BCloss_dist", index$defaultprob, dim(index$defaultprob)[1], dim(index$defaultprob)[2],
-            as.double(index$issuerweights), as.double(index$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( sum(index$cs$df * diff(c(0, as.numeric(ELvec)))) )
+    dist <- BClossdistC(index$defaultprob, index$issuerweights, index$recov, rho, index$Z, index$w, index$N)
+    return( -tranche.pl(dist$L, index$cs, 0, K) )
 }
 
 BCstrikes <- function(defaultprob, issuerweights, recov, cs, K, rho, Z, w, N=101) {
@@ -800,22 +816,11 @@ BCstrikes <- function(defaultprob, issuerweights, recov, cs, K, rho, Z, w, N=101
     el <- EL(index)
     ELvec <- c()
     for(i in 2:length(K)){
-        ELvec <- c(ELvec, ELtrunc(index, K[i], rho[i], Z, w, N))
+        ELvec <- c(ELvec, ELtrunc(index, K[i], rho[i]))
     }
     return( ELvec/el )
 }
 
-skewmapping <- function(index1, rhofun, index2, K2, Z, w, N=101) {
-    EL1 <- EL(index1)
-    EL2 <- EL(index2)
-    f <- function(x, ...){
-        return(abs(ELtrunc(index1, x, rhofun(x), Z, w, N)/EL1-
-                   ELtrunc(index2, K2, rhofun(x), Z, w, N)/EL2))
-    }
-    return( optimize(f, interval=c(0,1), index1, rhofun, index2, K2, Z, w, N) )
-}
-
-
 tranche.factor <- function(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