Arithmetic Brownian Motion Process and SDEs

We discuss various things related to the Arithmetic Brownian Motion Process- these include solution of the SDEs, derivation of its Characteristic Function and Moment Generating Function, derivation of the mean, variance, and covariance, and explanation of the calibration and simulation of the process.

Calibration and Simulation

We will now try to get some feel for the Arithmetic Brownian motion process by using it to model some real financial data - say FTSE100 index. We will need to calibrate the parameters of our process based on this index data. Note calibrating the model is the same thing as fitting the model.

Let’s say we decided to calibrate using FTSE 100 daily observed data over 17 business days. This is just to show how the model can be calibrated so the choice of the data, and the length of the observation is just arbitrary. In reality, the calibration approach will depend on what the model is used for. For example, calibration using historical observations of the process values may be appropriate when the model is used for simulating some function/quantity which depends on the observed values of the process. On the other hand, if the purpose were to use the model to price options for example, then one would have to calibrate the process parameters using the current option prices - this amounts to modelling an implied process whose values may be systematically different from the realised observations. This is a topic for another day, for now let’s get back to out tiny data series, which is shown below:

13-Sep-18, 7,281.6
14-Sep-18, 7,304.0
17-Sep-18, 7,302.1
18-Sep-18, 7,300.2
19-Sep-18, 7,331.1
20-Sep-18, 7,367.3
21-Sep-18, 7,490.2
24-Sep-18, 7,458.4
25-Sep-18, 7,507.6
26-Sep-18, 7,511.5
27-Sep-18, 7,545.4
28-Sep-18, 7,510.2
01-Oct-18, 7,495.7
02-Oct-18, 7,474.6
03-Oct-18, 7,510.3
04-Oct-18, 7,418.3
05-Oct-18, 7,318.5

The Arithmetic Brownian models the change in the value of the process, so we also generate the difference series- simply, today value minus the previous day value- which produces the following series in terms of Python array:

FTSE100_d1=[22.4,-1.9,-1.9,30.9,36.2,122.9,-31.8,49.2,3.9,33.9,-35.2,-14.5,-21.1,35.7,-92.0,-99.8]

This difference series is now going to be the main item of interest. Recall that the mean and variance of a variable can be calculated using the following formulae,

$${\rm Mean}=m_d= \frac{\sum_{i=1}^n{r_i}}{n}$$

$${\rm Variance}=v_d=\frac{\sum_{i=1}^n{(r_i}-m)^2}{n-1}$$

Where the subscript d stands for daily as these are the mean and variance of daily changes/returns. Applying these formulae to our difference series, we get the mean and variance of the series: \(m_d={\rm 2.3}\), \(v_d={\rm 2,932.2}\). The standard deviation, which is the square root of the variance, is 54.1. Here is sample Python code to generate these numbers:

import numpy as np

FTSE100_d1=[22.4,-1.9,-1.9,30.9,36.2,122.9,-31.8,49.2,3.9,33.9,-35.2,-14.5,-21.1,35.7,-92.0,-99.8]

mean=np.average(FTSE100_d1)

variance=np.var(FTSE100_d1,ddof=1)

print(mean)

print(variance)

So each of the daily change is normally distributed, with the above computed mean and variance,

$$R_i =X_i-X_{i-1}=\sim N \left[m={\rm 2.3},v={\rm 2,932.2}\right]$$

The question now is, what is the distribution of the yearly changes, as most often we are interested in the annual returns. If we assume that there are 252 trading days in a year - that is roughly 365 minus104 weekend days, and about 9/10 public holidays - and if we also assume that the daily changes are independent, we can quickly determine the distribution of the annual changes, as the sum of 252 independent normal variables.

$$S=R_1+R_2+ \dots +R_{252}=\sum_{i=1}^{n=252}{R_i}$$

Mean of this annual change is just the sum of the means, each of which is constant, so it is essentially, 252 times the mean of the daily returns,

$$m_a=E \left[S\right]==252 \, E \left[R_i \right]=252 \, m_d=581.2$$

Variance of the sum of independent variables again is just the sum of their variances. So essentially, we multiply it by 252,

$$v_a=V \left[S\right]==252 V \left[R_i \right]=252 \, v_d=738,912.2$$

So annualising the mean and variance of our daily returns is easy, we just multiply them by 252. Note standard deviation gets multiplied by the square root of 252 because it is equal to the square root of variance.

We have to annualise the parameters because in most cases time is measured in years, so t=1 means 1 year. In this regard the Arithmetic Brownian motion parameters that we had in the SDE can be interpreted as the annual drift, and the annual variance.

Now if t represents time in years,

$$t=1 \Rightarrow \text{1 year}$$

Then 1 day, which we denote by \(\Delta t\), would then be 1 divided by 252 (assuming 252 working days in a year),

$$\Delta t=1/252 \Rightarrow \text{1 Day}$$

And daily mean would be annual mean divided by 252, or annual mean times \(\Delta t\). So to compute annual mean, we divide the daily mean by the \(\Delta t\), which is 1 over 252.

$$m_d=m_a \; \Delta t \Rightarrow m_a=\frac{m_d}{\Delta t}$$

Dividing by 1 over 252 is the same thing as multiplying by 252, so the answer will be the same as before, but the relationship can look different when you encounter it for the first time. We will have the same relationship between the annual and daily variances.

And we have thus calibrated our Arithmetic Brownian motion parameters,

$$\mu=m_a={\rm 581}$$

$$\sigma=\sqrt{v_a}=\sqrt{\rm 738,912.2}=860$$

Thus our SDE is,

$$dX_t={\rm 581} \; dt +{\rm 860} \; dB_t$$

or in terms of the process value (which is the solution of the above SDE),

$$X_T=X_0 +{\rm 581} \; T+{\rm 860} \; B_T$$

Now that we have the parameters of the stochastic process, we can do many useful things with it. For example, we can simulate its paths. Lets recall the solution,

$$X_T=X_0+\mu \, T +\sigma \, B_T$$

It gives the value of the process at time T given the value of the process at time zero.

We can generalise this to an arbitrary start time t, and end time \(t +\Delta t\), meaning we are generating/ observing the value of the process after an interval of length \(\Delta t\),

$$X_{t+\Delta t}=X_t+\mu \, \Delta t +\sigma \, (B_{t+ \Delta t}-B_t)$$

We can then recursively generate/simulate the path of the process. Let’s make it more concrete by substituting for the annual drift and diffusion coefficient that we estimated based on FTSE data,

$$X_{t+\Delta t}=X_t+{\rm 581} \; \Delta t +{\rm 860} \; (B_{t+ \Delta t}-B_t)$$

We can replace the change in Brownian by \(\sqrt{\Delta t}\) times a standard normal,

$$X_{t+\Delta t}=X_t+581 \, \Delta t +860 \, \sqrt{\Delta t} \, N\left[ 0,1\right]$$

This is because the change in Brownian motion over the interval \(\Delta t\) is normally distributed, with mean zero and variance equal to \(\Delta t\). Square root of \(\Delta t\) times standard normal also has a mean of zero, and variance of \(\Delta t\), so the two represent the same thing in the distribution sense.

Now we can use this to simulate the daily values of the index, starting at \(X_0\), and recursively applying the above equation, with \(\Delta t\) equal to 1 over 252. In this case we start with 7,318.5, which was the latest value of the index,

$$X_0=7,318.5$$

We then generate a random value of the standard normal, and then plug it into this equation to get the simulated value one day ahead,

$$X_{\frac{1}{252}}=X_0+581 \, \frac{1}{252} +860 \, \sqrt{\frac{1}{252}} \, N\left[ 0,1\right]$$

You keep recursively generating the values of the process, one day at a time, until we reach our desired maturity/horizon. Assuming we are simulating the process over one year, here is simple Python code to simulate the daily observations:

import numpy as np
dt=1/252
N=252
X=np.zeros(N+1,dtype=np.float)
X[0]=7318.5
for i in range(N):

X[i+1]=X[i]+581*dt+860*np.sqrt(dt)*np.random.standard_normal()

You can then plot the path to visualise how it looks like. This is one way of simulating the process. Recall that the value of the process at any time is normally distributed with mean and variance that depend on the length of the time. So we can also generate the probability distribution of the process at any time using the normal distribution density.

$$f_X (x)=\frac{1}{\sqrt{2\pi v}}e^{-\frac{\left(x-m\right)^2}{2 v}} $$

If you plot the distribution of the process over time, you will see that the distribution drift to the right as the horizon increases. This is because the drift here is positive, so we would expect the process to drift upwards over time. You will also see that the distribution flattens as the horizon increases. This is because the variance increases with time, and you may recall that the higher the variance, the flatter the distribution becomes.