diff options
Diffstat (limited to 'build_SC.R')
| -rw-r--r-- | build_SC.R | 88 |
1 files changed, 50 insertions, 38 deletions
@@ -315,65 +315,76 @@ stonln1.data <- getdealdata("stonln1") A <- SPmatrix2(stonln1.portfolio$SC, getdealdata("stonln1"))
S <- 1 - sapply(stonln1.portfolio$SC, attr, "recov")
issuerweights <- stonln1.portfolio$notional/sum(stonln1.portfolio$notional)
-testC <- lossrecovdist(A$DP[,39], A$PP[,39], issuerweights, S, 100)
-test <- lossrecovdist(A$DP[,39], A$PP[,39], issuerweights, S, 100, useC=FALSE)
-poire <- recovdistC(A$DP[,39], A$PP[,39], issuerweights, S, 100)
-test2 <- lossdistribprepayC.joint(A$DP[,39], A$PP[,39], w, S, 100)
dp <- A$DP
pp <- A$PP
-dpmod <- MFupdate.prob(Z, w.mod, rho, dp)
-ppmod <- MFupdate.prob(-Z, w.mod, rho, pp)
-cl <- makeCluster(6)
-clusterExport(cl, list("shockprob", "rho", "Z", "issuerweights", "lossdistribprepayC.joint",
- "lossrecovdist.joint.term", "lossdistribC.joint", "Ngrid"))
-MFtest.joint <- MFlossrecovdist2(cl, w.mod, Z, rho, dp, dpmod, pp, ppmod, issuerweights, 1-S, Ngrid=201)
-MFtest <- MFlossrecovdist(w.mod, Z, rho, dp, dpmod, pp, ppmod, issuerweights, 1-S, Ngrid=201)
+Smat <- matrix(S, length(issuerweights), ncol(dp))
+#no correlation setup
+test <- lossrecovdist.joint.term(dp, pp, issuerweights, Smat, Ngrid)
-## build scenarios:
-## compute the joint density of (R/D, D)
-## u=y/(x+y)
-## v=x+y
-## colSums(distchv[t,,]) is the density of D
-## rowSums(distchv[t,,]) is the density of R/D, which is
-distchv <- array(0, dim=dim(MFtest.joint))
+dpmod <- MFupdate.prob(Z, w.mod, rho, dp)
+ppmod <- MFupdate.prob(-Z, w.mod, rho, pp)
+## cl <- makeCluster(6)
+## clusterExport(cl, list("shockprob", "rho", "Z", "issuerweights", "lossdistribprepayC.joint",
+## "lossrecovdist.joint.term", "lossdistribC.joint", "Ngrid"))
+dist <- MFlossrecovdist.prepay(w.mod, Z, rho, dp, dpmod, pp, ppmod, issuerweights, 1-S, Ngrid=201, TRUE)
+dist.joint <- MFlossdist.prepay.joint(cl, w.mod, Z, rho, dp, dpmod,
+ pp, ppmod, issuerweights, 1-S, Ngrid=201, FALSE)
+## two ways to compute the joint (D, R) distribution
+## first using the actual function
+dist.joint2 <- MFlossdist.prepay.joint(cl, w.mod, Z, rho, dp, dpmod,
+ pp, ppmod, issuerweights, 1-S, Ngrid=201, TRUE)
+## second, by doing a change of variable seems to work better for now
+distDR <- array(0, dim=dim(dist.joint))
+u <- seq(0, 1, length=Ngrid)
+v <- seq(0, 1, length=Ngrid)
for(t in 1:ncol(dp)){
- dist <- MFtest.joint[t,,]
- u <- seq(0, 1, length=Ngrid)
- v <- seq(0, 1, length=Ngrid)
- y <- u %o% v
- x <- sweep(-y, 2, v, "+")
- xgrid <- round(x*(Ngrid-1))
- ygrid <- round(y*(Ngrid-1))
for(i in 1:Ngrid){
for(j in 1:Ngrid){
- distchv[t,i,j] <- dist[xgrid[i,j]+1, ygrid[i,j]+1]*seq(0,1, length=Ngrid)[j]
+ index <- i+j
+ if(index<=Ngrid){
+ distDR[t,index,j] <- distDR[t,index,j] + dist.joint[t,i,j]
+ }
}
}
}
+## compute E(R|D)
+R <- matrix(0, Ngrid, ncol(dp))
+for(t in 1:ncol(dp)){
+ R[,t] <- (sweep(distDR[t,,], 1, rowSums(distDR[t,,]), "/") %*% support)/support
+}
+
n.scenarios <- 50
percentiles <- (seq(0, 1, length=n.scenarios+1)[-1]+
seq(0, 1, length=n.scenarios+1)[-(n.scenarios+1)])/2
+## compute scenariosd
scenariosd <- matrix(0, n.scenarios, ncol(dp))
-scenariosr <- matrix(0, n.scenarios, ncol(dp))
-for(t in 1:39){
- D <- colSums(distchv[t,,])/sum(colSums(distchv[t,,]))
- Dfun <- splinefun(c(0, cumsum(D)),c(0, seq(0, 1, length=Ngrid)), "monoH.FC")
+support <- seq(0, 1, length=Ngrid)
+for(t in 1:ncol(dp)){
+ D <- rowSums(distDR[t,,])
+ D <- D/sum(D)
+ Dfun <- splinefun(c(0, cumsum(D)), c(0, support), "monoH.FC")
dval <- round(Dfun(percentiles)*(Ngrid-1))
- dvallow <- floor(Dfun(percentiles)*(Ngrid-1))
- dvalup <- floor(Dfun(percentiles)*(Ngrid-1))
+ ## dvallow <- floor(Dfun(percentiles)*(Ngrid-1))
+ ## dvalup <- ceil(Dfun(percentiles)*(Ngrid-1))
scenariosd[,t] <- Dfun(percentiles)
- for(i in 1:n.scenarios){
- scenariosd[i,t]*
- scenariosr[i,t] <- seq(0,1,length=Ngrid)[which(cumsum(distchv[t,,dval[i]])/
- sum(distchv[t,,dval[i]])>=0.5)[1]]
- }
}
+## compute scenariosr
+scenariosr <- matrix(0, n.scenarios, ncol(dp))
+support <- seq(0, 1, length=Ngrid)
+for(t in 1:ncol(dp)){
+ scenariosr[,t] <- R[scenariosd[,t]*Ngrid+1,t]
+}
+
+
+
+
dates <- seq(stonln1.data$"Deal Next Pay Date", stonln1.data$maturity, by="3 months")
+dates <- dates[today()]
cdrfromscenarios <- function(scenario, dates){
cdr <- matrix(0, nrow(scenarios), ncol(scenarios))
@@ -383,10 +394,11 @@ cdrfromscenarios <- function(scenario, dates){ return( cdr )
}
-cdrmonthly <- matrix(0, n.scenarios, 38*3)
+cdrmonthly <- matrix(0, n.scenarios, *3)
for(i in 1:n.scenarios){
cdrmonthly[i,] <- cdr[i, as.numeric(rbind(1:38, 1:38, 1:38))]
}
+
scenariosr.temp <- scenariosr[,-1]
recoverymonthly <- matrix(0, n.scenarios, 38*3)
for(i in 1:n.scenarios){
|