0
$\begingroup$

I would like to price a n-th to default swap on a basket of 5 assets or reference names. I started to code in R and I put the routine hereby. my doubt is how to use the m = {m1,m2,m3,m4,m5} series which represents the simulated defaults for each reference. Many Thanks for any tips you can give me if you have already worked for a similar project.

Many thanks for the support,

Eugenio


library(moments) 
library(psych)
library(plyr)
library(copula)
library(MASS) 
# Importing the times series of the ref. names 
GS     = read.csv("GS.csv")
BAC    = read.csv("BAC.csv")
JPM    = read.csv("JPM.csv") 
CIT    = read.csv("CIT.csv") 
WFC    = read.csv("WFC.csv")
# Calculating the weekly returns and the standardised returns   
GS_r   = diff(log(GS[,6])) 
BAC_r  = diff(log(BAC[,6]))
JPM_r  = diff(log(JPM[,6])) 
CIT_r  = diff(log(CIT[,6]))
WFC_r  = diff(log(WFC[,6]))
Z_GS   = (GS_r -mean(GS_r))/sd(GS_r) 
Z_BAC  = (BAC_r-mean(BAC_r))/sd(BAC_r)
Z_JPM  = (JPM_r-mean(JPM_r))/sd(JPM_r)
Z_CIT  = (CIT_r-mean(CIT_r))/sd(CIT_r) 
Z_WFC  = (WFC_r-mean(WFC_r))/sd(WFC_r)
TAU    = cor(cbind(Z_GS, Z_BAC, Z_JPM, Z_CIT, Z_WFC), method = "kendall")
RHO    = sin((pi/2)*TAU)
# Credit spreads for the 5 Reference names from the market   
CS_GS        = c(150, 155, 160, 170, 190)
CS_BAC       = c(20, 25, 32, 45, 57)
CS_JPM       = c(30, 45, 59, 67, 89)
CS_CIT       = c(100, 150, 155, 175, 210)
CS_WFC       = c(50, 59, 80, 120, 180)
CS           = matrix(c(CS_GS,CS_BAC,CS_JPM,CS_CIT,CS_WFC), nrow = 5, ncol = 5)

# Calculate the survival probabilities Matrix 

ISP_mat  = function(CS) { 
RR       = 0.40 
df       = c(1,1.24,1.41,1.60,1.77,1.93) # Discounted Factors curve 
nrows    = length(df) 
ncols    = nrows-1 
DF       = matrix(df,nrows,ncols) 
ISP      = matrix(0,nrows,ncols)
h        = matrix(0,nrows,ncols)
terms    = matrix(0,nrows,ncols)
term     = matrix(0,nrows,ncols)
quotient = matrix(0,nrows,ncols)
firstTerm    = matrix(0,nrows,ncols)
lastTerm     = matrix(0,nrows,ncols)
for (i in 1:nrows) { 
 for (j in 1:ncols) { 
     if (i==1) { ISP[i,j] = 1 }
     if (i==2) { ISP[i,j] = (1-RR) / ((1-RR) + CS[i-1,j]/10000) } 
     if (i> 2) {     
        for (k in 1:(i-2) ) {
          term[i,j] = DF[k,j]*((1-RR)*ISP[k,j]-((1-RR)+(CS[i-1,j]/10000))*ISP[k+1,j])
          terms[i,j] = terms[i,j] + term[i,j]
        } 
    quotient[i,j]  = (DF[i-1,j]*((1-RR)+1*CS[i-1,j]/10000))
    firstTerm[i,j] = (terms[i,j] / quotient[i,j])
    lastTerm[i,j]  = ISP[i-1,j]*(1-RR)/(1-RR+1*CS[i-1,j]/10000)
    ISP[i,j]       = (firstTerm[i,j] + lastTerm[i,j])
}
}
}
return(ISP) 
} 

surv.prob = ISP_mat(CS)

# Calculate hazard matrix and cumulative hazard rates matrix 

nrows.h = nrow(surv.prob) - 1 
ncols.h = ncol(surv.prob)
hazard  = matrix(0,nrows.h, ncols.h) 
cumulativeHazard = matrix(0, nrows.h, ncols.h) 
for (i in 1:nrows.h) { 
 for (j in 1:ncols.h) { 
   cumulativeHazard[i,j] = - log(surv.prob[i+1,j]) 
    if (i==1) { 
      hazard[i,j] = cumulativeHazard[i,j] 
    }  
    if (i>1) { 
      hazard[i,j] = cumulativeHazard[i,j] - cumulativeHazard[i-1,j]
    }  
}
}

# Calculate the Default Times Matrix 

nSim   = 10000
nrows  = nSim 
ncols  = 5 
time   = c(0,1,2,3,4,5)
df     = c(1,1.24,1.41,1.60,1.77,1.93)
e      = matrix(c(time,df), nrow = length(df), ncol = 2)
t      = matrix(0,nrows,ncols)
tenor  = matrix(0,ncols,ncols)
dt           = matrix(0,nrows,ncols) 
defaultTenor = matrix(0,nrows,ncols) 
defaultTime  = matrix(0,nrows,ncols) 
u = pnorm(matrix(rnorm(nrows*ncols), ncol = ncols) %*% chol(RHO))
for (i in 1:nrows) {
  for (j in 1:ncols) {
    t[i,j] = abs(log(1-u[i,j]))
      for (k in 1:ncols) {
             if (cumulativeHazard[k,j] >= t[i,j] ) {  
             tenor[k,j] = k-1 
             } 
             if (tenor[k,j] >= 1) { 
             dt[i,j] = -(1/hazard[k,j])*log((1-u[i,j])/(exp(-cumulativeHazard[k-1,j])))
             defaultTenor[i,j] = tenor[k,j] + dt[i,j] 
             } 
             } 
             if (defaultTenor[i,j] > 5 ) { 
             defaultTenor[i,j] = 0
             } 
             if (defaultTenor[i,j] < 1 & defaultTenor[i,j] > 0) {
             defaultTenor[i,j] = 0.50
             }
}
}

defaultTime = defaultTenor
DL = matrix(0,nrow(defaultTime),ncol(defaultTime)) 
PL = matrix(0,nrow(defaultTime),ncol(defaultTime))

m  = defaultTime 
# I extract the Default Times for each asset removing the zeros and sorting 
m1 = sort(m[,1][!m[,1]==0])
m2 = sort(m[,2][!m[,2]==0])
m3 = sort(m[,3][!m[,3]==0])
m4 = sort(m[,4][!m[,4]==0])
m5 = sort(m[,5][!m[,5]==0])

# Linear Interpolation of the default times  
getDF1 = rep(NA, length(m1))
getDF2 = rep(NA, length(m2)) 
getDF3 = rep(NA, length(m3)) 
getDF4 = rep(NA, length(m4))
getDF5 = rep(NA, length(m5))
r1 = rep(NA, length(m1))
r2 = rep(NA, length(m2)) 
r3 = rep(NA, length(m3)) 
r4 = rep(NA, length(m4))
r5 = rep(NA, length(m5))

for (i in 1:length(m1)) { 
   for (k in 1:5) { 
       if (m1[i] >= e[k,1] & m1[i] <= e[k+1,1]) {  
            r1[i] = e[k,2] + (e[k+1,2]-e[k,2])*(m1[i]-k)
              getDF1[i] = exp(-r1[i]*m1[i])
} 
} 
} 
for (i in 1:length(m2)) { 
   for (k in 1:5) { 
       if (m2[i] >= e[k,1] & m2[i] <= e[k+1,1]) {  
            r2[i] = e[k,2] + (e[k+1,2]-e[k,2])*(m2[i]-k)
              getDF2[i] = exp(-r2[i]*m2[i])
} 
} 
} 
for (i in 1:length(m3)) { 
   for (k in 1:5) { 
       if (m3[i] >= e[k,1] & m3[i] <= e[k+1,1]) {  
            r3[i] = e[k,2] + (e[k+1,2]-e[k,2])*(m3[i]-k)
              getDF3[i] = exp(-r3[i]*m3[i])
} 
} 
} 
for (i in 1:length(m4)) { 
   for (k in 1:5) { 
       if (m4[i] >= e[k,1] & m4[i] <= e[k+1,1]) {  
            r4[i] = e[k,2] + (e[k+1,2]-e[k,2])*(m4[i]-k)
              getDF4[i] = exp(-r4[i]*m4[i])
} 
} 
} 
for (i in 1:length(m5)) { 
   for (k in 1:5) { 
       if (m5[i] >= e[k,1] & m5[i] <= e[k+1,1]) {  
            r5[i] = e[k,2] + (e[k+1,2]-e[k,2])*(m5[i]-k)
              getDF5[i] = exp(-r5[i]*m5[i])
} 
} 
}

# I calculate the Default Leg for each asset 

DL1 = 1/5*(1-RR)*getDF1
DL2 = 1/5*(1-RR)*getDF2
DL3 = 1/5*(1-RR)*getDF3
DL4 = 1/5*(1-RR)*getDF4
DL5 = 1/5*(1-RR)*getDF5

# now I calculate the Protection Leg or Premium Leg

PL1 = getDF1*m1*5/5 
PL2 = PL1 + getDF2*(m2-m1)*4/5 
PL3 = PL2 + getDF3*(m3-m2)*3/5
PL4 = PL3 + getDF4*(m4-m3)*2/5
PL5 = PL4 + getDF5*(m5-m4)*1/5

# I calculate the spread as follows: 

s1 = mean(DL1)/mean(PL1) 
s2 = mean(DL2)/mean(PL2) 
s3 = mean(DL3)/mean(PL3) 
s4 = mean(DL4)/mean(PL4) 
s5 = mean(DL5)/mean(PL5)

as you can see I extracted the default times for every single asset (in this case 5 names). Now I have to calculate the spread of each instrument and in a continuous time. the key concept is to calculate the Default Leg and the Protection Leg. Default Leg is calculated for each reference as following:

DL = 1/5 x (1-RR) x (Discounted Factors interpolated)

the protection leg is calculated in a different for each reference name.

For the 1st default, is

PL(1) = (DF_Interpolated_1 x m1 ) x (5/5).

For the 2nd Default I have:

PL(2) = PL(1) + (DF_Interpolated_2 x (m2 - m1) x (4/5)).

for the 3rd Default I have

PL(3) = PL(2) + (DF_Interpolated_3 x (m3 - m2) x (3/5)).

for the 4th Default I have

PL(4) = PL(3) + (DF_Interpolated_4 x (m4 - m3) x (2/5)).

for the 5th Default I have:

PL(5) = PL(4) + (DF_Interpolated_5 x (m5 - m4) x (1/5)).

What it is not clear is how to calculate the difference between two vectors having different lengths, since in my first round of 10,000 simulations, length(m1) = 1388, length(m2) = 406, length(m3) = 689, length(m4) = 1493, length(m5) = 1255.

if you need more details about the pricing model, please let me know.

$\endgroup$
3
  • $\begingroup$ Please ensure that all data necessary to run this code is provided. Are you sure that tenor[k, j] is not NA the moment the error occurs? $\endgroup$
    – Bob Jansen
    Apr 4 '20 at 17:38
  • $\begingroup$ Hi Bob, yes now the code works out correctly because as you stated tenor was set as a matrix of NA , whereas it has to be a matrix of Zeroes. Yes I can provide more data, but at this stage , I am pretty sure I can continue with the code. thanks a lot $\endgroup$ Apr 6 '20 at 21:41
  • $\begingroup$ I got the point. But my If condition is based on an inequality between a matrix c and another matrix t. Anyway, in the original code the tenor is not a matrix but it is defined as an Integer , hence A Number . I will post the new code with more data and comments so you can have a big picture. Many thanks in advance. $\endgroup$ Apr 8 '20 at 12:07
2
$\begingroup$

You have to initialize the tenor matrix to zero instead of NA. If you compare a number with NA it will yield an NA. Evaluating NA like this in an if-statement gives the error described.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.