josef teichmann :: personal homepage at ETH Zurich

Lecture Wahrscheinlichkeitstheorie und Statistik

We provide several commented implementations of important concepts to accompany the lecture. We shall use R for this purpose. Just paste the following code into your most favorite R editor and execute.

Calculation of the price of a European Call contract in a SABR model

S0 = 100; #Spot price
K = 98; #strike 
alpha0 = 0.852; #alpha parameter
beta= 0.5; #beta parameter
rho = 0.5; #rho parameter
nu = 0.42; #variance of volatility paramater
T = 1; #maturity price
m = 250; #no. of time step points
dt = T/m; #length of time discretization

M = 100; #no of repetitions of the MC experiment to see the normal distribution around the result
n = 10000; #no of repetitions of the MC simulations

payoff_mc_value<-rep(0,times=M); #vector with zeros of length M which stores the results of the M MC experiments
emp_var<-rep(0,times=M); #vector with zeros of length M which stores the empirical variance of the M MC experiments
cl<-rainbow(M)

trajectoryS<-rep(S0,times=m);
trajectoryA<-rep(alpha0,times=m);
par(mfrow=c(1,3)) 
plot(trajectoryS,type='l')
    for (i in 1:M)
    {
     S<-rep(S0,times=n);
     alpha<-rep(alpha0,times=n);
     trajectoryS<-rep(S0,times=m);
     trajectoryA<-rep(alpha0,times=m);
     for (j in 1:m)
     {
        B1=rnorm(n,0,1)*sqrt(dt);
        B2=rnorm(n,0,1)*sqrt(dt);
	S = S + alpha*S^beta*B1;
	alpha = alpha*exp(nu*(rho*B1+sqrt(1-rho^2)*B2)-nu^2*dt/2);
        trajectoryS[j]=S[1]
        trajectoryA[j]=alpha[1]
     }
     payoff_mc_value[i]=1/n*sum((S-K)*(S-K >=0));
     emp_var[i]=1/n*sum(((S-K)*(S-K>=0)-payoff_mc_value[i])^2)
     lines(trajectoryS,col=cl[i])
    }
plot(sqrt(emp_var)/sqrt(n),ylim=c(0,2*max(sqrt(emp_var)/sqrt(n))))
hist(payoff_mc_value)

Illustration zu Satz 7.1

N=10000 # number of experiments 
M=1000 # number of repetitions of N experiments
mu<-0.1 # expectation of one experiment
sigma<-0.33 # std of one experiment
m<-rep(0,times=M)
v<-rep(0,times=M)
t<-rep(0,times=M)
for (i in 1:M) # repeat N experiments M times to obtain a histogram of outcomes of estimator values
{
x<-rnorm(N,mu,sigma)
m[i]=(mean(x)-mu) # normalized difference between true expectation mu and estimated expectation normalized by the true standard deviation of the estimator, this should be N(0,1)
v[i]= var(x)/sigma^2 # normalized estimator variance by true variance, this value times (N-1) should be chi-square with df=N-1
t[i]=(mean(x)-mu)/(sqrt(v[i])*sigma/sqrt(N)) # difference between true and expectation mu and estimated expectation normalized by the estimated standard deviation
}
par(mfrow=c(1,3)) # draw three histograms with true probabilities together with the corresponding theoretical densities
hist(m/(sigma/sqrt(N)),prob=TRUE)
curve(dnorm(x,0,1),add=TRUE)
hist((N-1)*v,prob=TRUE)
curve(dchisq(x,df=N-1),add=TRUE)
hist(t,prob=TRUE)
curve(dt(x,df=N-1),add=TRUE)

Lecture Wahrscheinlichkeitstheorie und Statistik

Calculation of the price of a European Call contract in a SABR model

Illustration zu Satz 7.1

Bayesianisches Updating in einem Bernoulliexperiment