# PD and IR

Dec 10th, 2013
410
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. # visit http://reakkt.com for more
2.
3. # pd <- 0.01  # probability of default
4.
5. # cr <- 0     # (annual) coupon rate
6.
7. fa <- 100   # face value
8.
9. pd <- 0.05  # probability of default
10.
11. rv <- 0.37   # recovery rate
12.
13. rf <- 0.00429  # risk free rate
14.
15. rl <- 1     # recovery time
16.
17. n  <- 3     # years
18.
19. # expected (discounted) return
20. pv <- function(fa,n,cr,rf) {
21.
22.   -fa+sum(sapply(1:n, function(i) (fa*cr)/(1+rf)^i))+fa/(1+rf)^n
23.
24. }
25.
26. ## expected (discounted) return at default
27. # di - default time/year
28. # rv - recovery rate
29. # rl - recovery delay/time
30. de <- function(fa,di,cr,rv,rf,rl=1) {
31.
32.   rf.amount <- (fa*rv)/(1+rf)^(di+rl)
33.
34.   if (di==0)
35.
36.     # default in year 0
37.     return(-fa+rf.amount) else
38.
39.     # default after first year
40.     return(-fa+sum(sapply(1:di, function(i) (fa*cr)/(1+rf)^i))+rf.amount)
41.
42. }
43.
44. pv(fa=100,n=3,cr=0.05,rf)
45.
46. de(fa=100,di=0,cr=0.05,rv=0.4,rf,rl=1) # default in year one; no coupon payments
47. de(fa=100,di=1,cr=0.05,rv=0.4,rf,rl=1) # default after first coupon payment
48. de(fa=100,di=2,cr=0.05,rv=0.4,rf,rl=1) # default after second coupon payment
49. de(fa=100,di=3,cr=0.05,rv=0.4,rf,rl=1) # default after third coupon payment
50.
51. # k  - initial number of entities
52. # pd - annual/periodic PD
53. # n  - years/periods
54. PortfolioSurvival <- function(k,pd,n) {
55.
56.   defaulted  <- numeric(n+1)
57.   portfolio  <- rep(k,n+1)
58.
59.   for (i in 0:n) {
60.
61.     if (i==0) defaulted[1] <- ceiling(k*pd) else defaulted[i+1] <- ceiling(portfolio[i]*pd)
62.
63.     portfolio[i+1] <- ifelse(i==0,k,portfolio[i])-defaulted[i+1]
64.
65.   }
66.
67.   list(defaulted=defaulted,
68.        portfolio=portfolio)
69.
70. }
71.
72. # k  - no. of assets
73. # pd - probability of default
74. # n  - no. of years/periods
75. RateGivenPD <- function(k,pd,n,
76.                         fa,cr,rf,rv,rl) {
77.
78.   ps <- PortfolioSurvival(k,pd,n)
79.
80.   portfolio <- ps\$portfolio
81.   defaulted <- ps\$defaulted
82.
83.   cr   <- 0
84.   retS <- -1
85.
86.   while(retS<0) {
87.
88.     cr <- cr+0.001
89.
90.     ret <- numeric(n+1)
91.
92.     for (i in 0:n) {
93.
94.       if (i!=n) {
95.
96.         ret[i+1] <- de(fa,i,cr,rv,rf,rl)*defaulted[i+1]
97.
98.       } else {
99.
100.         ret[i+1] <- de(fa,i,cr,rv,rf,rl)*defaulted[i+1]+pv(fa,n,cr,rf)*portfolio[n+1]
101.
102.       }
103.
104.     }
105.
106.     ret
107.     k*fa
108.     retS <- sum(ret)
109.
110.   }
111.
112.   return( list(cr=cr,
113.                ps=ps,
114.                ret=ret,
115.                retS=retS) )
116.
117. }
118.
119. ##
120.
121. k  <- 1e6   # no. of assets
122.
123. pd.seq <- seq(0,0.99,0.01)
124.
125. cr.arr <- numeric(length(pd.seq))
126.
127. for (pd in pd.seq) cr.arr[which(pd.seq==pd)] <- RateGivenPD(k,pd,n,fa,cr,rf,rv,rl)\$cr
128.
129. plot(pd.seq,cr.arr,type="l",
130.      log="y",
131.      xlab="PD",ylab="IR (log)")
132. # abline(h=rf,col="Green",lty="dashed")
133.
134. cbind(pd.seq,cr.arr)