1 files changed, 23 insertions, 9 deletions
diff --git a/R/tranche_functions.R b/R/tranche_functions.R
index 2b7a6219..f22bdc68 100644
--- a/R/tranche_functions.R
+++ b/R/tranche_functions.R
@@ -450,17 +450,31 @@ dist.transform <- function(dist.joint){
 }
 
 shockprob <- function(p, rho, Z, log.p=F){
-    ## computes the shocked default probability as a function of the copula factor
+  ## computes the shocked default probability as a function of the copula factor
+  ## function is vectorized provided the below caveats:
+  ## p and rho are vectors of same length n, Z is a scalar, returns vector of length n
+  ## p and rho are scalars, Z is a vector of length n, returns vector of length n
+  if(length(p)==1){
     if(rho==1){
-        return(Z<=qnorm(p))
+      return(Z<=qnorm(p))
     }else{
-        pnorm((qnorm(p)-sqrt(rho)*Z)/sqrt(1-rho), log.p=log.p)
+      return(pnorm((qnorm(p)-sqrt(rho)*Z)/sqrt(1-rho), log.p=log.p))
     }
+  }else{
+    result <- double(length(p))
+    result[rho==1] <- Z<=qnorm(p[rho==1])
+    result[rho<1] <- pnorm((qnorm(p[rho<1])-sqrt(rho[rho<1])*Z)/sqrt(1-rho[rho<1]), log.p=log.p)
+    return( result )
+  }
 }
 
 shockseverity <- function(S, Stilde=1, Z, rho, p){
-    #computes the severity as a function of the copula factor Z
-    Stilde * exp( shockprob(S/Stilde*p, rho, Z, TRUE) - shockprob(p, rho, Z, TRUE))
+  ## computes the severity as a function of the copula factor Z
+  result <- double(length(S))
+  result[p==0] <- 0
+  result[p!=0] <- Stilde * exp( shockprob(S[p!=0]/Stilde*p[p!=0], rho[p!=0], Z, TRUE) -
+           shockprob(p[p!=0], rho[p!=0], Z, TRUE))
+  return(result)
 }
 
 dshockprob <- function(p,rho,Z){
@@ -657,21 +671,21 @@ tranche.pvvec <- function(K, L, R, cs){
 }
 
 BClossdist <- function(defaultprob, issuerweights, recov, rho, Z, w,
-                       N=length(recov)+1, defaultflag=FALSE){
+                       N=length(recov)+1, defaultflag=FALSE, n.int=500){
   if(missing(Z)){
     quadrature <- gauss.quad.prob(n.int, "normal")
     Z <- quadrature$nodes
     w <- quadrature$weights
   }
-  LZ <- matrix(0, N, n.int)
-  RZ <- matrix(0, N, n.int)
+  LZ <- matrix(0, N, length(Z))
+  RZ <- matrix(0, N, length(Z))
   L <- matrix(0, N, ncol(defaultprob))
   R <- matrix(0, N, ncol(defaultprob))
   for(t in 1:ncol(defaultprob)){
     for(i in 1:length(Z)){
       g.shocked <- shockprob(defaultprob[,t], rho, Z[i])
       S.shocked <- shockseverity(1-recov, 1, Z[i], rho, defaultprob[,t])
-      temp <- lossrecovdist(g.shocked, 0, issuerweights, S.shocked, N)
+      temp <- lossrecovdist(g.shocked, , issuerweights, S.shocked, N)
       LZ[,i] <- temp$L
       RZ[,i] <- temp$R
     }