M1.4 去噪就是取 Posterior Mean（Tweedie）

本篇重用3 篇

讀數是 72.3

全班量體重。用的體重計雜訊已知：每次讀數是真實體重加上一個標準差 1 kg 的 Gaussian 抖動，沒有系統偏差。全班的體重分佈也已知——老師去年量過，大致是平均 62 kg、標準差 8 kg 的鐘形。

某位同學站上去，讀數 72.3 kg。

你對他真實體重的最佳猜測是幾公斤？

雜訊對稱，只表示給定真實值時，讀數的平均沒有偏；看到讀數後反推真實值，還要看哪些來源本來就比較多。上一篇已經給了 posterior weights，平方誤差下的答案是它們加權的平均。

能不能不逐一算出所有 posterior weights，而直接從讀數分布算出這個平均？ 接下來的 Gaussian 微分，會把這個問題轉成密度的局部斜率。

課堂提問Q1

把情境翻成數學：隨機變數是什麼？「最佳猜測」的標準是什麼？「全班分佈已知」在這個問題裡扮演什麼角色？

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

最直覺的答案是「猜 72.3，因為雜訊平均是零」。這句話對了一半：它正確描述了 $Y$ 給定 $X$ 的分佈，卻回答錯了問題。

兩個隨機變數。 這位同學是從全班隨機挑出來的一個人，他的真實體重 $X$ 因此是隨機的，分佈 $p_X$ 就是全班的分佈（這是 prior）。讀數 $Y=X+\sigma\varepsilon$ ， $\varepsilon\sim\mathcal N(0,1)$ 、 $\sigma=1$ ——這是 likelihood：已知 $X=x$ ，讀數落在 $y$ 的密度是 $\varphi_\sigma(y-x)$ ，一個以 $x$ 為中心的鐘形。

最佳猜測的標準。 用平方誤差計分，M1.2 說最佳猜測是 $\mathbb E[X\mid Y=y]$ ：已知讀數是 $y$ ，真實體重的桶內平均。注意方向：直覺說的「雜訊平均是零」是 $\mathbb E[Y\mid X]=X$ ——已知真實值，讀數的平均；我們要的是反過來， $\mathbb E[X\mid Y]$ 。這兩個條件方向相反，正是 M1.3 那個「有病的人裡多少陽性」與「陽性的人裡多少有病」的差別。

全班分佈的角色。 反過來的條件需要 prior。「讀到 72.3」的人，有的是 71 kg 的人加了 $+1.3$ 的雜訊，有的是 74 kg 的人加了 $-1.7$ ——而全班裡 71 kg 的人比 74 kg 的人多（62 是中心，離得越遠越少），所以「讀到 72.3」的人裡偏輕的佔多數。桶內平均因此略低於 72.3。問題變成：這個修正量到底是多少、有沒有一條不必逐一算 posterior 的公式？

站在密度的斜坡上

想像許多人被隨機推向左右，現在站在讀數 72.3 的位置往回看。來自偏輕與偏重的人，不會只因雜訊對稱就各占一半：還要乘上兩邊原本的人數。這就是 posterior 重新加權的意思。

雜訊把每個人的真實體重往左右隨機推 1 kg 左右。站在 72.3 這個位置往回看：從左邊被推過來的人（本來更輕）比從右邊被推過來的人（本來更重）多，因為左邊本來就人多。所以「讀到 72.3」這一群人的真實體重，重心在 72.3 的左邊——往山頂方向。

Tweedie 會把這個平均回頭方向，寫成觀測密度 $p_Y$ 的 score，不是直接使用乾淨密度 $p_X$ 的 slope。兩個密度不要混在一起：noise 已經把原本的人群分布攤開了，下面的公式用的是攤開後的讀數。

讀數沿觀測密度的 score 往高密度方向修正。

圖：Tweedie 修正的方向由 score 決定，長度再乘上噪聲 variance。

三個工具接起來：Tweedie 公式

令 $Y=X+\sigma\varepsilon$ ，其中 $\sigma>0$ 、 $\varepsilon\sim\mathcal N(0,1)$ 與 $X$ 獨立，並假設相關條件平均存在、可以交換微分與積分。把讀數的密度寫出來：

p_Y(y)=\int p_X(x)\,\varphi_\sigma(y-x)\,dx,\qquad \varphi_\sigma(u)=\frac{1}{\sqrt{2\pi}\sigma}e^{-u^2/(2\sigma^2)} .

現在對 $y$ 微分，三步。

第一步（微分穿過積分，M0.3）。 積分變數是 $x$ 、微分變數是 $y$ ，Leibniz 法則讓 $\frac{d}{dy}$ 直接進到被積函數：

p_Y'(y)=\int p_X(x)\,\partial_y\varphi_\sigma(y-x)\,dx .

第二步（Gaussian 的 $\nabla\log$ ，M0.0）。 $\partial_y\varphi_\sigma(y-x)=-\dfrac{y-x}{\sigma^2}\,\varphi_\sigma(y-x)$ 。代入並把 $y$ 提出：

\sigma^2\,p_Y'(y)=\int(x-y)\,p_X(x)\varphi_\sigma(y-x)\,dx=\int x\,p_X(x)\varphi_\sigma(y-x)\,dx\;-\;y\,p_Y(y).

第三步（Bayes，M1.3）。 posterior 是 $p(x\mid y)=\dfrac{p_X(x)\varphi_\sigma(y-x)}{p_Y(y)}$ 。兩邊除以 $p_Y(y)$ ，第一項就變成用 posterior 加權的平均：

\sigma^2\,\frac{p_Y'(y)}{p_Y(y)}=\int x\,p(x\mid y)\,dx-y=\mathbb E[X\mid Y=y]-y .

整理，就是 Tweedie 公式：

\boxed{\ \mathbb E[X\mid Y=y]=y+\sigma^2\,\frac{d}{dy}\log p_Y(y)\ }

多維時一字不改， $\frac{d}{dy}$ 換成 gradient： $\mathbb E[X\mid Y=y]=y+\sigma^2\nabla_y\log p_Y(y)$ 。

讀法：最佳猜測＝讀數＋（雜訊的 variance）×（讀數密度的相對斜率）。 「往密度高處修正」是 $\nabla\log p_Y$ 的方向；「坡越陡修越多」是它的大小；「雜訊越大修越多」是前面的 $\sigma^2$ 。三個直覺各對應公式的一個因子。

展開細節兩個關鍵：只需要 p_Y；以及 ∇log 不在乎正規化

公式需要的是讀數 density $p_Y$ ，不必先逐點求出 prior。這和 Robbins [1] 的 empirical Bayes 思路相連：由一批 noisy observations 估計 marginal，再借它做推論。Efron [2] 討論 Tweedie 公式及其 selection bias 應用；實際由資料估 score 時，還會有 density estimation 的誤差，不能只畫 histogram 就當成精確答案。

$\nabla\log$ 不在乎正規化。 M0.0 說過 $\nabla\log(cf)=\nabla\log f$ 。所以估 $p_Y$ 時不需要它的積分等於 1，只要形狀對就夠——這在實務上省掉了一個最難算的量。

72.14，還是 71.0

代數字。 全班 $X\sim\mathcal N(62,8^2)$ 、雜訊 $\sigma=1$ 。讀數 $Y=X+\varepsilon$ 是兩個獨立 Gaussian 的和，由 M1.0， $Y\sim\mathcal N(62,\,64+1)$ 。它的 $\frac{d}{dy}\log p_Y(y)=-\dfrac{y-62}{65}$ 。代 $y=72.3$ ：

\mathbb E[X\mid Y=72.3]=72.3-1\times\frac{10.3}{65}=72.3-0.16=72.14\ \text{kg}.

換一台雜訊 $\sigma=3$ 的體重計： $Y\sim\mathcal N(62,73)$ ，修正 $=-9\times\dfrac{10.3}{73}=-1.27$ ，最佳猜測 71.03 kg。

Tweedie 告訴我們，posterior mean 恰好等於 y 的條件是該點的 marginal score 為零；它不要求整個 prior 都平坦。例如對稱分布的中心也會零修正。在這個 Gaussian 例子中， $\sigma=1$ 時修正約 $-0.16$ kg， $\sigma=3$ 時約 $-1.03$ kg；因為 marginal score 也隨 noise 改變，不能只看 $\sigma^2$ 增加九倍就認定修正也增加九倍。

這個結論依賴什麼。 第一，雜訊是 Gaussian——第二步用了 Gaussian 的 $\nabla\log$ 是線性的；換別種雜訊，公式的形狀會變。第二， $\sigma$ 要對——修正乘的是 $\sigma^2$ ，所以 $\sigma$ 估成一半，修正只剩四分之一。第三，prior 要是這個人所屬的母體——如果站上去的是來訪的老師，全班的分佈就不是他的 prior，公式會把他往學生的體重拉。第四，這是平方誤差下的最佳猜測；換計分方式，答案不再是 posterior mean（M1.2）。

模擬全班兩百萬次

import numpy as np
rng = np.random.default_rng(4); n = 2_000_000
def run(noise, sigma=1.0, ys=(66.0, 72.3, 78.0)):
    X = rng.normal(62, 8, n)                                          # 全班體重（prior）
    if noise == "gauss": eps = rng.normal(0, 1, n)
    else: eps = np.where(rng.random(n) < 0.95, rng.normal(0, 1, n), rng.normal(0, 10, n))  # 5% 機率跳格
    Y = X + sigma*eps
    emp = [X[np.abs(Y - y) < 0.1].mean() for y in ys]                 # 標準答案：讀數≈y 的人的平均真實體重
    hist, edges = np.histogram(Y, bins=350, range=(30, 100), density=True)
    c = 0.5*(edges[1:] + edges[:-1]); dlogp = np.gradient(np.log(hist + 1e-12), c)
    tweedie = np.array(ys) + sigma**2*np.interp(ys, c, dlogp)         # Tweedie：只用讀數的 histogram
    print(f"{noise:6s} σ={sigma}: 分桶平均 {np.round(emp, 2)}  Tweedie {np.round(tweedie, 2)}")
run("gauss", 1.0); run("gauss", 3.0); run("glitch", 1.0)

第一行：讀數 72.3 的分桶平均約 72.14、Tweedie 也給 72.14——公式在沒看過 $p_X$ 的情況下，只靠讀數的 histogram 就算出了正確的修正；讀數 78 修得更多（ $-0.25$ ），66 修得少（ $-0.06$ ），「坡越陡修越多」兌現。第二行 $\sigma=3$ ：72.3 修到約 71.0。第三行是刻意違反「雜訊是 Gaussian」：體重計有 5% 的機率跳出 $\pm10$ kg 的離譜值。這時分桶平均（標準答案）與 Tweedie 分岔——讀數 78 的人有相當比例其實是跳格的普通同學，標準答案把他們大幅往回拉，而 Tweedie 用 $\sigma^2=1$ 只修一點點。公式沒有錯，是它的前提沒了。

檢查理解

想一想

一個線上遊戲用玩家最近 20 場的平均分數當「實力估計」。排行榜前十名的玩家，下一個月的平均分數幾乎都會下降。用本篇的工具，最好的解釋是：

想一想

Tweedie 公式的推導裡，哪一步必須用到「雜訊是 Gaussian」？

想一想

承起點的例子（全班 $\mathcal N(62,8^2)$ 、 $\sigma=1$ ）。若站上體重計的是一位讀數 72.3 的來訪老師，而你仍用全班分佈套公式，會發生什麼事？

參考文獻

Robbins, H. An Empirical Bayes Approach to Statistics. Proceedings of the Third Berkeley Symposium on Mathematical Statistics and Probability, Vol. 1, 1956.（公式的出處之一；Robbins 將其歸功於 M. C. K. Tweedie 的私人通訊。全文在 Project Euclid 的 Berkeley Symposium 專輯裡，這裡不附連結——那個站的網址結構常變。）
Efron, B. Tweedie’s Formula and Selection Bias. Journal of the American Statistical Association 2011.（Tweedie 公式、empirical Bayes 與 selection bias。）
Stein, C. M. Estimation of the Mean of a Multivariate Normal Distribution. Annals of Statistics 1981.（Stein’s lemma——Gaussian 下「微分穿過期望」的同一個代數——與多維的收縮估計。）