# Chapter 7. Simulation Methods for VaR for Options and Bonds

##### Listing 7.1: Transformation in R Last edited: 2011

x=seq(-3,3,by=0.1)
plot(x,pnorm(x),type="l")

##### Listing 7.2: Transformation in Matlab Last edited: 2011

x=-3:0.1:3;
plot(x,normcdf(x))


##### Listing 7.3: Various RNs in R Last edited: 2011

set.seed(12)
S=10
runif(S)
rnorm(S)
rt(S,4)

##### Listing 7.4: Various RNs in Matlab Last edited: December 2016

rng default;
S=10
rand(S,1)
randn(S,1)
trnd(4,S,1)


##### Listing 7.5: Price bond in R Last edited: Augua 2016

yield=c(5.00, 5.69, 6.09, 6.38, 6.61, 6.79, 6.94, 7.07, 7.19, 7.30)
TT=length(yield)
r=0.07
Par=10
coupon=r*Par
cc=1:10*0+coupon
cc[10]=cc[10]+Par
P=sum(cc/((1+yield/100)^(1:TT)))

##### Listing 7.6: Price bond in Matlab Last edited: 2011

yield = [5.00 5.69 6.09 6.38 6.61 6.79 6.94 7.07 7.19 7.30];
T = length(yield);
r=0.07;
Par=10;
coupon=r*Par;
cc=zeros(1,T)+coupon;
cc(10)=cc(10)+Par;
P=sum(cc./((1+yield./100).^(1:T)))


##### Listing 7.7: Simulate yields in R Last edited: August 2016

set.seed(12)
sigma = 1.5
S = 8
r = rnorm(S,0,sigma)
ysim = matrix(nrow=TT,ncol=S)
for (i in 1:S) ysim[,i]=yield+r[i]
matplot(ysim,type='l')

##### Listing 7.8: Simulate yields in Matlab Last edited: 2011

S = 8
sigma = 1.5
randn('state',123);
r = randn(1,S)*sigma;
ysim=nan(T,S);
for i=1:S
ysim(:,i)=yield+r(i);
end
ysim = repmat(yield',1,S)+repmat(r,T,1);
plot(ysim)


##### Listing 7.9: Simulate bond prices in R Last edited: August 2016

SP = vector(length=S)
for (i in 1:S){
SP[i] = sum(cc/((1+ysim[,i]/100)^(TT)))
}
SP = SP-(mean(SP) - P)
barplot(SP)
hist(SP,probability=TRUE)
x=seq(6,16,length=100)
lines(x, dnorm(x, mean = mean(SP), sd = sd(SP)))

##### Listing 7.10: Simulate bond prices in Matlab Last edited: 2011

SP = nan(S,1);
for s = 1:S
SP(s) = sum(cc./(1+ysim(:,s)'./100).^((1:T))) ;
end
SP = SP-(mean(SP) - P);
bar(SP)
histfit(SP)


##### Listing 7.11: Black-Scholes valuation in R Last edited: August 2016

P0 = 50
sigma = 0.2
r = 0.05
TT = 0.5
X = 40
f = bs(X,P0,r,sigma,TT)

##### Listing 7.12: Black-Scholes valuation in Matlab Last edited: 2011

P0 = 50;
sigma = 0.2;
r = 0.05;
T = 0.5;
X = 40;
f = bs(X,P0,r,sigma,T);


##### Listing 7.13: Black-Scholes simulation in R, continued Last edited: August 2016

S = 1e6
set.seed(12)
F = P0*exp(r*TT)
ysim = rnorm(S,-0.5*sigma*sigma*TT,sigma*sqrt(TT))
F=F*exp(ysim)
SP = F-X
SP[SP<0] = 0
fsim = SP*exp(-r*TT)

##### Listing 7.14: Black-Scholes simulation in Matlab, continued Last edited: 2011

randn('state',0);
S = 1e6;
F = P0*exp(r*T);
ysim = randn(S,1)*sigma*sqrt(T)-0.5*T*sigma^2;
F=F*exp(ysim);
SP = F-X;
SP(find(SP < 0)) = 0;
fsim =SP*exp(-r*T) ;


##### Listing 7.15: Option density plots in R, continued Last edited: 2011

hist(F,probability=TRUE,ylim=c(0,0.06))
x=seq(min(F),max(F),length=100)
lines(x, dnorm(x, mean = mean(F), sd = sd(SP)))
hist(fsim,nclass=100,probability=TRUE)

##### Listing 7.16: Option density plots in Matlab, continued Last edited: 2011

histfit(F);
hist(fsim,100);


##### Listing 7.17: Simulate VaR in R Last edited: 2011

set.seed(1)
S = 1e7
s2 = 0.01^2
p = 0.01
r = 0.05
P = 100
ysim = rnorm(S,r/365-0.5*s2,sqrt(s2))
Psim = P*exp(ysim)
q = sort(Psim-P)
VaR1 = -q[p*S]

##### Listing 7.18: Simulate VaR in Matlab Last edited: 2011

randn('state',0);
S = 1e7;
s2 = 0.01^2;
p = 0.01;
r = 0.05;
P = 100;
ysim = randn(S,1)*sqrt(s2)+r/365-0.5*s2;
Psim = P*exp(ysim);
q = sort(Psim-P);
VaR1 = -q(S*p)


##### Listing 7.19: Simulate option VaR in R Last edited: August 2016

TT = 0.25;
X = 100;
sigma = sqrt(s2*250);
f = bs(X,P,r,sigma,TT)
fsim = bs(X,Psim,r,sigma,TT-(1/365))
q = sort(fsim$Call-f$Call)
VaR2 = -q[p*S]

##### Listing 7.20: Simulate option VaR in Matlab Last edited: 2011

T = 0.25;
X = 100;
sigma = sqrt(s2*250);
f = bs(X,P,r,sigma,T);
fsim=bs(X,Psim,r,sigma,T-(1/365));
q = sort(fsim.Call-f.Call);
VaR2 = -q(p*S)


##### Listing 7.21: 2 Options Last edited: August 2016

X1 = 100
X2 = 110
f1 = bs(X1,P,r,sigma,TT)
f2 = bs(X2,P,r,sigma,TT)
f2sim = bs(X2,Psim,r,sigma,TT-(1/365))
f1sim = bs(X1,Psim,r,sigma,TT-(1/365))
q = sort(f1sim$Call+f2sim$Put+Psim-f1$Call-f2$Put-P);
VaR3 = -q[p*S]

##### Listing 7.22: 2 Options Last edited: 2011

X1 = 100;
X2 = 110;
f1 = bs(X1,P,r,sigma,T);
f2 = bs(X2,P,r,sigma,T);
f1sim=bs(X1,Psim,r,sigma,T-(1/365));
f2sim=bs(X2,Psim,r,sigma,T-(1/365));
q = sort(f1sim.Call+f2sim.Put+Psim-f1.Call-f2.Put-P);
VaR3 = -q(p*S)


##### Listing 7.23: Simulated two asset returns in R Last edited: 2011

library (MASS)
mu = c(r/365,r/365)
Sigma = matrix(c(0.01, 0.0005, 0.0005, 0.02),ncol=2)
set.seed(12)
y = mvrnorm(S,mu,Sigma)

##### Listing 7.24: Simulated two asset returns in Matlab Last edited: 2011

mu = [r/365 r/365]';
Sigma = [0.01 0.0005; 0.0005 0.02];
randn('state',12)
y = mvnrnd(mu,Sigma,S);


##### Listing 7.25: Two asset VaR in R Last edited: 2011

K=2
P = c(100,50)
x = c(1,1)
Port = P %*% x
Psim = matrix(t(matrix(P,K,S)),ncol=K)*exp(y)
PortSim = Psim %*% x
q = sort(PortSim-Port[1,1])
VaR4 = -q[S*p]

##### Listing 7.26: Two asset VaR Matlab Last edited: 2011

K = 2;
P = [100 50]';
x = [1 1]';
Port = P'*x;
Psim = repmat(P,1,S)' .*exp(y);
PortSim=Psim * x;
q = sort(PortSim-Port);
VaR4 = -q(S*p)


##### Listing 7.27: A 2 asset case in R with an option Last edited: August 2016

f = bs(P[2],P[2],r,sigma,TT)
fsim = bs(P[2],Psim[,2],r,sigma,TT-(1/365))
q = sort(fsim$Call+Psim[,1]-f$Call-P[1]);
VaR5 = -q[p*S]

##### Listing 7.28: A 2 asset case in Matlab with an option Last edited: 2011

f = bs(P(2),P(2),r,sigma,T);
fsim=bs(P(2),Psim(:,2),r,sigma,T-(1/365));
q = sort(fsim.Call+Psim(:,1)-f.Call-P(1));
VaR5 = -q(p*S)