M3.3 山谷裡隨機走的人群：Langevin Dynamics 與 Stationary Distribution

本篇重用2 篇

• M0.0霧中下山：Gradient 與方向
• M3.2一千片葉子在亂流裡：Fokker–Planck 方程式

一座山谷，一千個迷路的人

一座碗狀的山谷，起霧了。一千個人被放在山谷的同一側、同一個地方。每個人每一秒鐘隨機挪動一步——但因為腳下有坡，每一步都稍微傾向往低處：坡越陡，往下偏的成分越大；在谷底附近坡近乎零，就純粹亂走。

一小時後，這一千個人分佈在哪？如果換一座整體寬兩倍的山谷呢？谷底平坦、邊坡陡峭的山谷呢？

只有往低處走，大家可能聚到谷底；一直隨機挪動，又讓人散開。這兩種移動能不能維持一個固定的 density，而不是全部停在同一點？ 我們先找一個「流進與流出恰好抵消」的 distribution，再問是否真的會走到它。

課堂提問Q1

翻成數學。「山谷」是什麼？「傾向往低處」與「亂走」各寫成什麼？我們想知道的量——「很久以後大家在哪」——是哪個物件？

先想一想，再展開看整理後的答案

課堂上會先把山谷寫成一個函數：位置 $x$ 的高度 $V(x)$，谷底是 $V$ 的最小值。「往低處最陡的方向」是 $-\nabla V(x)$（M0.0），坡越陡向量越長——這正是題目說的「坡越陡往下偏越多」。所以 drift 是 $f(x)=-\nabla V(x)$。

「亂走」是 M3.1 的噪聲項。把兩者放進 SDE：

$$dx_t=-\nabla V(x_t)\,dt+\sqrt2\,dW_t .$$

這叫 Langevin dynamics（過阻尼版本）。$\sqrt2$ 是慣例，讓後面的答案乾淨；一般寫 $\sqrt{2\tau}$，$\tau$ 是「亂走的強度」，物理上是溫度。

Stationary distribution（平穩分布）的定義是：從這個 distribution 出發，之後的 distribution 不變，即 $\partial_t p=0$。這不表示每個人停止移動，也不自動保證任意 $p_0$ 都會收斂到它；存在、唯一與吸引其他初始分布，是不同的問題。

拉回與攤開的拔河

M3.2 說分佈的演化是兩股力量：drift 搬、噪聲攤。這裡 drift 永遠指向谷底，所以「搬」永遠是往中間拉回；噪聲則永遠把人群往外攤。

一開始一千人堆在谷的一側：這裡坡陡，拉回很強，整團被快速拉向谷底。到了谷底附近，坡變緩、拉回變弱，噪聲開始佔上風、把人群攤開。攤到某個寬度，被攤到邊坡上的人又因為坡陡被拉回。拉回與攤開在某個寬度上達成平衡——這個寬度由山谷的形狀決定，與出發點無關（出發點的記憶早就被亂走洗掉了）。

把同一座山谷橫向拉寬兩倍：每一處的坡都變緩、拉回變弱，要攤得更開才會碰到足夠陡的坡，所以人群也散開兩倍。反過來，窄而陡的山谷聚得緊。谷底特別平、但邊坡特別陡的山谷，則會給出頂部平坦、尾巴很薄的分佈——形狀跟著山谷走。

圖：−∇V 往谷底拉、噪聲往外攤；兩者在 e⁻ⱽ 的形狀上達到平衡。

把 $e^{-V}$ 代進 Fokker–Planck

Langevin dynamics 的 $f=-\nabla V$、$g=\sqrt2$，所以 $g^2/2=1$，Fokker–Planck 方程式（M3.2）是

$$\partial_tp=\nabla\cdot(p\,\nabla V)+\Delta p=\nabla\cdot\big(p\,\nabla V+\nabla p\big).$$

Stationary 的意思是右邊為零。猜 $p_\infty=e^{-V}/Z$，$Z=\int e^{-V}dx$（假設有限），然後代入驗證（一維寫法，多維逐字相同）：

$$\partial_xp_\infty=-V'\,p_\infty \quad\Longrightarrow\quad p_\infty V'+\partial_xp_\infty=p_\infty V'-V'p_\infty=0 .$$

括號裡的量——它是 flux——逐點為零，所以它的 divergence 當然為零，$\partial_tp_\infty=0$。✓ $\square$

注意符號：真正的 probability flux 是 $J=-p\nabla V-\nabla p$，因為守恆式是 $\partial_t p=-\nabla\cdot J$。上面括號內的量是 $-J$；代入 $p\propto e^{-V}$ 時兩項逐點抵消，所以確實沒有淨 flux。

會不會走到它？ 這需要另加 non-explosion、ergodicity 等條件；收斂速率也需要額外分析。$Z<\infty$ 讓候選 density 可以正規化，卻不足以保證指數收斂。Roberts 與 Tweedie [2] 分析了哪些尾部條件能或不能給出指數收斂；本篇的二次山谷則是可以直接求解的例子。

山谷形狀怎麼改答案。 $V=x^2/2$：$p_\infty=\mathcal N(0,1)$，標準差 1。橫向拉寬兩倍 $V=x^2/8$：$\mathcal N(0,4)$，標準差 2——散開兩倍。$V=x^2$：標準差 $0.71$——窄谷聚得緊。$V=x^4/4$（谷底平、邊坡陡）：$p_\infty\propto e^{-x^4/4}$，頂部平、尾巴比 Gaussian 薄，標準差約 $0.82$；這個數字沒有直覺可以猜——只有公式能給。加上溫度 $\tau$：$dx=-\nabla V\,dt+\sqrt{2\tau}\,dW$ 的 stationary distribution 是 $e^{-V/\tau}/Z_\tau$，$\tau$ 越大攤得越開。

展開細節離散步長不為零時，答案會偏：V = x²/2 的精確計算

用 Euler–Maruyama 模擬時步長 h 有限：xₙ₊₁ = (1−h)xₙ + √(2h)·εₙ。這是一條線性遞迴，它自己的 stationary variance 可以精確算：設 Var = s，則 s = (1−h)²s + 2h，解得 s = 2h / (1 − (1−h)²) = 2h / (2h − h²) = 1 / (1 − h/2)。 h → 0 時 s → 1；h = 0.1 時 variance 約偏 5%，h = 0.5 時是 1.33，h = 1 時是 2（標準差是 √2 倍，不是兩倍）。上述 stationary variance 需要 0 < h < 2；h ≥ 2 時噪聲持續累積，沒有有限的 stationary variance。這個二次例子的 variance bias 是 O(h)，不代表所有非線性 drift 的 Euler 鏈都穩定或有相同的保證 [2]。

直覺對了，但「很久」有多久

回答起點問題。一小時後（若一小時算「很久」），一千人的密度 $\propto e^{-V}$：與出發點無關、由山谷形狀完全決定。寬兩倍的谷讓人散開兩倍——直覺對，而且現在有數字；平底陡坡的谷給出平頂薄尾的分佈——這一個直覺猜不出來。

但這個判斷有兩個假設，各對應一種現實中的失效：

• 時間要夠久。 「指數速率收斂」的速率取決於 $V$ 的形狀。單一碗狀的谷收斂很快；若山谷有兩個谷底、中間隔一道高 $\Delta V$ 的山脊，理論上 $p_\infty\propto e^{-V}$ 仍然對（兩個谷底各占一團、比例由 $e^{-V}$ 決定），但從一側出發的人要翻過山脊才能填滿另一側，翻越的平均等待時間約 $e^{\Delta V/\tau}$ 的量級（Arrhenius 定律 [3]）。$\Delta V/\tau=10$ 時是兩萬多倍的時間尺度——「一小時」可能連一個人都沒翻過去。這時模擬的直方圖只填了一個谷，而理論的 $e^{-V}$ 說兩個谷都有人：理論沒錯，是你等得不夠久。
• 步長要夠小。 電腦模擬用有限的 $h$，離散鏈的 stationary distribution 與 $e^{-V}$ 差 $O(h)$（展開框有精確的算法：$V=x^2/2$ 時 variance 是 $1/(1-h/2)$）。$h=0.5$ 偏 33%，$h\ge2$ 直接發散。

跑一次，等不夠久一次，步長太大一次

import numpy as np
rng = np.random.default_rng(0)
def langevin(gradV, x0, T, h, n=20000):
    x = np.full(n, x0, float)
    for _ in range(int(T/h)):
        x = x - h*gradV(x) + np.sqrt(2*h)*rng.standard_normal(n)
    return x
x = langevin(lambda x: x, x0=3.0, T=20.0, h=0.01)          # 碗狀谷 V = x²/2，從 x=3 出發
print(x.mean().round(2), x.var().round(2))                  # ≈ 0.00, 1.00：與起點無關
for h in [0.1, 0.5, 1.0]:                                   # 步長太大
    print(h, langevin(lambda x: x, 3.0, 20.0, h).var().round(2), 1/(1-h/2))
# 雙井 V = k(x²−1)²/4，k=32（山脊高 ΔV = k/4 = 8），從 x=−1 出發
xb = langevin(lambda x: 32*x*(x**2-1), -1.0, 20.0, 0.01)
print((xb > 0).mean())                                      # ≈ 0.05，理論說 0.5；改 T=2000 得 ≈ 0.50

碗狀谷：從 $x=3$ 出發，$T=20$ 後均值 0、variance 1，與 $\mathcal N(0,1)$ 一致，起點的記憶消失。步長：$h=0.1,0.5,1.0$ 的 variance 分別約 $1.05$、$1.33$、$2.0$，與 $1/(1-h/2)$ 吻合——理論不只預告會偏，還算出偏多少。雙井（山脊高 8）：$T=20$ 時右井只有約 5% 的人，理論卻說各占一半；把 $T$ 拉長一百倍到 2000，右井才填到 50%。

檢查理解

想一想

一群氣體分子在重力場裡上下亂撞，高度 $z$ 的位能是 $V(z)=mgz$（$z\ge0$，地面是牆）。很久以後分子密度隨高度的分佈是： 均勻分佈，因為碰撞把一切攤平。 全部落在地面，因為重力一直往下拉。 $p(z)\propto e^{-mgz/\tau}$，指數遞減——這就是 Langevin 的 stationary distribution，也是大氣壓力隨高度指數下降的來源。

看解答

想一想

假設 continuous Langevin 確實收斂到 $e^{-V}/Z$，numerical approximation 也已檢查穩定，但從雙井的左側出發，有限時間內仍幾乎沒人到右側。合理的解釋是： 理論錯了，stationary distribution 取決於起點。 理論沒錯，但翻越山脊的等待時間約 $e^{\Delta V/\tau}$，你等的時間遠不夠；收斂的「指數速率」在這種地形上慢得驚人。 步長太大造成的偏差。

看解答

想一想

對二次山谷 $V(x)=x^2/2$，拿掉噪聲後只剩 $\frac{dx_t}{dt}=-x_t$。起點有有限二階矩，則 $t\to\infty$ 時： $x_t=e^{-t}x_0$ 趨近零，distribution 趨近 $\delta_0$。 還是 $e^{-V}$。 隨機分佈在整個山谷。

看解答

參考文獻

1. Pavliotis, G. A. Stochastic Processes and Applications: Diffusion Processes, the Fokker–Planck and Langevin Equations. Springer, 2014.（第 4.5 節與第 6 章：Langevin dynamics 的 stationary distribution、收斂到平衡的速率。）
2. Roberts, G. O., Tweedie, R. L. Exponential Convergence of Langevin Distributions and Their Discrete Approximations. Bernoulli 2(4), 1996.（連續 Langevin 的指數收斂；離散 Euler 版本的偏差與可能發散——展開框的現象。）
3. Gardiner, C. Stochastic Methods: A Handbook for the Natural and Social Sciences. Springer, 4th ed., 2009.（第 5 章：一維 Fokker–Planck 的 stationary solution；第 5.2.7 節：跨越位能障礙的 Kramers／Arrhenius 逃逸時間。）
4. Särkkä, S., Solin, A. Applied Stochastic Differential Equations. Cambridge University Press, 2019.（第 5.5 節：stationary Fokker–Planck 與線性 SDE 的 closed form。）