Conditional Expectation and the Tower Property

Companion notebook to the Conditional Expectation and the Tower Property article (Probability Theory, Module 3). This notebook verifies every quantitative claim made in the article using pure Python — no external libraries are required.

§0 — Imports and shared utilities
§1 — Discrete CE from scratch: defining-property verification on Ω = {a,b,c,d}
§2 — Tower property on a 3-level filtration (Ω = {1,…,8})
§3 — Pulling-out-known-factors property
§4 — Jensen's inequality for conditional expectation
§5 — Bivariate Gaussian CE formula + validation summary

💡 Click any code cell to edit it. Press Shift+Enter or click ▶ Run to execute. Your edits are auto-saved in your browser. Click ↻ Reset to restore the original code.

Cells §1–§4 use fractions.Fraction for exact rational arithmetic. Cell §5 uses math floats (continuous distribution). Run all cells in order.

§0 — Imports and shared utilities

Standard library only. fractions.Fraction gives exact rational arithmetic for the discrete cells. math provides standard math functions for the Gaussian cell.

from __future__ import annotations
from fractions import Fraction
from itertools import combinations
import math

# Type alias: an outcome is a string label; an atom is a frozenset of outcomes.
Outcome = str
Atom = frozenset

def uniform_prob(omega: frozenset) -> dict[str, Fraction]:
    """Return uniform probability {outcome: 1/|omega|} using exact Fractions."""
    n = len(omega)
    return {w: Fraction(1, n) for w in omega}

def expectation(X: dict[str, Fraction], prob: dict[str, Fraction]) -> Fraction:
    """Compute E[X] = sum_w X(w) * P({w})."""
    return sum(X[w] * prob[w] for w in X)

def cond_expectation_on_partition(
    X: dict[str, Fraction],
    prob: dict[str, Fraction],
    partition: list[frozenset],
) -> dict[str, Fraction]:
    """
    Compute E[X | G] where G is the sigma-algebra generated by `partition`.
    On each atom A in the partition, E[X|G](omega) = E[X * 1_A] / P(A)
    for all omega in A (constant on atoms).
    Returns a dict omega -> conditional expectation value.
    """
    result: dict[str, Fraction] = {}
    for atom in partition:
        p_atom = sum(prob[w] for w in atom)
        if p_atom == 0:
            raise ValueError(f"Atom {set(atom)} has probability 0.")
        e_x_on_atom = sum(X[w] * prob[w] for w in atom)
        # E[X|G](omega) = E[X * 1_atom] / P(atom)  for omega in atom
        ce_value = e_x_on_atom / p_atom
        for w in atom:
            result[w] = ce_value
    return result

print("✓ Imports and utilities ready.")
print("  - uniform_prob: uniform measure on finite Omega")
print("  - expectation: E[X] via sum formula")
print("  - cond_expectation_on_partition: E[X|G] for partition-generated G")

§1 — Discrete conditional expectation: defining-property verification

The article's hand exercise: Ω = {a, b, c, d}, P uniform (each weight 1/4), G = σ({a,b}, {c,d}), X(a)=1, X(b)=3, X(c)=0, X(d)=4.

We compute E[X|G] from the defining property: on atom {a,b}, the conditional expectation is constant and equals the average of X over {a,b} under P. We then verify the integral identity (CE) on each G-atom.

✦ Key point: E[X|G] is a random variable (a function on Ω), not a single number. It takes the value 2 on both {a,b} and {c,d} here — an illustrative coincidence; in general the values on different atoms differ.

# ── Sample space and measure ─────────────────────────────────────────────────
omega4 = frozenset({"a", "b", "c", "d"})
prob4  = uniform_prob(omega4)   # each outcome: 1/4

# ── Random variable X ────────────────────────────────────────────────────────
X4 = {"a": Fraction(1), "b": Fraction(3), "c": Fraction(0), "d": Fraction(4)}

# ── Sub-sigma-algebra G = sigma({a,b}, {c,d}) ────────────────────────────────
# Atoms of G (the minimal non-empty measurable sets):
atom_ab = frozenset({"a", "b"})
atom_cd = frozenset({"c", "d"})
partition4 = [atom_ab, atom_cd]

# ── Compute E[X|G] ────────────────────────────────────────────────────────────
ce4 = cond_expectation_on_partition(X4, prob4, partition4)

print("Conditional expectation E[X | G]:")
for w in sorted(omega4):
    print(f"  E[X|G]({w}) = {ce4[w]}  (X({w}) = {X4[w]})")

# ── Verify the defining property (CE) on each atom ───────────────────────────
# (CE) requires: integral of E[X|G] over G == integral of X over G, for all G in sigma-algebra
print()
print("Verification of defining property (CE):")

for atom, label in [(atom_ab, "{a,b}"), (atom_cd, "{c,d}")]:
    lhs = sum(ce4[w] * prob4[w] for w in atom)   # int_G E[X|G] dP
    rhs = sum(X4[w]  * prob4[w] for w in atom)   # int_G X dP
    match = "✓" if lhs == rhs else "✗"
    print(f"  G = {label}: int_G E[X|G] dP = {lhs},  int_G X dP = {rhs}  {match}")

# Whole omega is also in G (trivially)
lhs_omega = sum(ce4[w] * prob4[w] for w in omega4)
rhs_omega  = sum(X4[w]  * prob4[w] for w in omega4)
match_omega = "✓" if lhs_omega == rhs_omega else "✗"
print(f"  G = Omega:  int_Omega E[X|G] dP = {lhs_omega},  int_Omega X dP = {rhs_omega}  {match_omega}")

§2 — Tower property on a 3-level filtration

We verify: if H ⊆ G ⊆ F, then E[E[X|G]|H] = E[X|H] a.s.

Setup: Ω = {1,2,3,4,5,6,7,8}, P uniform (each weight 1/8).

H = σ({1,2,3,4}, {5,6,7,8}) — coarsest non-trivial level (2 atoms)
G = σ({1,2}, {3,4}, {5,6}, {7,8}) — finer level (4 atoms)
F = 2^Ω — full information (8 atoms, one per outcome)

We check that conditioning E[X|G] further on H gives the same result as conditioning X on H directly.

# ── Omega = {1,...,8}, uniform ────────────────────────────────────────────────
omega8 = frozenset(range(1, 9))
prob8  = uniform_prob(omega8)   # each outcome: 1/8

# ── Random variable: X(k) = k^2 (arbitrary, non-constant across atoms) ────────
X8 = {k: Fraction(k * k) for k in range(1, 9)}

# ── H = sigma({1,2,3,4}, {5,6,7,8})  — 2 atoms ───────────────────────────────
H_partition = [frozenset(range(1, 5)), frozenset(range(5, 9))]

# ── G = sigma({1,2},{3,4},{5,6},{7,8})  — 4 atoms ────────────────────────────
G_partition = [
    frozenset({1, 2}), frozenset({3, 4}),
    frozenset({5, 6}), frozenset({7, 8}),
]

# ── Compute E[X|G] then E[E[X|G]|H] ─────────────────────────────────────────
ceG  = cond_expectation_on_partition(X8, prob8, G_partition)
# Now treat E[X|G] as a new random variable and condition it on H
ceGH = cond_expectation_on_partition(ceG, prob8, H_partition)

# ── Compute E[X|H] directly ───────────────────────────────────────────────────
ceH  = cond_expectation_on_partition(X8, prob8, H_partition)

# ── Compare ───────────────────────────────────────────────────────────────────
print("Tower property verification: E[E[X|G]|H] == E[X|H]?")
print()
print(f"  X(k) = k^2 for k in {{1,...,8}}")
print()
all_match = True
for k in sorted(omega8):
    lhs = ceGH[k]
    rhs = ceH[k]
    match = "✓" if lhs == rhs else "✗"
    if lhs != rhs:
        all_match = False
    print(f"  k={k}: E[E[X|G]|H]={lhs}  E[X|H]={rhs}  {match}")

print()
print(f"  Tower property holds for all outcomes: {'✓ YES' if all_match else '✗ NO'}")

# Also verify that conditioning E[X|H] on G returns E[X|H] (since H ⊆ G)
ceHG = cond_expectation_on_partition(ceH, prob8, G_partition)
all_match2 = all(ceHG[k] == ceH[k] for k in omega8)
print(f"  E[E[X|H]|G] == E[X|H] (since E[X|H] is G-measurable): {'✓ YES' if all_match2 else '✗ NO'}")

§3 — Pulling out known factors: E[XY|G] = Y · E[X|G]

If Y is G-measurable (constant on each G-atom), then Y is "already known" given G. The theorem states that E[XY|G] = Y · E[X|G] a.s.

We verify this on Ω = {a,b,c,d} from §1, using Y defined to be G-measurable: Y takes one value on {a,b} and another on {c,d}.

✦ Finance application: In risk-neutral pricing, the discount factor e^-r(T-t) is known at time t (it is F_t-measurable), so we can pull it out: E[e^-r(T-t)·Φ(S_T)|F_t] = e^-r(T-t)·E[Φ(S_T)|F_t]. This is used constantly to move discount factors in and out of expectations.

# Reuse Omega = {a,b,c,d}, prob4, X4, partition4 from §1
# Y is G-measurable: constant on each atom
# Y({a,b}) = 2, Y({c,d}) = 5
Y4 = {"a": Fraction(2), "b": Fraction(2), "c": Fraction(5), "d": Fraction(5)}

# XY pointwise product
XY4 = {w: X4[w] * Y4[w] for w in omega4}

# E[XY | G]  — direct computation
ceXY = cond_expectation_on_partition(XY4, prob4, partition4)

# Y * E[X|G]  — pulling-out-known-factors formula
# E[X|G] was computed as ce4 in §1 — recompute here to be self-contained
ce4_local = cond_expectation_on_partition(X4, prob4, partition4)
Y_times_ceX = {w: Y4[w] * ce4_local[w] for w in omega4}

print("Pulling-out-known-factors: E[XY|G] vs Y * E[X|G]")
print()
print(f"  X:  {' '.join(f'{w}={X4[w]}' for w in sorted(omega4))}")
print(f"  Y:  {' '.join(f'{w}={Y4[w]}' for w in sorted(omega4))}  (G-measurable: const on atoms)")
print(f"  XY: {' '.join(f'{w}={XY4[w]}' for w in sorted(omega4))}")
print()

all_match = True
for w in sorted(omega4):
    lhs = ceXY[w]
    rhs = Y_times_ceX[w]
    match = "✓" if lhs == rhs else "✗"
    if lhs != rhs:
        all_match = False
    print(f"  {w}: E[XY|G] = {lhs},   Y*E[X|G] = {rhs}   {match}")

print()
print(f"  Property E[XY|G] = Y·E[X|G] holds for all outcomes: {'✓ YES' if all_match else '✗ NO'}")

# Contrast: if Y is NOT G-measurable, the equality breaks
Y_nonG = {"a": Fraction(1), "b": Fraction(3), "c": Fraction(1), "d": Fraction(3)}
XY_nonG = {w: X4[w] * Y_nonG[w] for w in omega4}
ceXY_nonG = cond_expectation_on_partition(XY_nonG, prob4, partition4)
ce_X_local = cond_expectation_on_partition(X4, prob4, partition4)
Y_nonG_times_ceX = {w: Y_nonG[w] * ce_X_local[w] for w in omega4}
breaks = any(ceXY_nonG[w] != Y_nonG_times_ceX[w] for w in omega4)
print(f"\n  With Y NOT G-measurable: property fails = {'✓ YES (expected)' if breaks else 'NO (unexpected)'}")

§4 — Jensen's inequality: φ(E[X|G]) ≤ E[φ(X)|G]

For a convex function φ and X ∈ L¹ with φ(X) ∈ L¹:

φ(𝔼[X|G]) ≤ 𝔼[φ(X)|G] a.s.

We verify this for φ(x) = x² (convex) on Ω = {a,b,c,d}. Note: for linear φ, equality holds — verifying that equality in Jensen characterises linearity.

⚠ Finance implication: By Jensen, for a convex payoff Φ: E[Φ(S_T)|F_t] ≥ Φ(E[S_T|F_t]). This is the mathematical reason why an at-the-money call option (convex payoff) has positive time value — the expectation of a convex function exceeds the convex function of the expectation. Options have value precisely because of convexity.

# Reuse Omega = {a,b,c,d}, prob4, X4, partition4 from §1
# phi(x) = x^2  (convex)
phi = lambda x: x * x

# phi(X): apply phi pointwise
phi_X4 = {w: phi(X4[w]) for w in omega4}

# E[X|G]  (the conditional expectation of X)
ceX = cond_expectation_on_partition(X4, prob4, partition4)

# phi(E[X|G]): apply phi to conditional expectation
phi_ceX = {w: phi(ceX[w]) for w in omega4}

# E[phi(X)|G]: conditional expectation of phi(X)
ce_phi_X = cond_expectation_on_partition(phi_X4, prob4, partition4)

print("Jensen's inequality: phi(E[X|G]) <= E[phi(X)|G]   for phi(x) = x^2")
print()
print(f"  X values:   { {w: float(X4[w]) for w in sorted(omega4)} }")
print(f"  phi(X):     { {w: float(phi_X4[w]) for w in sorted(omega4)} }")
print()

all_ok = True
for w in sorted(omega4):
    lhs = phi_ceX[w]         # phi(E[X|G](w))
    rhs = ce_phi_X[w]        # E[phi(X)|G](w)
    ok  = lhs <= rhs
    if not ok:
        all_ok = False
    diff = rhs - lhs
    print(f"  {w}: phi(E[X|G]) = {lhs},   E[phi(X)|G] = {rhs},   gap = {diff}   {'✓' if ok else '✗'}")

print()
print(f"  Jensen's inequality satisfied for all outcomes: {'✓ YES' if all_ok else '✗ NO'}")

# Check: equality holds on atom {a,b} because X takes values 1 and 3 (non-constant),
# strict inequality is expected. Confirm.
print()
print("  Atom {a,b}: E[X|G] = 2, phi(E[X|G]) = 4")
print(f"  E[phi(X)|G] on {{a,b}} = (phi(1)+phi(3))/2 = (1+9)/2 = 5  =>  4 < 5  (strict inequality)")

# Also verify equality for linear phi (should be exact equality)
phi_linear = lambda x: x * Fraction(3) + Fraction(1)
phi_lin_X4 = {w: phi_linear(X4[w]) for w in omega4}
ceX_l = cond_expectation_on_partition(X4, prob4, partition4)
phi_ceX_l = {w: phi_linear(ceX_l[w]) for w in omega4}
ce_phi_X_l = cond_expectation_on_partition(phi_lin_X4, prob4, partition4)
equality = all(phi_ceX_l[w] == ce_phi_X_l[w] for w in omega4)
print(f"\n  Linear phi(x)=3x+1 gives equality phi(E[X|G])=E[phi(X)|G]: {'✓ YES' if equality else '✗ NO'}")

§5 — Bivariate Gaussian CE formula and validation summary

The article derives: for jointly Gaussian (X, Y) with parameters (μ_X, μ_Y, σ_X, σ_Y, ρ):

𝔼[X | Y = y] = μ_X + ρ · (σ_X / σ_Y) · (y − μ_Y)

We implement this formula and verify it numerically against the regression coefficient formula α = Cov(X,Y)/Var(Y). Cell uses math module floats (no numpy).

✦ Key structural insight: The conditional expectation is linear in y for the Gaussian case — it is the best linear predictor. For non-Gaussian distributions, E[X|Y=y] may be non-linear, and the linear projection (regression) is only an approximation. In Monte Carlo LSMC, we exploit this by using polynomial regression — which is exact for Gaussian state variables and approximate otherwise.

# ── Bivariate Gaussian CE formula ─────────────────────────────────────────────
def gaussian_ce(y: float, mu_X: float, mu_Y: float,
                sigma_X: float, sigma_Y: float, rho: float) -> float:
    """
    E[X | Y = y] for jointly Gaussian (X, Y).
    Formula: mu_X + rho * (sigma_X / sigma_Y) * (y - mu_Y)
    Assumptions: sigma_X > 0, sigma_Y > 0, rho in (-1, 1).
    """
    if sigma_X <= 0 or sigma_Y <= 0:
        raise ValueError("Standard deviations must be positive.")
    if abs(rho) >= 1:
        raise ValueError("Correlation rho must be in (-1, 1).")
    alpha = rho * (sigma_X / sigma_Y)
    return mu_X + alpha * (y - mu_Y)

# ── Test case 1: standard bivariate (mu_X=mu_Y=0, sigma_X=sigma_Y=1) ─────────
mu_X, mu_Y = 0.0, 0.0
sigma_X, sigma_Y = 1.0, 1.0
rho = 0.6

print("Test 1: mu_X=mu_Y=0, sigma_X=sigma_Y=1, rho=0.6")
print("  E[X|Y=y] = rho * y = 0.6 * y  (symmetric standardised case)")
for y_val in [-2.0, -1.0, 0.0, 1.0, 2.0]:
    result = gaussian_ce(y_val, mu_X, mu_Y, sigma_X, sigma_Y, rho)
    expected = rho * y_val
    diff = abs(result - expected)
    print(f"  y={y_val:+.1f}: E[X|Y=y] = {result:.4f},  rho*y = {expected:.4f},  error = {diff:.2e}")

# ── Test case 2: shifted parameters ─────────────────────────────────────────
print()
mu_X2, mu_Y2 = 2.0, 5.0
sigma_X2, sigma_Y2 = 3.0, 2.0
rho2 = -0.4
alpha2 = rho2 * sigma_X2 / sigma_Y2

print(f"Test 2: mu_X={mu_X2}, mu_Y={mu_Y2}, sigma_X={sigma_X2}, sigma_Y={sigma_Y2}, rho={rho2}")
print(f"  alpha = rho * sigma_X / sigma_Y = {rho2} * {sigma_X2}/{sigma_Y2} = {alpha2:.4f}")
print(f"  E[X|Y=y] = {mu_X2} + {alpha2:.4f} * (y - {mu_Y2})")
for y_val in [3.0, 5.0, 7.0]:
    result = gaussian_ce(y_val, mu_X2, mu_Y2, sigma_X2, sigma_Y2, rho2)
    manual  = mu_X2 + alpha2 * (y_val - mu_Y2)
    diff    = abs(result - manual)
    print(f"  y={y_val:.1f}: formula={result:.6f}, manual={manual:.6f}, error={diff:.2e}")

# ── Verify regression coefficient formula Cov(X,Y)/Var(Y) ────────────────────
# For bivariate Gaussian: Cov(X,Y) = rho * sigma_X * sigma_Y
# Var(Y) = sigma_Y^2
# alpha = Cov/Var = rho * sigma_X * sigma_Y / sigma_Y^2 = rho * sigma_X / sigma_Y
cov_XY = rho2 * sigma_X2 * sigma_Y2
var_Y  = sigma_Y2 ** 2
alpha_regress = cov_XY / var_Y
print()
print(f"Regression formula check:")
print(f"  Cov(X,Y) = rho*sigma_X*sigma_Y = {cov_XY:.4f}")
print(f"  Var(Y)   = sigma_Y^2 = {var_Y:.4f}")
print(f"  alpha = Cov(X,Y)/Var(Y) = {alpha_regress:.6f}")
print(f"  Direct formula: rho*sigma_X/sigma_Y = {alpha2:.6f}")
print(f"  Match: {'✓' if abs(alpha_regress - alpha2) < 1e-12 else '✗'}")

# ── FULL VALIDATION SUMMARY ────────────────────────────────────────────────────
print()
print("=" * 60)
print("  VALIDATION SUMMARY — Module 3: Conditional Expectation")
print("=" * 60)

checks = []

# Check 1: Discrete CE defining property on {a,b,c,d}
from fractions import Fraction as F
prob4c  = {"a": F(1,4), "b": F(1,4), "c": F(1,4), "d": F(1,4)}
X4c     = {"a": F(1), "b": F(3), "c": F(0), "d": F(4)}
part4c  = [frozenset({"a","b"}), frozenset({"c","d"})]
ce4c    = cond_expectation_on_partition(X4c, prob4c, part4c)
def_ok  = all(
    sum(ce4c[w]*prob4c[w] for w in atom) == sum(X4c[w]*prob4c[w] for w in atom)
    for atom in part4c
)
checks.append(("Discrete CE: defining property (CE) holds on both atoms", def_ok))

# Check 2: Tower property on {1,...,8}
prob8c = {k: F(1,8) for k in range(1,9)}
X8c    = {k: F(k*k) for k in range(1,9)}
Hpart  = [frozenset(range(1,5)), frozenset(range(5,9))]
Gpart  = [frozenset({1,2}), frozenset({3,4}), frozenset({5,6}), frozenset({7,8})]
ceGc   = cond_expectation_on_partition(X8c, prob8c, Gpart)
ceGHc  = cond_expectation_on_partition(ceGc, prob8c, Hpart)
ceHc   = cond_expectation_on_partition(X8c, prob8c, Hpart)
tower_ok = all(ceGHc[k] == ceHc[k] for k in range(1,9))
checks.append(("Tower property: E[E[X|G]|H] = E[X|H] on Omega={1..8}", tower_ok))

# Check 3: Pulling out known factors
Y4c    = {"a": F(2), "b": F(2), "c": F(5), "d": F(5)}
XY4c   = {w: X4c[w]*Y4c[w] for w in ["a","b","c","d"]}
ceXYc  = cond_expectation_on_partition(XY4c, prob4c, part4c)
ceXc   = cond_expectation_on_partition(X4c,  prob4c, part4c)
pof_ok = all(ceXYc[w] == Y4c[w]*ceXc[w] for w in ["a","b","c","d"])
checks.append(("Pulling-out-known-factors: E[XY|G]=Y*E[X|G] (G-meas Y)", pof_ok))

# Check 4: Jensen's inequality phi(E[X|G]) <= E[phi(X)|G] for phi=x^2
phi_X4c   = {w: X4c[w]**2 for w in ["a","b","c","d"]}
phi_ceXc  = {w: ceXc[w]**2 for w in ["a","b","c","d"]}
ce_phi_Xc = cond_expectation_on_partition(phi_X4c, prob4c, part4c)
jensen_ok = all(phi_ceXc[w] <= ce_phi_Xc[w] for w in ["a","b","c","d"])
checks.append(("Jensen's inequality: phi(E[X|G]) <= E[phi(X)|G] for phi=x^2", jensen_ok))

# Check 5: Bivariate Gaussian CE formula
def _ce_ok(mu_X, mu_Y, sx, sy, rho):
    for y in [-1.5, 0.0, 2.3]:
        expected = mu_X + rho*(sx/sy)*(y - mu_Y)
        got      = gaussian_ce(y, mu_X, mu_Y, sx, sy, rho)
        if abs(got - expected) > 1e-12:
            return False
    return True
gauss_ok = _ce_ok(0,0,1,1,0.6) and _ce_ok(2,5,3,2,-0.4)
checks.append(("Bivariate Gaussian CE formula matches closed-form", gauss_ok))

all_pass = True
for label, passed in checks:
    icon   = "✓" if passed else "✗"
    status = "PASS" if passed else "FAIL"
    print(f"  {icon} [{status}]  {label}")
    if not passed:
        all_pass = False
print("=" * 60)
print(f"\n  {'All checks passed ✓' if all_pass else 'Some checks FAILED — review above'}")