Monte Carlo simulation (Matlab)

This program is under development.
is complete, remove the {{develop}} tag.

Default MATLAB functions

In the MATLAB Finance Toolbox, a function called portsim.m allows to perform Monte Carlo simulations. It uses the mvnrnd.m function (from the Statistics Toolbox) which returns a random vector given a chosen variance and expectation.

Those two functions can be used to perform Monte Carlo simulations of simple ITo diffusion like:

where r, σ and λ are constants of time and S.

The simplest Monte Carlo Simulation with Matlab

Here is the simplest Monte Carlo possible using the basic functions of MATLAB (to understand the mechanisms of this simulation).

The parameters of a Monte Carlo simulation in finance

The goal is the simulation of trajectories of several stochastic process defined in (D.1). S is the price of the stocks, S0 is their prices at the begining of the simulation, , σ their volatilities (annual), T is the length of the simulations (step is the time step), and nb_traj is the number of trajectories to simulate.

<<test_simplest_montecarlo0.m>>=
S0      = 1000
mu      = .05
sigma   = .12
T       =   5
nb_traj = 100
step    = 1/255
[S, t]  = simplest_montecarlo0( sigma, T, nb_traj, S0, mu, step);
figure('Color',[0.9412 0.9412 0.9412 ]);
plot(t,S, 'linewidth',2);
axis([t(1) t(end) min(S(:)) max(S(:))]);
xlabel('time');ylabel('value');
title(sprintf('%d trajectories', size(S,2)));


Simplest Monte Carlo from scratch

Obtained simulation
We will use the fact that a solution of (D.1) is:

Of course we will have to take into account the size of the time steps.

<<simplest_montecarlo0.m>>=
function [S, t] = simplest_montecarlo0( sigma, T, nb_traj, S0, mu, step)
% SIMPLEST_MONTECARLO - the simplest Monte Carlo simulation with Matlab
% use:
%  S = simplest_montecarlo( sigma, T, nb_traj, S0, mu, step)
% example:
%  [S, t] = simplest_montecarlo( .12, 5, 50, 100, .05, 1/255);
simplest_monovariate_montecarlo

<<simplest_monovariate_montecarlo>>=
nT = ceil(T/step);
W  = sigma * sqrt(step) * cumsum(randn(nT, nb_traj));
c  = repmat((mu - sigma^2/2) *step * (1:nT)',1,nb_traj);
S  = [repmat(S0,1,nb_traj); S0 * exp( c + W)];
if nargout > 1
t = [0;step * (1:nT)'];
end


Mutlivariate Monte Carlo from scratch

 Obtained mutli assets simulation Now we will simulate this kind of diffusion (for N correlated assets): Here we won't solve the diffusion process, but rather simulate a discretized version of it. That is we will write each on the N equation of (D.1') as: In this expression, δt is the time step, is a correlated i.i.d. gaussian noise.
<<simplest_montecarlo.m>>=
function [S, t] = simplest_montecarlo( mode, sigma, T, nb_traj, S0, mu, step)
% SIMPLEST_MONTECARLO - the simplest Monte Carlo simulation with Matlab
% use:
%  S = simplest_montecarlo( mode, sigma, T, nb_traj, S0, mu, step)
% example:
%  [S, t] = simplest_montecarlo( 'simplest', .12, 5, 50, 100, .05, 1/255);
switch lower(mode)
case 'test-simple'
test_simplest_montecarlo
case 'test-multi'
test_simplest_multiple_montecarlo
case {'simple', 'simplest'}
simplest_monovariate_montecarlo
case {'multi', 'multiple'}
simplest_multivariate_montecarlo
end


Simple Multivariate Monte Carlo

• At first we use the corr2cov function, which returns the standard deviation and the correlation matrix associated with a covariance matrix. The only goal (before reverting it back to the covariance matrix), is to take into account the used step.
• Then the Cholesky decomposition (chol.m Matlab function) is used to be able to build a correlated multivariate gaussian noise.
• The randn function can generate multivariate independant gaussian variables, here we need nT sample of realisations.
• Then we use a typical Matlab trick in two lines:
Sd  = repmat({Sd},1,nb_traj);
Sdb = blkdiag(Sd{:});

the first one build a cellarray (a list) of nb_traj replications of the Sd matrix ; the second one give them as arguments to the blkdiag function. those two lines can be read as:
• Thanks to that trick the product of our generated gaussians by Sdb is:
id est:
so we went from nbtraj n-dimensionnal i.i.d. gaussian noises dX to nbtraj n-dimensionnal correlated i.i.d. normal noises dW.
<<simplest_multivariate_montecarlo>>=
nT  = ceil(T/step);
[s,c] = cov2corr(sigma);
s   = s * sqrt(step);
sigma = corr2cov(s,c);
Sd  = chol(sigma);
dW  = randn(nT,length(S0)*nb_traj);
Sd  = repmat({Sd},1,nb_traj);
Sdb = blkdiag(Sd{:});
dWc = dW * Sdb;
c   = repmat(1+mu*step, nT, nb_traj);
S   = [repmat(S0, 1, nb_traj); repmat(S0, nT, nb_traj) .* cumprod(c + dWc)];
S   = reshape(S, [nT+1, size(S0,2), nb_traj]);
if nargout > 1
t = [0;step * (1:nT)'];
end


Tests simplest simulations

Each of those two blocks only call the main function. They:

• initialize the parameters
• call the main function in the correct mode
• plot the obtained result
<<test_simplest_motecarlo>>=
mode    = 'simplest';
S0      = 100;
mu      = .05;
sigma   = .12;
T       = 5;
nb_traj = 100;
step    = 1/255;
[S, t]  = simplest_montecarlo( mode, sigma, T, nb_traj, S0, mu, step);
figure('Color',[0.9412 0.9412 0.9412 ]);
plot(t,S, 'linewidth',2);
axis([t(1) t(end) min(S(:)) max(S(:))]);
xlabel('time');ylabel('value');
title(sprintf('%d trajectories', size(S,2)));

<<test_simplest_multiple_montecarlo>>=
mode    = 'multiple';
T       = 1/255*10;
nb_traj = 2;
step    = 1/255;
S0      = [100, 150];
s       = [.12 .18];
sigma   = corr2cov(s, [1 .3; .3 1]);
mu      = [0.05 0.04];
T       = 5;
nb_traj = 100;
step    = 1/255;
[S, t]  = simplest_montecarlo( mode, sigma, T, nb_traj, S0, mu, step);
figure('Color',[0.9412 0.9412 0.9412 ]);
subplot(2,1,1);
plot(t,reshape(S(:,1,:),size(S,1),size(S,3)), 'linewidth',2);
axis([t(1) t(end) min(S(:)) max(S(:))]);
xlabel('time');ylabel('value');
title(sprintf('%d trajectories - stock 1 (\\sigma=%5.2f%%)', size(S,3), s(1)*100));
subplot(2,1,2);
plot(t,reshape(S(:,2,:),size(S,1),size(S,3)), 'linewidth',2);
axis([t(1) t(end) min(S(:)) max(S(:))]);
xlabel('time');ylabel('value');
title(sprintf('%d trajectories - stock 2 (\\sigma=%5.2f%%)', size(S,3), s(2)*100));


A more generic Monte Carlo simulation in MATLAB

We often need to simulate diffusion process more generic than (D.1) like:

To be able to perform a Monte Carlo simulation of diffusion (D.2), we will need:

1. to simulate the 2-dimensional brownian diffusion process (W1,W2);
2. to simulate the diffusion of σ;
3. to finally obtain a simulation of S.