Geometrische Brownsche Bewegung

Preise ändern sich in Prozent

Die arithmetische Brownsche Bewegung addiert dasselbe Dollar-Rauschen, ob die Aktie 5 oder 500 kostet, und kann unter null fallen. Eine Aktie verhält sich anders: Eine Bewegung von 1 % ist bei jedem Preis 1 %. Die geometrische Brownsche Bewegung legt den Zufall in die prozentuale Änderung: dS = mu * S * dt + sigma * S * dW. Sie ist das Modell von Black–Scholes. mu ist die erwartete Jahresrendite und sigma die Volatilität, die am Markt gestellte Größe.

import math

def gbm_step(s, mu, sigma, dt, z):
    return s * math.exp((mu - 0.5 * sigma ** 2) * dt + sigma * math.sqrt(dt) * z)

print(round(gbm_step(100.0, 0.08, 0.2, 1 / 252, 0.0), 4))     # a calm day: only the (mu - sigma^2/2) drift acts
print(round(gbm_step(100.0, 0.08, 0.2, 1 / 252, 1.0), 4))     # a one-sigma up day
print(round(gbm_step(100.0, 0.08, 0.2, 1 / 252, -50.0), 4))   # an absurd fall: still positive

Die exakte Lösung und das minus ein halbes Sigma zum Quadrat

Die Lösung der Gleichung lautet S(t) = S(0) * exp((mu - sigma**2 / 2) * t + sigma * W(t)). Warum das seltsame - sigma**2 / 2? Weil die Exponentialfunktion konvex ist: Nach der Jensenschen Ungleichung hebt Rauschen den Mittelwert von exp(Rauschen) an, und das Abziehen von sigma**2 / 2 hebt das genau auf, sodass der *erwartete* Preis mit Rate mu wächst. Typische *Pfade* wachsen langsamer, mit mu - sigma**2 / 2, weshalb der Median des Aktienkurses unter dem Mittel liegt. Das nächste Thema, das Itô-Lemma, leitet diesen Term her.

import numpy as np

def gbm_path(s0, mu, sigma, dt, z):
    log_moves = (mu - 0.5 * sigma ** 2) * dt + sigma * np.sqrt(dt) * z
    return s0 * np.exp(np.concatenate([[0.0], np.cumsum(log_moves)]))

rng = np.random.default_rng(1)
z = rng.standard_normal(252)
path = gbm_path(100.0, 0.08, 0.2, 1 / 252, z)
print(path[:4].round(3), round(float(path[-1]), 3))

Viele Pfade: Mittelwert und Median

Simulieren Sie Tausende Einjahrespfade und sehen Sie, wo sie enden. Der *mittlere* Endpreis liegt nahe S0 * exp(mu * T), der Formel für den erwarteten Preis. Der *Median* ist niedriger, S0 * exp((mu - sigma**2/2) * T), und bei hoher Volatilität ist der Abstand groß: Die meisten Pfade enden unter dem Mittel, das einige spektakuläre hochziehen. Diese Rechtsschiefe, eine Lognormalverteilung, ist die Form von Aktienkursen über lange Horizonte.

import numpy as np

rng = np.random.default_rng(5)
s0, mu, sigma, T, steps, paths = 100.0, 0.08, 0.4, 1.0, 252, 20000
dt = T / steps
z = rng.standard_normal((paths, steps))
final = s0 * np.exp(np.sum((mu - 0.5 * sigma ** 2) * dt + sigma * np.sqrt(dt) * z, axis=1))

print(round(float(final.mean()), 1), round(s0 * np.exp(mu * T), 1))                      # close to each other: the expected price
print(round(float(np.median(final)), 1), round(s0 * np.exp((mu - 0.5 * sigma ** 2) * T), 1))   # the median is lower
print(round(float(np.mean(final < s0)), 3))                                              # the share of paths that lose money

mu und sigma schätzen

Aus alle dt Jahre erfassten Preisen sind die Log-Renditen r = diff(log(prices)) unabhängige Normalverteilungen mit Mittel (mu - sigma**2/2) * dt und Standardabweichung sigma * sqrt(dt). Also sigma = std(r) / sqrt(dt) und mu = mean(r) / dt + sigma**2 / 2. Die Volatilitätsschätzung ist mit einigen hundert Punkten gut; der Drift ist notorisch verrauscht, weil sein Standardfehler nur mit der gesamten *Zeit* schrumpft, nicht mit der Zahl der Beobachtungen.

import numpy as np

def estimate_gbm(prices, dt):
    r = np.diff(np.log(prices))
    sigma = r.std(ddof=1) / np.sqrt(dt)
    return r.mean() / dt + 0.5 * sigma ** 2, sigma

rng = np.random.default_rng(8)
dt = 1 / 252
for years in (1, 10, 100):
    n = int(years * 252)
    log_moves = (0.10 - 0.5 * 0.25 ** 2) * dt + 0.25 * np.sqrt(dt) * rng.standard_normal(n)
    prices = 100 * np.exp(np.concatenate([[0.0], np.cumsum(log_moves)]))
    mu, sigma = estimate_gbm(prices, dt)
    print(years, "years:", round(float(mu), 3), round(float(sigma), 3))           # true mu 0.10 and sigma 0.25

• GBM: S(t) = S(0) * exp((mu - sigma**2 / 2) * t + sigma * W(t)); prozentuale Änderungen sind zufällig, Preise bleiben positiv.
• Das - sigma**2 / 2 hebt die Konvexität von exp auf: Der erwartete Preis ist S0 * exp(mu * t), der Median-Pfad wächst aber mit mu - sigma**2 / 2.
• Schätzen Sie aus Log-Renditen: sigma = std(r) / sqrt(dt), mu = mean(r) / dt + sigma**2 / 2; Volatilität lässt sich gut schätzen, Drift schlecht.

Als Nächstes: das Itô-Lemma, die Kettenregel, die erklärt, woher das - sigma**2 / 2 kommt.

Lab: Was bewertet wird

• Mit gbm_step einen exakten Schritt machen25 Pkt.
• Mit gbm_path einen Kurspfad simulieren25 Pkt.
• Mit gbm_expected über alle Zukünfte mitteln25 Pkt.
• Mit estimate_gbm die Parameter zurücklesen25 Pkt.