# Simulation of a Geometric Brownian Motion in R

The geometric Brownian motion (GBM) is the most basic processes in financial modelling.  Consider a stockprice S(t) with dynamics

.
It can be shown (just use Ito) that the solution to this stochastic differential equation is

where W(t) is a  Brownian Motion.

In other words, a geometric Brownian motion is nothing else than a transformation of a Brownian motion. For this, we sample the Brownian W(t) (this is "f" in the code, and the red line in the graph). This is being illustrated in the following example, where we simulate a trajectory of a Brownian Motion and then plug the values of W(t) into our stock price S(t).
The exponential lines are the expected value of the GBM's due to

An interesting observation is that around timepoint 12 the green GBM is above its expected value, whereas the  blue line is below its expected value, even when they follow the same underlying trajectory of a BM. This may seem odd at first, and in fact I thought I had made a mistake in my program. Whenever I simulated a trajectory, in most cases the blue line was below its expected value. Should the stock price not be above its expectation value when the BM is positive, and below it is negative?
The answer is no. Even though the increment of a BM from 0 to 15 is normally distributed (with mean zero and variance 15), wrapped in an exponential function,  positive deviations of the BM from 0 have an exponential effect on the mean value, whereas negative trajectories are bound from below by 0 and thus have limited downside effect.
The conclusion of this is that the longer the timeline graph goes, the more likely it is for the GBM to be significantly below its expected value.

Finally the R code for the simulation of the Geometric Brownian Motions on the graph:

# 01/05/2013

maturity <- 15
simulation.length <- 10001
dt <-  maturity/(simulation.length-1)

timeline <- seq(0,maturity, dt)

S0<-1
r<-0.05
mu<-0.1
mu0<-0.2
sigma<-0.2
sigma0<-0.375

f <- g <- g0 <- h <- h0 <- rep(0, times=simulation.length)
g0[1] <- h0[1] <- g[1] <-  h[1] <- S0

for(i in 2:simulation.length){
f[i] <- f[i-1]+sqrt(dt)*rnorm(1)
g[i] <- g[1]*exp((mu-(sigma^2)/2)*(i-1)*dt+sigma*f[i])
g0[i] <- g0[1]*exp(mu*(i-1)*dt)
h[i] <- h[1]*exp((mu0-sigma0^2/2)*(i-1)*dt+sigma0*f[i])
h0[i] <- h0[1]*exp(mu0*(i-1)*dt)
}

o_range <- range(f,g,g0,h,h0)

plot(timeline,f, ylim=o_range, type="l", col="coral1")
lines(timeline,g0, col="chartreuse3")
lines(timeline,g, col="chartreuse2")
lines(timeline,h, col="deepskyblue1")
lines(timeline,h0, col="deepskyblue3")

title(main="Geometric Brownian Motion trajectories", col.main="red", font.main=4)

legend(1, o_range[2], c("mu = 0.2,  sigma = 0.375","mu = 0.1,  sigma = 0.2","Brownian motion"), cex=0.8,
col=c("deepskyblue1","chartreuse2","coral1"), pch=1, lty=1);