Simulating Stock Prices with Monte Carlo and Brownian Motion in Python

Brownian Motion

Brownian motion refers to the incessant random movement of microscopic particles suspended in a liquid or gas, named after the British botanist Robert Brown. The particles typically have diameters ranging from nanometers to micrometers. Their motion results from collisions with fluid molecules, and the intensity increases with temperature.

Monte Carlo Simulation

The Monte Carlo method comprises a class of computational algorithms that rely on repeated random sampling to obtain numerical results. The basic idea is to use randomness to solve problems that are deterministic in principle. It is commonly applied to physics and mathematics, especially for optimization, numerical integration, and sampling from probability distributions.

Monte Carlo Simulation for Stock Prices

We use Monte Carlo simulation to observe potential asset price changes over time, assuming daily returns follow a normal distribution. This price evolution is also called a random walk, similar to Brownian motion.

Random walk theory suggests stock price changes are random without a fixed pattern. Prices reflect supply‑demand and intrinsic value determined by fundamentals such as assets, P/E ratio, and dividend yield. Market participants collectively keep prices near intrinsic value, causing seemingly arbitrary fluctuations.

Assuming stock prices follow geometric Brownian motion, we derive the asset price dynamics for a given period.

Code

Return Simulation

Generate returns

<code>import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
import pandas as pd
import math
'''
s: current stock price
t: time horizon (years)
r: annualized return
sigma: annualized volatility
nper_per_year: number of periods per year
'''
def generate_simulated_stock_returns(t, r, sigma, nper_per_year):
    simulated_returns = []  # list of simulated returns
    dt = 1 / nper_per_year  # time step
    term = int(t * nper_per_year)  # total steps
    for i in range(1, term + 1):  # loop
        z = np.random.normal()  # random number
        simulated_return = (r - (sigma**2 / 2)) * dt + z * sigma * (dt**(1/2))  # return each step
        simulated_returns.append(simulated_return)
    array_return = np.array(simulated_returns)
    return array_return  # return list of returns
</code>

Return simulation example

<code># initial stock price s=100; expected return r=0.1; sigma=0.3
s=100; r=0.1; sigma=0.3
# 1‑year, 2 periods per year
t=1; nper_per_year=2
array_return = generate_simulated_stock_returns(t, r, sigma, nper_per_year)
print(array_return)
# 2‑year, 24 periods per year
t=2; nper_per_year=24
array_return = generate_simulated_stock_returns(t, r, sigma, nper_per_year)
print(array_return)
</code>

Simulation results:

<code>[ 0.2036485  -0.03690405]
[-5.67879573e-02 -1.50695772e-01 ...]
</code>

Stock Price Simulation

Generate stock prices

<code>def generate_simulated_stock_values(s, t, r, sigma, nper_per_year):
    rate = generate_simulated_stock_returns(t, r, sigma, nper_per_year)
    stock_price = [s]
    term = int(t * nper_per_year)
    for i in range(1, term + 1):
        values = stock_price[i-1] * math.e**(rate[i-1])
        stock_price.append(values)
    array_price = np.array(stock_price)
    return array_price
</code>

Stock price simulation example

<code># 1‑year, 2 periods per year
t=1; nper_per_year=2
array_price = generate_simulated_stock_values(s, t, r, sigma, nper_per_year)
print(array_price)
# 2‑year, 24 periods per year
t=2; nper_per_year=24
array_price = generate_simulated_stock_values(s, t, r, sigma, nper_per_year)
print(array_price)
</code>

Plotting Stock Prices

Plotting function

<code>def plot_simulated_stock_values(s, t, r, sigma, nper_per_year, num_trials=1):
    term = int(t * nper_per_year) + 1
    x_axis = np.linspace(0, t, term)
    for i in range(num_trials):
        price = generate_simulated_stock_values(s, t, r, sigma, nper_per_year)
        plt.plot(x_axis, price)
    plt.title(str(num_trials) + " simulated trials")
    plt.xlabel("years")
    plt.ylabel("value")
    plt.show()
</code>

Simulate 5 trials (2‑year, 250 periods per year)

<code># 2‑year, 250 periods per year, simulate 5 times
t=2; nper_per_year=250; num_trials=5
plot_simulated_stock_values(s, t, r, sigma, nper_per_year, num_trials)
</code>

Simulate 1000 trials (2‑year, 250 periods per year)

<code># 2‑year, 250 periods per year, simulate 10000 times
t=2; nper_per_year=250; num_trials=10000
plot_simulated_stock_values(s, t, r, sigma, nper_per_year, num_trials)
</code>

References

Python Quantitative: Monte Carlo Simulation to Predict Stock Trends https://mp.weixin.qq.com/s/43KQgH-BArop29uJ3Fe7MA