This program is in Matlab to simulate Brownian motion
Generating GBM:
T = 1; N = 500; dt = T/N;
t = 0:dt:T;
dW = sqrt(dt)*randn(1,N);
mu = 0.1; sigma = 0.01;
plot(t,exp(mu*t + sigma*[0,cumsum(dW)]))
the command "randn(1,N)
" what is it for is this generate array from 1 to N of of independent Gaussian random numbers, all of which are $N(0;1)$.
because I know that the command The command will be randn(M,N)
, which generates an $M \times N
$ matrix of independent Gaussian random numbers, all of which are $N(0;1)$.
and please if anyone can explain to me what the line
"dW = sqrt(dt)*randn(1,N);
" means?
Thanks
Best Answer
First, this code returns the analytical solution to the Stratonovich stochastic differential equation for geometric Brownian motion:
$$dY = \mu Y dt+\sigma Y \circ dW$$
The full solution to this is:
$$Y(t) = Y_0 \exp(\mu t+\sigma W_t)$$
Your code assumes $Y_0=1$. The solution is slightly different if you use the Itô formulation.
The function
randn
returns random variates (values) independently drawn from the normal distribution. So,randn(1,N)
returns a sequence of independent normally-distributed values (500 in this case).N
in this case is just the number of time-steps. These are used to produce the independent Wiener increments, $dW_t$. Thus,dW = sqrt(dt)*randn(1,N)
producesN
independent Wiener increments where the standard deviation is equal to the square root of the time span,dt
, between them $\dagger$. These correspond to theN
points in the vectort
. The code[0,cumsum(dW)]
integrates these increments to produce $W_t$, a Wiener process, which is also called standard Brownian motionFor further details on SDEs, Brownian motion, and simulating them with Matlab I recommend this excellent paper:
The URL to the Matlab files in the paper won't work, use this one: http://personal.strath.ac.uk/d.j.higham/algfiles.html
$\dagger$ Why is
dW
scaled bysqrt(dt)
instead of justdt
like when numerically integrating ODEs? A discrete Wiener process by definition has a variance equal to the time step between the increments, in this casedt
. However, to scale output ofrandn
properly, we need to express this in terms of a standard deviation.