Stochastic Differential Equations: Black–Scholes, Heat Kernels, and Mean Reversion

A stochastic differential equation (SDE) specifies how a state changes through a predictable drift and unpredictable shocks. Solving it gives an entire random trajectory, rather than one deterministic curve. A related partial differential equation (PDE) can describe the distribution of that trajectory or the value of a claim on it.

Learning goals

• Derive the Black–Scholes PDE using a self-financing hedge.
• Connect a terminal payoff to Gaussian smoothing through the heat equation.
• State sufficient conditions for a well-defined, unique SDE solution.
• Solve geometric Brownian motion and the Ornstein–Uhlenbeck process.
• Distinguish expected growth, typical path growth, and stationary distributions.
• Interpret a matrix of diffusion coefficients in a multidimensional system.

1. Reading an SDE correctly

On a filtered probability space, let $B_t$ be standard Brownian motion relative to the information filtration $\mathcal F_t$. The equation

$$dX_t=b(t,X_t)\,dt+a(t,X_t)\,dB_t$$

means the integral identity

$$X_t=X_0+\int_0^t b(s,X_s)\,ds+\int_0^t a(s,X_s)\,dB_s.$$

The last term is an Itô integral, not an ordinary derivative multiplied by a small number. Locally, conditional on today’s state, the drift contributes approximately $b\,\Delta t$ and the shock has variance approximately $a^2\Delta t$. Brownian shocks scale like $\sqrt{\Delta t}$.

The process must use information already available. Future Brownian increments cannot be used to choose today’s coefficient or position.

2. Delta hedging and the Black–Scholes PDE

Assume a non-dividend-paying asset follows

$$dP_t=\mu P_t\,dt+\sigma P_t\,dB_t,$$

with constant volatility $\sigma>0$, constant risk-free rate $r$, frictionless continuous trading, and borrowing and lending at that rate. For a sufficiently smooth derivative value $G(t,P)$, Itô’s formula gives

$$dG=\left(G_t+\mu P G_P+\frac12\sigma^2P^2G_{PP}\right)dt+\sigma P G_P\,dB_t.$$

Here $G_t$ denotes the time partial derivative; $G(t,P_t)$ denotes the value process.

Use a self-financing replicating strategy

Let a portfolio hold $\Delta_t$ shares and a cash balance $C_t$, with value $G=\Delta P+C$. Self-financing means that rebalancing is paid for from the portfolio:

$$dG=\Delta\,dP+rC\,dt=\Delta\,dP+r(G-\Delta P)dt.$$

Matching the Brownian coefficients requires $\Delta=G_P$. Matching the drift terms then yields

$$G_t+\mu PG_P+\frac12\sigma^2P^2G_{PP}=\mu PG_P+r(G-PG_P).$$

Therefore

$$\boxed{G_t+rPG_P+\frac12\sigma^2P^2G_{PP}-rG=0.}$$

The real-world drift $\mu$ cancels because the hedge removes exposure to the asset’s instantaneous random move.

Equivalently, short one derivative and hold $G_P$ shares. Over an instant with the holdings fixed, the hedged gain is

$$-dG+G_P\,dP=-\left(G_t+\frac12\sigma^2P^2G_{PP}\right)dt.$$

Its value is $-G+PG_P$, so equating its gain to interest on that value gives the same PDE. Do not differentiate $-G+PG_P$ as though rebalancing shares were free: the cash financing must be included.

Terminal and spatial conditions

The contract supplies the terminal payoff:

$$G(T,P)=(P-K)^+\quad\text{or}\quad G(T,P)=(K-P)^+.$$

For a European call, $G(t,0)=0$ and $G(t,P)\sim P-Ke^{-r(T-t)}$ as $P\to\infty$. A numerical PDE solver also needs appropriate spatial boundary treatment. American exercise and path-dependent claims require additional structure; the displayed terminal-value problem describes European claims depending on the terminal asset price.

3. Why a heat equation appears

The heat equation spreads initial data through a Gaussian kernel. For $\lambda>0$,

$$u_\tau=\lambda u_{xx},\qquad u(x,0)=h(x),$$

has the convolution solution, for suitable initial data,

$$u(x,\tau)=\int_{-\infty}^{\infty}\frac{e^{-(y-x)^2/(4\lambda\tau)}}{\sqrt{4\pi\lambda\tau}}h(y)\,dy.$$

The fundamental solution for an initial unit point mass is the kernel itself. Its Gaussian variance is $2\lambda\tau$, and its integral over space is one.

For $h(x)=\mathbf1_{\{a<x<b\}}$,

$$u(x,\tau)=\Phi\!\left(\frac{b-x}{\sqrt{2\lambda\tau}}\right)-\Phi\!\left(\frac{a-x}{\sqrt{2\lambda\tau}}\right).$$

The denominator is $\sqrt{\tau}$ only when $\lambda=1/2$. For example, $\lambda=1/2$, $\tau=1$, $x=0$, and $(a,b)=(-1,1)$ give $u=2\Phi(1)-1\approx0.6827$.

An explicit Black–Scholes transformation

Set remaining time $\tau=T-t$, log price $x=\log P$, and

$$v(\tau,x)=G(T-\tau,e^x),\qquad c=r-\frac12\sigma^2.$$

The pricing PDE becomes

$$v_\tau=\frac12\sigma^2v_{xx}+cv_x-rv.$$

Discount removal and a moving coordinate eliminate the last two terms: write

$$v(\tau,x)=e^{-r\tau}u(\tau,x+c\tau).$$

Substitution gives $u_\tau=\tfrac12\sigma^2u_{xx}$ with $u(0,z)=g(e^z)$, where $g$ is the terminal payoff. Thus

$$G(t,P)=e^{-r\tau}E\left[g\left(Pe^{(r-\sigma^2/2)\tau+\sigma\sqrt{\tau}Z}\right)\right],\qquad Z\sim N(0,1).$$

This is the discounted risk-neutral expectation. The SDE describes paths, the PDE describes their valuation consequences, and the Gaussian kernel connects both views. See Robert Kohn’s NYU notes on heat equations and financial PDEs for the kernel and transformation approach.

4. Existence and uniqueness: when is the model well-defined?

A useful sufficient theorem assumes measurable deterministic coefficients, an $\mathcal F_0$-measurable initial value with $E|X_0|^2<\infty$, and, uniformly over a finite time interval,

$$|b(t,x)-b(t,y)|^2+|a(t,x)-a(t,y)|^2\le K|x-y|^2,$$ $$|b(t,x)|^2+|a(t,x)|^2\le K(1+|x|^2).$$

The first condition is global spatial Lipschitz continuity: nearby starting states cannot have arbitrarily different coefficients. The second is linear growth: coefficient magnitudes grow at most linearly in the state.

Under these conditions there is a continuous adapted strong solution with finite second-moment bounds on each finite horizon. Pathwise uniqueness means that solutions driven by the same Brownian motion and initial value satisfy

$$P(X_t=Y_t\text{ for every }t\in[0,T])=1.$$

“Strong” means the solution is constructed using the given driving noise and initial information. “Finite-horizon bounded” does not mean that its moments are bounded uniformly as $T\to\infty$. GBM can have rapidly growing moments while satisfying this theorem. These conditions are sufficient, not necessary: some useful models need different arguments.

5. Geometric Brownian motion by coefficient matching

Consider $dX_t=\mu X_tdt+\sigma X_tdB_t$, with $X_0=x_0>0$. Try $X_t=f(t,B_t)$. Itô’s formula requires

$$f_x=\sigma f,\qquad f_t+\frac12f_{xx}=\mu f.$$

The first equation gives $f(t,x)=A(t)e^{\sigma x}$. The second gives $A'(t)=(\mu-\sigma^2/2)A(t)$. Hence

$$\boxed{X_t=x_0e^{(\mu-\sigma^2/2)t+\sigma B_t}.}$$

In particular,

$$\log(X_t/x_0)\sim N((\mu-\sigma^2/2)t,\sigma^2t),$$ $$E[X_t]=x_0e^{\mu t},\qquad\operatorname{median}(X_t)=x_0e^{(\mu-\sigma^2/2)t}.$$

Coefficient matching is an ansatz, not a universal solution algorithm: many SDE solutions depend on the Brownian history and cannot be expressed as a function of $(t,B_t)$ alone.

Expected growth versus almost-sure decay

Brownian motion satisfies $B_t/t\to0$ almost surely. Therefore

$$\frac{1}{t}\log(X_t/x_0)\longrightarrow\mu-\frac12\sigma^2\quad\text{almost surely}.$$

If $0<\mu<\sigma^2/2$, then $X_t\to0$ almost surely even though $E[X_t]\to\infty$. Consequently, for every fixed $c>0$, $P(X_t<c)\to1$.

There is no contradiction: an increasingly rare and extreme upper tail dominates the expectation. Interchanging limit and expectation is invalid here; the family is not uniformly integrable.

For $x_0=100$, $\mu=0.02$, and $\sigma=0.30$, the log-growth rate is $-0.025$. At ten years, the expected value is about 122.14, while the median is about 77.88. These are distribution summaries, not forecasts of a particular realized path. If $\mu=\sigma^2/2$, the zero limiting log-growth rate alone does not imply convergence of $X_t$ to zero.

6. Ornstein–Uhlenbeck: solving mean reversion

Use $m$ for the long-run level to distinguish it from GBM’s drift:

$$dX_t=-\alpha(X_t-m)dt+\sigma dB_t,\qquad\alpha>0.$$

Multiplying $X_t-m$ by $e^{\alpha t}$ and applying the product rule yields

$$d\big(e^{\alpha t}(X_t-m)\big)=\sigma e^{\alpha t}dB_t.$$

Integrating gives

$$\boxed{X_t=m+(x_0-m)e^{-\alpha t}+\sigma\int_0^t e^{-\alpha(t-s)}dB_s.}$$

The deterministic kernel makes the integral Gaussian. Its mean is zero and its variance follows from Itô’s isometry:

$$E[X_t]=m+(x_0-m)e^{-\alpha t},\qquad\operatorname{Var}(X_t)=\frac{\sigma^2}{2\alpha}(1-e^{-2\alpha t}).$$

Stationarity and the Vasicek interpretation

From a fixed initial value, the distribution approaches $N(m,\sigma^2/(2\alpha))$. Paths continue fluctuating; they do not settle at $m$. Starting with this invariant distribution, independently of future Brownian increments, gives a stationary OU process with covariance

$$\operatorname{Cov}(X_t,X_{t+h})=\frac{\sigma^2}{2\alpha}e^{-\alpha|h|}.$$

The Vasicek short-rate model uses this form for interest rates. The fitted long-run level is a model parameter, not automatically a central bank’s announced target. Physical-measure parameters describe historical dynamics; risk-neutral parameters determine bond pricing. Gaussian rates can become negative.

For $x_0=0.05$, $m=0.03$, $\alpha=0.5$, and $\sigma=0.01$, the two-year mean is 3.7358% and standard deviation is 0.9299 percentage points. The expected deviation halves after $\log(2)/\alpha\approx1.3863$ years.

7. Multidimensional systems and correlated shocks

For $d$ state variables and $k$ independent Brownian drivers,

$$d\mathbf X_t=\mathbf b(t,\mathbf X_t)dt+A(t,\mathbf X_t)d\mathbf W_t,$$

where $A$ is a $d\times k$ matrix. The dimensions need not be equal. Instantaneous covariance is

$$d[X_i,X_j]_t=(AA^T)_{ij}dt.$$

Independent drivers can therefore produce correlated states through shared matrix loadings. For example,

$$dX_1=b_1dt+\sigma_1dW_1,\qquad dX_2=b_2dt+\sigma_2\left(\rho dW_1+\sqrt{1-\rho^2}dW_2\right)$$

has instantaneous covariance $\rho\sigma_1\sigma_2dt$. This construction supports cross-asset and multi-factor interest-rate models.

A position–velocity system illustrates unequal state and noise dimensions:

$$dQ_t=V_tdt,\qquad dV_t=-\alpha V_tdt+\sigma dB_t.$$

Two state variables share one noise driver. Position depends on the integrated velocity history.

For a smooth function, the generator is

$$\mathcal L f=\mathbf b\cdot\nabla f+\frac12\operatorname{tr}(AA^T\nabla^2f).$$

Under suitable regularity and integrability, a discounted terminal expectation solves the backward PDE $f_t+\mathcal L f-rf=0$. In pricing, use the risk-neutral generator; a physical-measure generator generally gives a different expectation.

8. Exact transitions before numerical approximations

For GBM, simulate exactly over a step $h$ using

$$X_{t+h}=X_t\exp((\mu-\sigma^2/2)h+\sigma\sqrt h Z).$$

For OU, use

$$X_{t+h}=m+(X_t-m)e^{-\alpha h}+\sigma\sqrt{\frac{1-e^{-2\alpha h}}{2\alpha}}Z,$$

with fresh independent standard normals across successive steps. These transitions avoid discretization error at sampled times. They do not reveal every intervening barrier crossing.

For a general SDE, Euler–Maruyama uses $X_{n+1}=X_n+b(t_n,X_n)h+a(t_n,X_n)\sqrt h Z_n$. Check convergence and stability: a coarse Euler step can make a positive GBM approximation negative or distort OU mean reversion.

Practice questions

1. Why does $\mu$ disappear from the Black–Scholes PDE?
2. What variance corresponds to $u_t=3u_{xx}$ after two time units?
3. What extra scaling is needed in the step-payoff CDF formula for general $\lambda$?
4. For GBM with $\mu=0.03$ and $\sigma=0.40$, what is the almost-sure long-run log-growth rate?
5. If OU has $\alpha=2$ and $\sigma=3$, what is its invariant variance?
6. Does convergence to an OU stationary distribution imply paths converge to a constant?
7. If $A$ is $5\times2$, how many states and independent Brownian drivers are present?
8. Why is differentiating $-G+PG_P$ alone insufficient to justify a dynamic hedge?

Answers

1. The replicating stock position cancels the asset drift on both sides of the self-financing equation, together with matching its Brownian exposure.
2. $2\lambda t=12$.
3. The denominator is $\sqrt{2\lambda t}$.
4. $0.03-0.40^2/2=-0.05$ per year; paths tend toward zero almost surely in this model.
5. $\sigma^2/(2\alpha)=9/4=2.25$.
6. No. The limiting distribution retains positive variance when $\sigma>0$.
7. Five states and two independent drivers; covariance has rank at most two at each instant.
8. Changing the share holding requires cash financing. The self-financing gains equation accounts for that rebalancing.

Further study

The MIT lecture on stochastic differential equations provides a companion treatment of hedging, diffusion equations, and OU dynamics. For the underlying calculus and SDE framework, consult NYU’s stochastic calculus course resources. Continue with Volatility Modeling to examine how constant-volatility diffusion assumptions are extended.

Feynman study lab

Explain it simply

An SDE specifies both the predictable part of motion and the random shocks.

Concrete example. For dX = 2 dt + 3 dB, a step of length 0.01 has drift increment 0.02 and shock standard deviation 0.3.

Teach it back without looking

Close the lesson and explain the idea aloud to a beginner in two minutes. Define each technical term you use, walk through the example, and state one assumption or limitation. Then answer:

Why is the shock scale 0.3 rather than 0.03?

Check your explanation after you try

Brownian increments have standard deviation equal to the square root of the time step. Multiply 3 by the square root of 0.01.

Find the gap and try again

If you hesitated, relied on a formula you could not explain, or missed the self-check, write down the exact gap. Revisit the relevant section of this lesson and work through its example. Close the notes, explain it again in simpler language, then change one number or assumption and predict what changes. Recalculate or check the lesson to test your prediction. Try once more tomorrow from memory.