5.6.5 Case study: the lognormal from the normal

The model: a sum of normal log returns

The inputs are those of case studies 2.4.6 and 2.6.5. Today's price is $S_0 = 100$ dollars, the horizon $T = 0.5$ years, the interest rate $r = 0.03$ and the volatility $\sigma = 0.2$ , both annual decimals with $r$ continuously compounded. $Z$ is a standard normal random variable, with the density $\varphi$ that case study 2.6.5 used as a weight.

The log return $Y = \ln(S_T/S_0)$ (case study 1.4.7) is

\begin{gathered} Y = m + vZ \\[4pt] m = \left(r - \tfrac{\sigma^2}{2}\right)T, \qquad v = \sigma\sqrt{T} \end{gathered}

a linear function of $Z$ , so it is normal (lesson 5.6.1), with mean $m = 0.005$ and standard deviation $v \approx 0.141421$ .

Why the log return is normal

Case study 5.4.6 split the six months into 126 trading days and added up the variances of the daily log returns. Lesson 5.6.3 finishes the argument. Each day's log return is independent and $\mathcal{N}\left(\frac{m}{126}, \frac{v^2}{126}\right)$ , and the MGF of their sum is the 126th power of one day's MGF:

\left(e^{\frac{m}{126}t + \frac{v^2}{126}\cdot\frac{t^2}{2}}\right)^{126} = e^{mt + v^2t^2/2}

That is the MGF of $\mathcal{N}(m, v^2)$ , so by uniqueness (lesson 5.6.2) the six-month log return is exactly normal. Its standard deviation is $v = \sigma\sqrt{T}$ : four times the horizon gives twice the spread.

m

import math import numpy as np from scipy.stats import kurtosis, skew S0, r, sigma, T = 100, 0.03, 0.2, 0.5 m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) print(round(m, 6), round(v, 6), round(sigma / math.sqrt(252), 6)) t = 0.7 print(round(math.exp(m / 126 * t + v**2 / 126 * t**2 / 2) ** 126, 10), round(math.exp(m * t + v**2 * t**2 / 2), 10)) days = 126 dt = T / days rng = np.random.default_rng(5651) daily = rng.normal((r - sigma**2 / 2) * dt, sigma * math.sqrt(dt), size=(20_000, days)) Y = daily.sum(axis=1) print(round(float(Y.mean()), 4), round(float(Y.std()), 4)) print(round(float(skew(Y)), 2) + 0.0, round(float(kurtosis(Y)), 2) + 0.0)

⌘/Ctrl + Enter runs the cell

The density of the price

The price $S_T = S_0e^Y$ rises with $Y$ . By the distribution-function method of lesson 5.6.1, $S_T \le s$ exactly when $Y \le \ln(s/S_0)$ , so for $s > 0$

F_{S_T}(s) = N\!\left(\frac{\ln(s/S_0) - m}{v}\right)

with $N$ the standard normal distribution function. For the density, use the monotone rule with $y = \ln(s/S_0)$ and $\frac{dy}{ds} = \frac{1}{s}$ .

Rule: Density of the future price

\begin{gathered} f_{S_T}(s) = \frac{1}{v\,s}\,\varphi\!\left(\frac{\ln(s/S_0) - m}{v}\right) \\[4pt] s > 0 \end{gathered}

The lognormal density of lesson 5.4.4, with the factor $\frac{1}{s}$ from the stretch of $e^y$ . It is 0 for $s \le 0$ : a price $S_0e^Y$ is always positive.

The density is 0.009598 at $80, 0.028192 at $100 and 0.010711 at $120, per dollar. It is not symmetric about $100: the price cannot fall below 0 but can rise without limit, so the right tail is longer.

From $F$ come the probabilities that case study 5.4.6 tabulated. $P(S_T \le 90) \approx 0.217588$ , one minus its 0.782412, and $P(S_T > 100) = N\!\left(\frac{m}{v}\right) \approx 0.514102$ , the $N(d_2)$ of case study 2.6.5.

v

import math from scipy.integrate import quad from scipy.stats import lognorm, norm S0, r, sigma, T = 100, 0.03, 0.2, 0.5 m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) def f(s): return norm.pdf((math.log(s / S0) - m) / v) / (v * s) dist = lognorm(s=v, scale=S0 * math.exp(m)) for s in [80, 100, 120]: print(s, round(f(s), 6), round(dist.pdf(s), 6)) print(round(quad(f, 0, math.inf)[0], 8)) print(round(norm.cdf((math.log(90 / S0) - m) / v), 6), round(norm.cdf(m / v), 6))

⌘/Ctrl + Enter runs the cell

The mean, from the normal's MGF

The mean needs one line with lesson 5.6.2. Since $S_T = S_0e^m e^{vZ}$ , its mean is $S_0e^m\,\mathbb{E}[e^{vZ}] = S_0e^m M_Z(v)$ , and $M_Z(v) = e^{v^2/2}$ :

\begin{aligned} \mathbb{E}[S_T] &= S_0\,e^{m + v^2/2} \\ &= S_0\,e^{(r - \sigma^2/2)T + \sigma^2 T/2} = S_0\,e^{rT} \end{aligned}

That is $101.5113. Case study 2.6.5 found it by completing the square inside an integral. The normal's MGF is that step done once, for every $t$ .

Rule: Mean of a lognormal

\begin{gathered} Y \sim \mathcal{N}(m, v^2) \\[4pt] \mathbb{E}[e^Y] = M_Y(1) = e^{m + v^2/2} \end{gathered}

The mean that lesson 5.4.4 stated: the normal's MGF at $t = 1$ . It is larger than $e^{\mathbb{E}[Y]} = e^m$ by the factor $e^{v^2/2}$ .

The same function gives the spread. $\mathbb{E}[S_T^2] = S_0^2\,\mathbb{E}[e^{2Y}] = S_0^2 M_Y(2) = S_0^2e^{2m + 2v^2} \approx 10512.711$ . Subtracting the squared mean leaves 208.1656, the variance lesson 5.4.4 stated, so the standard deviation is $14.4279, or $S_0e^{rT}\sqrt{e^{\sigma^2 T} - 1}$ .

e^{v^2/2}

import math from scipy.integrate import quad from scipy.stats import norm S0, r, sigma, T = 100, 0.03, 0.2, 0.5 m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) def f(s): return norm.pdf((math.log(s / S0) - m) / v) / (v * s) print(round(math.exp(v**2 / 2), 6), round(math.exp(m), 6), round(S0 * math.exp(m + v**2 / 2), 4)) print(round(quad(lambda s: s * f(s), 0, math.inf)[0], 4), round(S0 * math.exp(r * T), 4)) second = S0**2 * math.exp(2 * m + 2 * v**2) var = second - (S0 * math.exp(r * T))**2 print(round(second, 4), round(var, 4), round(math.sqrt(var), 4), round(S0 * math.exp(r * T) * math.sqrt(math.exp(sigma**2 * T) - 1), 4))

⌘/Ctrl + Enter runs the cell

Why the $-\frac{\sigma^2}{2}$ is there

Prices in case study 2.6.5 are averages in a risk-neutral model: every asset grows on average like money in the bank. So the model must give $\mathbb{E}[S_T] = S_0e^{rT}$ . Then the discounted average, $e^{-rT}\,\mathbb{E}[S_T]$ , is today's price $100.

A log return of $rT + vZ$ , with no correction, would overshoot. The average of $e^{vZ}$ is $e^{v^2/2}$ , above $e^0 = 1$ , though $vZ$ averages 0. The mean would be $S_0e^{rT + \sigma^2T/2} \approx 102.5315$ , which discounts to 101.005: more than the $100 the stock costs today.

Subtracting half the variance from the log drift cancels the factor $e^{v^2/2}$ exactly. Lesson 5.7.1 states the general fact behind it, Jensen's inequality: for a convex $g$ , $\mathbb{E}[g(Y)] \ge g(\mathbb{E}[Y])$ .

The cost of the correction: a median below the mean

Case study 5.4.6 watched the median fall behind the mean as the horizon grew. In the pricing model the mean is pinned at $S_0e^{rT}$ whatever $\sigma$ is, while the median $S_0e^m$ falls as $\sigma$ grows. Six months out:

volatility $\sigma$	0.1	0.2	0.4	0.8
median price	101.26	100.50	97.53	86.50
mean price	101.51	101.51	101.51	101.51
$P(S_T < \text{mean})$	0.5141	0.5282	0.5562	0.6114

At 80% volatility the typical outcome is a loss, to $86.50, while the mean is still $101.51. A few very large outcomes hold the mean up, and 61% of outcomes fall below it.

The last row is $N\!\left(\frac{v}{2}\right)$ , as in lesson 5.4.4. The reason: the log of the mean, $\ln S_0 + m + \frac{v^2}{2}$ , is $\frac{v^2}{2}$ above the median log price, which is $\frac{v}{2}$ standard deviations.

Try it The first line is the model without the correction: its mean and its discounted mean. The loop rebuilds the table.

import math from scipy.stats import norm S0, r, T = 100, 0.03, 0.5 v = 0.2 * math.sqrt(T) wrong = S0 * math.exp(r * T + v**2 / 2) print(round(wrong, 4), round(wrong * math.exp(-r * T), 4)) for sigma in [0.1, 0.2, 0.4, 0.8]: m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) print(sigma, round(S0 * math.exp(m), 2), round(S0 * math.exp(r * T), 2), round(norm.cdf(v / 2), 4))

⌘/Ctrl + Enter runs the cell

The median and the mode

$S_T$ is an increasing function of $Y$ , so quantiles pass straight through (lesson 5.6.1). The median of $Y$ is $m$ , so the median price is $S_0e^m \approx 100.5013$ . In the same way the 5% and 95% quantiles are $S_0e^{m \mp 1.644854v}$ , the $79.64 and $126.82 of case study 5.4.6.

The mode, where the density peaks, is $S_0e^{m - v^2} \approx 98.5112$ ; it is stated here and checked on a grid below. So the mode, the median and the mean come in that order, all within $3.

Five years out they spread apart: the mode is $86.07, the median $105.13 and the mean $116.18. The price then ends below its mean with probability 0.5885.

The density of $S_T$ , six months and five years out. The price cannot fall below 0 but has a long right tail. So the mean sits above the median, the median above the mode, and the gaps grow with $\sigma\sqrt{T}$ .

P(S_T < \text{mean})

import math import numpy as np from scipy.stats import norm S0, r, sigma = 100, 0.03, 0.2 for T, lo, hi in [(0.5, 50, 150), (5, 20, 300)]: m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) print(T, round(S0 * math.exp(m - v**2), 4), round(S0 * math.exp(m), 4), round(S0 * math.exp(r * T), 4), round(norm.cdf(v / 2), 4)) s = np.arange(lo, hi, 0.001) f = norm.pdf((np.log(s / S0) - m) / v) / (v * s) print(round(float(s[np.argmax(f)]), 2))

⌘/Ctrl + Enter runs the cell

A simulation check

Draw $Z$ , compute $S_T = S_0e^{m + vZ}$ , and every claim above can be checked by counting. In one seeded run of 200,000 draws, the mean is $101.53, the median $100.52 and the standard deviation $14.43.

The share at or below $90 is 0.2174, against 0.217588. The share below the mean is 0.5285, against 0.528186. The lowest of the 200,000 prices is $51.62: never negative.

Try it The same seed draws the histogram in the figure above. Case study 5.7.5 prices the call by averaging discounted payoffs over draws like these, and measures how far such an average can stray.

import math import numpy as np S0, r, sigma, T = 100, 0.03, 0.2, 0.5 m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) mean = S0 * math.exp(r * T) rng = np.random.default_rng(5652) S = S0 * np.exp(m + v * rng.standard_normal(200_000)) print(round(float(S.mean()), 2), round(float(np.median(S)), 2), round(float(S.std()), 2)) print(round(float(np.mean(S <= 90)), 4), round(float(np.mean(S < mean)), 4)) print(round(float(S.min()), 2))

⌘/Ctrl + Enter runs the cell

Check yourself

For this stock six months out, which is largest?

The mode of

S_T

The median of

S_T

The mean of

S_T

Check yourself

A colleague simulates

S_T = S_0e^{rT + \sigma\sqrt{T}Z}

, leaving out the

-\frac{\sigma^2}{2}

. What happens?

Nothing:

Z

has mean 0, so the mean price is still

S_0e^{rT}

The mean price rises to

S_0e^{rT + \sigma^2T/2}

, about $102.53, and the stock is worth more than $100 today in the modelThe prices can now be negative

On the desk

Conventions and traps

SciPy's parameters. lognorm(s=v, scale=S0 * exp(m)): the shape is $\sigma\sqrt{T}$ and the scale is the median. With scale=S0, forgetting the drift, the mean comes out $101.005. With s=sigma, forgetting $\sqrt{T}$ , the standard deviation is $20.71 instead of $14.43.
Which drift. Every formula here holds with a forecast drift $\mu$ in place of $r$ , as case study 5.4.6 used. But only $r$ makes the discounted mean today's price: with $\mu = 8\%$ , $e^{-rT}\,\mathbb{E}[S_T]$ is $102.53, a forecast and not a price.
Units. $\sigma$ and $r$ are annual and $T$ is in years. A daily volatility of 0.0126 is $\frac{0.2}{\sqrt{252}}$ , not 0.2 divided by 252.
Thin tails. The normal log return brings the normal's thin tails with it. Case study 5.4.6 found the invented daily returns of Level 5 fatter-tailed, so the model understates large moves.
Edge cases. With $\sigma = 0$ or $T = 0$ the density divides by $v = 0$ . Then $S_T$ is the certain $S_0e^{rT}$ : treat it separately.

\mu = 8\%

import math from scipy.stats import lognorm S0, r, sigma, T = 100, 0.03, 0.2, 0.5 m, v = (r - sigma**2 / 2) * T, sigma * math.sqrt(T) right = lognorm(s=v, scale=S0 * math.exp(m)) print(round(right.mean(), 4), round(lognorm(s=v, scale=S0).mean(), 4)) print(round(right.std(), 4), round(lognorm(s=sigma, scale=S0 * math.exp(m)).std(), 4)) print(round(math.exp(-r * T) * S0 * math.exp(0.08 * T), 4))

⌘/Ctrl + Enter runs the cell

Exercise

The distribution of a future price

Write price_summary(S0, r, sigma, T) and prob_below(level, S0, r, sigma, T) for the model of this case study.

price_summary returns a dict with keys "mean", "median", "mode" and "sd" of $S_T$ : $S_0e^{rT}$ , $S_0e^m$ , $S_0e^{m - v^2}$ and $S_0e^{rT}\sqrt{e^{\sigma^2T} - 1}$ .
prob_below returns $P(S_T \le \text{level})$ , which is $N\!\left(\frac{\ln(\text{level}/S_0) - m}{v}\right)$ . A level of 0 or below gives 0.
With $\sigma = 0$ or $T = 0$ , $S_T$ is certain to be $S_0e^{rT}$ . Then all three centres equal it and sd is 0; prob_below is 1 at or above it and 0 below.
Raise ValueError if $S_0 \le 0$ , $\sigma < 0$ or $T < 0$ .

The tests check the six-month and five-year numbers above, the probabilities 0.217588 and 0.528186, and the edge cases.

See the tests

import math


def test_six_months():
    """Six months out: mean 101.5113, median 100.5013, mode 98.5112, sd 14.4279"""
    got = price_summary(100, 0.03, 0.2, 0.5)
    expect_close(got["mean"], 101.51130646157189, what="mean", rel=1e-9)
    expect_close(got["median"], 100.5012520859401, what="median", rel=1e-9)
    expect_close(got["mode"], 98.51119396030626, what="mode", rel=1e-9)
    expect_close(got["sd"], 14.427945946151569, what="sd", rel=1e-7)


def test_five_years():
    """Five years out: mode 86.07, median 105.13, mean 116.18"""
    got = price_summary(100, 0.03, 0.2, 5)
    expect_close(got["mode"], 100 * math.exp(-0.15), what="mode", rel=1e-9)
    expect_close(got["median"], 100 * math.exp(0.05), what="median", rel=1e-9)
    expect_close(got["mean"], 100 * math.exp(0.15), what="mean", rel=1e-9)


def test_probabilities():
    """P(S_T <= 90) is 0.217588, and P(S_T <= mean) is N(v/2) = 0.528186"""
    expect_close(prob_below(90, 100, 0.03, 0.2, 0.5), 0.21758753431780564, what="prob_below(90)", rel=1e-6)
    expect_close(prob_below(100 * math.exp(0.015), 100, 0.03, 0.2, 0.5), 0.5281859888985083, what="prob_below(mean)", rel=1e-6)
    expect_close(prob_below(100, 100, 0.03, 0.2, 0.5), 1 - 0.514101801652164, what="prob_below(100) = 1 - N(d2)", rel=1e-6)


def test_order():
    """With volatility, mode < median < mean, and the median has probability 1/2 below it"""
    got = price_summary(50, 0.01, 0.35, 2)
    assert got["mode"] < got["median"] < got["mean"], f"expected mode < median < mean, got {got}"
    expect_close(prob_below(got["median"], 50, 0.01, 0.35, 2), 0.5, what="prob_below(median)", rel=1e-9)


def test_no_volatility():
    """With sigma = 0 the price is certain: S0 e^(rT) = 101.5113, sd 0, probabilities 0 or 1"""
    got = price_summary(100, 0.03, 0.0, 0.5)
    for key in ("mean", "median", "mode"):
        expect_close(got[key], 100 * math.exp(0.015), what=key, rel=1e-12)
    expect_close(got["sd"], 0.0, what="sd")
    expect_close(prob_below(101.6, 100, 0.03, 0.0, 0.5), 1.0, what="prob_below just above the certain price")
    expect_close(prob_below(101.4, 100, 0.03, 0.0, 0.5), 0.0, what="prob_below just below the certain price")


def test_at_expiry():
    """With T = 0 the price is today's, 100"""
    got = price_summary(100, 0.03, 0.2, 0)
    expect_close(got["median"], 100.0, what="median at T = 0", rel=1e-12)
    expect_close(got["sd"], 0.0, what="sd at T = 0")


def test_non_positive_level():
    """A level of 0 or below has probability 0"""
    expect_close(prob_below(0, 100, 0.03, 0.2, 0.5), 0.0, what="prob_below(0)")
    expect_close(prob_below(-5, 100, 0.03, 0.2, 0.5), 0.0, what="prob_below(-5)")


def test_bad_inputs():
    """S0 <= 0, a negative volatility or a negative horizon raises ValueError"""
    cases = [
        ("S0 = 0", lambda: price_summary(0, 0.03, 0.2, 0.5)),
        ("sigma = -0.2", lambda: price_summary(100, 0.03, -0.2, 0.5)),
        ("T = -1", lambda: prob_below(90, 100, 0.03, 0.2, -1)),
    ]
    for name, call in cases:
        try:
            call()
        except ValueError:
            continue
        raise AssertionError(f"expected ValueError for {name}")

Z

Case study: the lognormal from the normal

The model: a sum of normal log returns

Why the log return is normal

The density of the price

Rule: Density of the future price

The mean, from the normal's MGF

Rule: Mean of a lognormal

Why the $-\frac{\sigma^2}{2}$ is there

The cost of the correction: a median below the mean

The median and the mode

A simulation check

On the desk

The distribution of a future price

Summary card

Free further reading

The model: a sum of normal log returns

Why the log return is normal

The density of the price

Rule: Density of the future price

The mean, from the normal's MGF

Rule: Mean of a lognormal

Why the −σ22-\frac{\sigma^2}{2}−2σ2​ is there

The cost of the correction: a median below the mean

The median and the mode

A simulation check

On the desk

Summary card

Free further reading

Why the $-\frac{\sigma^2}{2}$ is there