:從貝葉斯推斷到不確定性量化)
不確定性這個東西做數(shù)據(jù)建模的人天天掛在嘴邊但真正把它當回事、拿進模型里去算的其實少得可憐。大多數(shù)業(yè)務(wù)分析還是停留在“點估計”的層面——均值多少、預(yù)測值多少、系數(shù)多大頂多再給個置信區(qū)間應(yīng)付一下。問題是現(xiàn)實世界從來不給你一個確定答案。渠道投放效果會有波動用戶對價格的反應(yīng)有隨機性供應(yīng)鏈缺貨風(fēng)險更是層層疊加。你要是只用單一數(shù)值去描述這些現(xiàn)象等于把所有風(fēng)險都壓在一個毫無彈性的點上一旦真實值偏離整個決策鏈跟著遭殃。概率編程解決的就是這件事。它讓你把不確定性當成模型的一等公民不是事后貼上去的誤差條而是從數(shù)據(jù)生成機制里就帶出來的東西。用概率編程做不確定性推理本質(zhì)上是在回答一個更誠實的問題基于當前掌握的數(shù)據(jù)和先驗知識某個未知量到底可能落在哪些值上各自的可能性有多大。這個思路在營銷歸因、推薦系統(tǒng)、風(fēng)控建模、實驗分析這些場景里都適用而且一旦跑通你會發(fā)現(xiàn)它能帶來的信息量遠不止一個預(yù)測數(shù)字那么簡單。這篇文章寫給正在用傳統(tǒng)統(tǒng)計或機器學(xué)習(xí)方法做建模、想嘗試概率編程但不知道從哪下手的人。我會從概念拆解講到工具選型再用一個完整的實操案例帶你走一遍從數(shù)據(jù)到后驗分布的建模流程最后把我踩過的坑和排查經(jīng)驗一并整理出來。適合作為你入坑概率編程的第一篇參考。1. 為什么常規(guī)建模方法在不確定性面前不夠用1.1 你平時用的模型其實是在“假裝確定”先做一個思想實驗。你在做某產(chǎn)品的付費轉(zhuǎn)化率預(yù)估用歷史數(shù)據(jù)訓(xùn)練出一個模型算出轉(zhuǎn)化率是 2.3%。這個數(shù)字看起來精確但它背后省略了太多信息。當你說轉(zhuǎn)化率是 2.3% 的時候你其實是在說“我手頭所有數(shù)據(jù)指向的最佳猜測是 2.3%?!笨蓡栴}是這個猜測的置信度有多高如果樣本量只有 500 個用戶這個 2.3% 的可信區(qū)間可能橫跨 1.5% 到 3.5%如果樣本量是 50 萬區(qū)間可能收窄到 2.25% 到 2.35%。同樣的點估計背后代表的不確定性完全不同。傳統(tǒng)頻率派方法不是不能描述這種差異它用置信區(qū)間、標準誤來刻畫。但這里有幾個脫節(jié)的地方。第一置信區(qū)間的解釋方式極其繞它的嚴格含義是“重復(fù)抽樣下區(qū)間包含真實值的頻率”而不是“真實值落在這個區(qū)間的概率是 95%”——后者才是業(yè)務(wù)決策者真正想問的問題。第二當你引入先驗信息的時候頻率派沒有一套自然的框架來做這件事。第三復(fù)雜的層級結(jié)構(gòu)、缺失數(shù)據(jù)、測量誤差這些實際建模中常見的問題用傳統(tǒng)方法處理起來要么得做大量近似假設(shè)要么計算復(fù)雜到難以下手。概率編程換了一個問法。它直接問你你對某個未知量先驗上有什么判斷觀察到的數(shù)據(jù)長什么樣然后通過貝葉斯公式把這兩部分融合起來輸出一個完整的后驗分布。這個分布就是你要的答案全貌——不光告訴你最大概率的值在哪還告訴你其他值有多大可能。1.2 偶然不確定性和認知不確定性兩個必須分清的概念進入概率編程之前先要把兩種性質(zhì)完全不同的不確定性區(qū)分清楚因為它們對模型設(shè)計的影響是決定性的。偶然不確定性Aleatoric Uncertainty這類不確定性來自系統(tǒng)本身的隨機性是數(shù)據(jù)生成機制自帶的。比如拋硬幣的結(jié)果、用戶訪問時長的波動、傳感器讀數(shù)的噪聲。就算你收集一億個樣本這類不確定性不會消失它只是被刻畫得更精確而已。在概率編程里它通常體現(xiàn)在似然函數(shù)的噪聲項上。比如你假設(shè)觀測值服從正態(tài)分布那個方差參數(shù)就是在描述偶然不確定性的大小。認知不確定性Epistemic Uncertainty這類不確定性來自你認知的局限是“不知道”造成的。比如你只有 100 個樣本對模型參數(shù)的真實值心里沒底再比如你完全不確定某個渠道的投放效果究竟是正還是負。認知不確定性的特點是可以隨著數(shù)據(jù)積累而減小——樣本越多你對參數(shù)的認識越清晰后驗分布越集中。這兩種不確定性在業(yè)務(wù)上有著完全不同的含義。拿庫存預(yù)測來說偶然不確定性對應(yīng)的是需求本身的隨機波動你再怎么優(yōu)化模型也改變不了消費者行為的隨機性能做的是在安全庫存里給它留出緩沖認知不確定性對應(yīng)的是你對“需求均值到底是多少”這個參數(shù)的認識不清這個可以通過增加數(shù)據(jù)量、做更細致的分層來縮小。如果你在建模時把這兩種不確定性混為一談模型給出的風(fēng)險評估就是失真的。概率編程天然適合區(qū)分這兩者。通過為模型參數(shù)設(shè)置先驗分布并做后驗推斷你可以直觀地看到哪些不確定性在數(shù)據(jù)進來后大幅縮小了——那是認知不確定性在消退哪些沒有縮小——那些就是系統(tǒng)本身的偶然不確定性需要你接受并適應(yīng)它。這個視角上的轉(zhuǎn)變對實際決策的影響非常大。2. 概率編程里的核心概念從先驗到后驗的完整鏈路2.1 貝葉斯公式概率編程的地基概率編程的本質(zhì)是貝葉斯推斷在計算層面的工程化實現(xiàn)。貝葉斯公式表達起來非常簡單$$P(\theta | D) \frac{P(D | \theta) P(\theta)}{P(D)}$$翻譯成大白話$P(\theta)$ 是先驗分布在看到數(shù)據(jù)之前你對參數(shù) $\theta$ 的認知。它可以是信息性的比如“根據(jù)過往經(jīng)驗?zāi)城赖霓D(zhuǎn)化率大概率在 1%-3% 之間”也可以是弱信息性的比如“轉(zhuǎn)化率不會小于 0 也不會大于 100%”就行。$P(D | \theta)$ 是似然函數(shù)在給定參數(shù) $\theta$ 的情況下觀測到當前數(shù)據(jù)集 $D$ 的概率。它是模型與數(shù)據(jù)的連接點描述了數(shù)據(jù)生成機制。$P(\theta | D)$ 是后驗分布在看到數(shù)據(jù)之后你更新過的、對 $\theta$ 的最終認知。這是你做推斷、做預(yù)測、做決策的依據(jù)。理解貝葉斯公式的關(guān)鍵在于“更新”這個詞。先驗是你最初的判斷數(shù)據(jù)像秤砣一樣改變著判斷的天平后驗是了這個過程之后的結(jié)果。數(shù)據(jù)量越大、數(shù)據(jù)越“硬”后驗就越被數(shù)據(jù)主導(dǎo)先驗的影響越弱數(shù)據(jù)量小的時候后驗則主要是先驗在支撐——這也是貝葉斯方法在小樣本場景下依然能給出合理推斷的原因。概率編程做的事情用一句話概括就是你定義先驗和似然它幫你算后驗。這里的核心難點是$P(D)$ 這個歸一化常數(shù)在絕大多數(shù)現(xiàn)實模型里是沒有解析解的你需要用數(shù)值方法來逼近后驗分布。這就是 MCMC 和變分推斷登場的時刻。2.2 后驗推斷的兩條路MCMC 與變分推斷**MCMC馬爾可夫鏈蒙特卡洛**是目前概率編程中最主流的后驗推斷方法。它的思路很巧妙既然無法直接算出后驗的解析形式那就構(gòu)造一種采樣策略讓采樣出來的樣本服從的分布逐漸逼近后驗分布。具體實現(xiàn)的算法很多比較常用的是 No-U-Turn SamplerNUTS它是 Hamilton Monte CarloHMC的一種自適應(yīng)變體。NUTS 能在參數(shù)空間里利用梯度信息高效地探索后驗分布在參數(shù)維度不高、后驗形態(tài)不算太扭曲的情況下收斂速度相當理想。使用 MCMC 得到的是一堆從后驗分布中抽取的樣本。通過對這些樣本做統(tǒng)計你可以得到后驗均值、分位數(shù)、可信區(qū)間等一切你關(guān)心的量。樣本量越大對后驗的刻畫就越精細。它的優(yōu)點是準確——沒有系統(tǒng)性偏差只要鏈長足夠、診斷通過樣本分布就趨近于真實后驗。缺點是需要大量計算尤其在高維參數(shù)空間或大數(shù)據(jù)集場景下采樣成本可能高到不可接受。**變分推斷Variational Inference**是另一個方向。它不是采樣而是把后驗推斷轉(zhuǎn)化成一個最優(yōu)化問題設(shè)定一個參數(shù)化的分布族 $q(\theta; \phi)$然后調(diào)整變分參數(shù) $\phi$讓 $q$ 盡可能接近真實后驗 $P(\theta | D)$通常用 KL 散度來衡量兩個分布的距離。優(yōu)化完成后$q$ 就作為后驗的近似輸出。變分推斷的優(yōu)勢是速度快能處理大規(guī)模數(shù)據(jù)和復(fù)雜模型。缺點是有偏——$q$ 的分布族選擇如果表達能力不夠近似結(jié)果就和真實后驗有偏差?,F(xiàn)實使用中一般是這樣取舍的先期探索模型結(jié)構(gòu)、需要快速迭代的時候用變分推斷最終關(guān)鍵結(jié)論的確定、當模型維度不算太高的時候用 MCMC。我在實際項目里經(jīng)常兩個一起跑變分推斷給出的結(jié)果當“預(yù)報員”MCMC 的結(jié)果當“裁判員”。2.3 概率編程語言到底“編程”了什么很多人第一次看到 PyMC 或 Stan 的代碼時會有個困惑這不就是寫個模型公式然后調(diào) sample() 嗎這算什么編程這個問題的答案在于概率編程語言提供的是一種生成式模型的表達框架。你寫的模型代碼實際上是在描述“數(shù)據(jù)是如何被生成出來的”——先定義參數(shù)的先驗分布再定義中間變量之間的確定性關(guān)系最后定義觀測變量的似然分布??蚣茇撠?zé)在后臺把這個描述轉(zhuǎn)換成對應(yīng)的計算圖然后自動推導(dǎo)采樣所需的梯度信息、自動完成推斷過程。這種“聲明式建?!钡姆绞綆韮蓚€巨大的好處。第一模型和推斷解耦你可以把精力全部放在模型結(jié)構(gòu)的合理性上而不需要自己實現(xiàn)采樣算法、推導(dǎo)梯度公式。第二模型靈活度極高你可以輕松地構(gòu)建多層級的回歸模型、帶有測量誤差的模型、混合分布模型等復(fù)雜結(jié)構(gòu)這些用傳統(tǒng)統(tǒng)計軟件做會非常痛苦但在概率編程框架里只是多寫幾行代碼的事。當然這并不意味著你可以完全不懂推斷原理。實際上對 MCMC 的收斂性、采樣器行為、先驗敏感性有基本了解是排查模型問題的前提——這個問題我在后面的章節(jié)會詳細展開。3. 工具選型主流概率編程框架怎么選3.1 PyMC、Stan、NumPyro 三足鼎立概率編程的生態(tài)經(jīng)過這么多年的發(fā)展基本形成了三足鼎立的局面老牌的 PyMC前身是 PyMC2/PyMC3、工業(yè)級口碑的 Stan尤其它的 R 接口、以及 JAX 生態(tài)下的新銳 NumPyro。選型沒有絕對的正確答案關(guān)鍵看你的技術(shù)棧和需求。PyMCPython 生態(tài)中最成熟、社區(qū)最活躍的選擇。API 設(shè)計對新手非常友好內(nèi)置的pm.Model()上下文管理器讓模型定義像寫公式一樣自然。它封裝了多種采樣器和變分推斷方法pm.sample()一把梭就能完成推斷。配合arviz庫做后驗診斷和可視化體驗非常順滑。如果你的主力語言是 Python想快速上手、快速驗證想法PyMC 基本是首選。Stan以 C 為基礎(chǔ)實現(xiàn)通過cmdstanpy或 R 的rstan接口調(diào)用。Stan 的建模語言是獨立的 Stan Language需要在.stan文件中寫模型然后編譯運行。它的 HMC/NUTS 實現(xiàn)非常扎實在某些復(fù)雜模型和后驗幾何相對崎嶇的情況下采樣效率往往優(yōu)于其他框架。代價是建模和調(diào)試的迭代周期比 Python 原生環(huán)境要長一些。適合需要極致采樣性能、模型復(fù)雜度高的場景。NumPyro基于 JAX 構(gòu)建用函數(shù)式編程風(fēng)格來描述模型。它的最大優(yōu)勢是支持 JIT 編譯和 GPU 加速在對性能極度敏感的大規(guī)模任務(wù)上有明顯優(yōu)勢。如果你是 JAX 生態(tài)的玩家或者模型大到采樣不得不靠 GPU 加速NumPyro 非常值得嘗試。學(xué)習(xí)曲線相對陡——函數(shù)式風(fēng)格、numpyro.sample的調(diào)用方式和plate與numpyro.handlers的交互邏輯都需要一些適應(yīng)時間。3.2 我的選型建議按場景對號入座這里給你一個比較實際的參考場景推薦工具理由Python 新手 / 快速原型PyMC上手快、社區(qū)資料多、可視化配套完善復(fù)雜模型 / 需要最高采樣質(zhì)量StanNUTS 實現(xiàn)成熟對崎嶇后驗的適應(yīng)性強大數(shù)據(jù)量 / GPU 加速 / JAX 生態(tài)NumPyroJIT 向量化性能天花板最高R 用戶做統(tǒng)計建模rstan與 R 的統(tǒng)計生態(tài)銜接自然生產(chǎn)環(huán)境部署NumPyro / PyMC 序列化Python 棧部署路徑最順還有兩個輕量級選手值得一提。sampyl和emcee都是輕量級 MCMC 庫體量小、依賴少適合教學(xué)和小型項目。但它們功能覆蓋面有限生產(chǎn)級建模我還是建議回到三大件上來。我的個人經(jīng)驗是不要貪多主攻一個熟悉另一個。比如你日常用 PyMC 建模遇到采樣性能瓶頸時能切到 Stan 去對比驗證或者你在 PyMC 里原型驗證沒問題之后用 NumPyro 在生產(chǎn)環(huán)境里重寫一遍上 GPU??蚣苤g模型表達的思路是相通的重點是把核心概念吃透。4. 完整實操用概率編程做一個轉(zhuǎn)化率不確定性分析4.1 業(yè)務(wù)場景與數(shù)據(jù)構(gòu)造我們來走一個完整的建模流程。場景設(shè)定某產(chǎn)品線做了一個新版本的功能產(chǎn)品經(jīng)理想知道新版上線后購買轉(zhuǎn)化率到底比舊版高多少。收集到的數(shù)據(jù)形式很簡單舊版觀察了 2000 個用戶其中 96 人完成購買新版觀察了 2140 個用戶其中 132 人完成購買。傳統(tǒng)做法是算轉(zhuǎn)換率、做個卡方檢驗或者兩比例 z 檢驗給個 p 值。但概率編程能給出的信息豐富得多你能直接得到轉(zhuǎn)化率差異的完整后驗分布并基于概率做出決策比如“新版轉(zhuǎn)化率比舊版高的概率是 93.7%”。為了演示我們用 Python 的 PyMC 框架。核心代碼如下import pymc as pm import arviz as az import numpy as np # 觀測數(shù)據(jù) n_old, x_old 2000, 96 n_new, x_new 2140, 132 with pm.Model() as conversion_model: # 先驗舊版本轉(zhuǎn)化率 ~ Beta(1, 1) 等價于 Uniform(0, 1) p_old pm.Beta(p_old, alpha1, beta1) p_new pm.Beta(p_new, alpha1, beta1) # 似然觀測到的購買人數(shù)服從二項分布 obs_old pm.Binomial(obs_old, nn_old, pp_old, observedx_old) obs_new pm.Binomial(obs_new, nn_new, pp_new, observedx_new) # 派生量轉(zhuǎn)化率差異 diff pm.Deterministic(diff, p_new - p_old) # 執(zhí)行 MCMC 采樣 trace pm.sample(draws2000, tune2000, chains4, random_seed42)這段代碼一共做了四件事定義轉(zhuǎn)化率的先驗Beta(1,1) 表示完全沒有任何先驗傾向、定義觀測數(shù)據(jù)的似然二項分布、定義關(guān)心的派生量轉(zhuǎn)化率差異、執(zhí)行 NUTS 采樣。跑完之后trace里存的就是所有參數(shù)的 MCMC 樣本。4.2 后驗分析從樣本到業(yè)務(wù)決策采樣完成后我們把后驗樣本展示出來az.summary(trace, var_names[p_old, p_new, diff])你會得到一張表格每一行一個變量包含均值、標準差、HDI最高密度區(qū)間、有效樣本量等指標。假設(shè)輸出結(jié)果是mean sd hdi_3% hdi_97% ess_bulk r_hat p_old 0.048 0.003 0.042 0.054 6234 1.00 p_new 0.062 0.004 0.055 0.069 6157 1.00 diff 0.013 0.005 0.004 0.022 6128 1.00這份結(jié)果信息的密度遠超一個 p 值。首先p_old和p_new的 HDI 區(qū)間都是兩側(cè)遠離 0 的表明數(shù)據(jù)充分支持這兩個轉(zhuǎn)化率確實在 4%-7% 這個量級。其次diff的 94% HDI 區(qū)間是 [0.004, 0.022]不包含 0——這說明新版轉(zhuǎn)化率顯著高于舊版等等先別急著說“顯著”在貝葉斯框架里我們換一種更直白的表述。更直觀的做法是直接計算概率prob_improvement (trace.posterior[diff] 0).mean().item() print(f新版轉(zhuǎn)化率高于舊版的概率: {prob_improvement:.1%})假設(shè)輸出為 98.3%那這個信息可以直接扔給業(yè)務(wù)方我們不需要糾結(jié)“顯著性水平 0.05”或“p0.03”這種解釋易產(chǎn)生歧義的表述而是可以說“根據(jù)數(shù)據(jù)表現(xiàn)新版轉(zhuǎn)化率高于舊版的概率為 98.3%”。這個概率表述業(yè)務(wù)方一聽就懂決策者可以直接基于它來評估是否全量上線。4.3 先驗敏感性分析給結(jié)論上雙保險有經(jīng)驗的建模者不會只看一組先驗跑出來的結(jié)果就下結(jié)論。貝葉斯方法的先驗選擇在一定程度上是主觀判斷所以要做敏感性分析換不同的先驗設(shè)置看結(jié)論會不會被顯著改變。我們把代碼改成用信息性強的先驗再跑一遍with pm.Model() as conversion_model_prior: # 更強的先驗根據(jù)歷史經(jīng)驗轉(zhuǎn)化率集中在 3%-7% 之間 p_old pm.Beta(p_old, alpha10, beta170) p_new pm.Beta(p_new, alpha10, beta170) obs_old pm.Binomial(obs_old, nn_old, pp_old, observedx_old) obs_new pm.Binomial(obs_new, nn_new, pp_new, observedx_new) diff pm.Deterministic(diff, p_new - p_old) trace_prior pm.sample(draws2000, tune2000, chains4, random_seed42)對比兩次運行的diff后驗均值如果結(jié)論基本一致說明你的結(jié)論受先驗影響不大、是數(shù)據(jù)本身在驅(qū)動如果差異明顯說明當前數(shù)據(jù)量還不足以克服先驗的影響。這時候就要慎重對結(jié)論的表述要加上“基于當前先驗判斷”的限定條件。這是貝葉斯建模的好習(xí)慣也是規(guī)避“先驗挑數(shù)據(jù)”這個爭議的有效方式。4.4 層級模型的擴展不滿足于單一體轉(zhuǎn)化率上述模型是極其簡化的——它假設(shè)所有用戶同質(zhì)忽略用戶來源渠道、訪問設(shè)備、時間區(qū)間等因素。一個更接近真實業(yè)務(wù)的建模方式是把這些變量吸納進模型。設(shè)想數(shù)據(jù)升級為用戶來源渠道有 A、B、C 三個每個渠道既有舊版也有新版訪問。你可以在模型里加入渠道作為分組變量構(gòu)建一個層級模型# 假設(shè) channels 是 0,1,2 編碼的渠道標簽 import pandas as pd # 這里簡化為構(gòu)造一個示例 DataFrame np.random.seed(42) df pd.DataFrame({ channel: np.repeat([0, 1, 2], [800, 700, 640]), version: 0 # 占位 }) # 實際使用時數(shù)據(jù)是每個 [channel, version] 組合的訪問量和購買量 with pm.Model() as hierarchical_model: # 渠道間共享的超先驗 mu pm.Normal(mu, mu0, sigma1) sigma pm.HalfCauchy(sigma, beta0.5) # 渠道效應(yīng)收縮到共同均值 z_channel pm.Normal(z_channel, mu0, sigma1, shape3) beta_channel pm.Deterministic(beta_channel, mu sigma * z_channel) # 版本效應(yīng)按渠道建模時也能分層處理 p_logit mu beta_channel * ... # 根據(jù)具體編碼設(shè)計 # 似然觀測銷售數(shù)據(jù)服從伯努利/二項分布 ...層級模型的好處在于部分池化樣本量少的渠道會參考整體分布的信息避免因樣本稀疏產(chǎn)生極端的估計樣本量大的渠道則保持自己的特性。這在轉(zhuǎn)化率、客單價等業(yè)務(wù)指標的多維度分析中非常實用也是傳統(tǒng)頻率派方法處理起來最棘手、而概率編程最舒展的場景之一。5. 建模實操中的常見問題與排查經(jīng)驗5.1 采樣不收斂最讓人抓狂的 MCMC 問題跑pm.sample()之后第一件事永遠是看r_hat和有效樣本量ess_bulk。r_hat大于 1.01 基本就意味著采樣鏈沒有收斂。常見原因和排查手段我整理成了一張速查表癥狀可能原因排查方向r_hat 恒大于 1.01鏈初始值差異太大burn-in 不夠增加 tune 長度提供更合理的初始值鏈呈“毛蟲”狀、游走緩慢后驗分布相關(guān)性高 / 參數(shù)化不佳嘗試標準化參數(shù)、非中心化參數(shù)化高維模型頻繁發(fā)散后驗幾何太崎嶇調(diào)整步長、換用pm.Metropolis對比、簡化模型結(jié)構(gòu)eff 樣本量極低500自相關(guān)性過高增加chains或draws檢查模型是否有冗余參數(shù)警告出現(xiàn)很多divergencesNUTS 步長自適應(yīng)遇到問題通常因先驗與數(shù)據(jù)沖突檢查先驗是否有不合理取值并調(diào)整或降低target_accept到 0.9有參數(shù)后驗分布“掛在邊界”先驗約束不夠或者數(shù)據(jù)對參數(shù)的識別力弱考慮更換先驗分布族例如在比例參數(shù)上使用 Beta 而不是 Uniform我遇到最多的還是參數(shù)化不當導(dǎo)致的無效采樣。一個典型的例子是方差參數(shù)用了HalfNormal(sigma10)這樣過寬的弱先驗。表面上看起來“不預(yù)設(shè)立場很客觀”但實際上過寬的先驗會在高參數(shù)空間分配大量概率質(zhì)量導(dǎo)致采樣器在極端值之間掙扎表現(xiàn)為大量發(fā)散。改用HalfCauchy或適度窄的先驗分布之后采樣立刻恢復(fù)順暢。這里想多說一句弱先驗不等于打死不設(shè)先驗。在你對參數(shù)的合理范圍多少有數(shù)的情況下給一個符合業(yè)務(wù)常識的約束或弱信息先驗往往既不會扭曲后驗結(jié)論又能顯著改善采樣性能。這也是貝葉斯工程師的經(jīng)驗之談。5.2 模型結(jié)果怪得離譜先檢查似然和數(shù)據(jù)的對應(yīng)關(guān)系有次我做一個客單價預(yù)估模型后驗均值明顯偏離業(yè)務(wù)直覺。排查了半天最后發(fā)現(xiàn)是“數(shù)據(jù)對齊”的問題訓(xùn)練數(shù)據(jù)里的一列應(yīng)是當日訂單量結(jié)果因為某次處理腳本的 bug混入了前一天的訂單數(shù)據(jù)。這類“模型以外”的錯誤在統(tǒng)計推斷里比模型本身的錯誤更隱蔽。排查順序建議為檢查數(shù)據(jù)描述性統(tǒng)計——均值、分位數(shù)、缺失率是否合理。畫出關(guān)鍵變量在模型下的生成模擬數(shù)據(jù)posterior_predictive對比真實數(shù)據(jù)的分布。如果模擬數(shù)據(jù)和觀測數(shù)據(jù)差別巨大說明模型的一個或多個核心假設(shè)有問題。逐一簡化模型先跑一個只含主要效應(yīng)、不含交互項和層級結(jié)構(gòu)的模型看看能否得到大概合理的后驗再逐步增加復(fù)雜度。對連續(xù)型變量檢查數(shù)據(jù)縮放是否均勻。數(shù)值范圍跨度過大的特征會讓 NUTS 的步長自適應(yīng)失靈建議對輸入數(shù)據(jù)做標準化同時可以在模型的確定性關(guān)系中使用相應(yīng)的偏差。曾有很長一段時間我在做這類分析時從不做第 2 步——posterior_predictive檢查——后來一個實際案例教育了我模型對某個分組變量的層級方差估得很小把分組間差異全都“解釋”成了偶然噪聲完全不符合業(yè)務(wù)預(yù)期。當時的模型層設(shè)計沒有包含組內(nèi)樣本量對先前方差收縮作用的權(quán)衡簡單來說就是層級結(jié)構(gòu)設(shè)置不當。如果沒有事后預(yù)測檢查這個問題會一直潛伏到部署階段才爆發(fā)。5.3 變分推斷的結(jié)果不靠譜別急著怪算法使用變分推斷時最容易掉進的坑是q的分布族表達能力不足導(dǎo)致求出的近似后驗和真實后驗差異巨大但收斂過程看起來完全正常。一種驗證方法是把變分推斷結(jié)果和 MCMC 結(jié)果做對比。在數(shù)據(jù)量不大、模型不太復(fù)雜的時候這個對比成本有限卻能一針見血地發(fā)現(xiàn)問題。如果兩者差異在可接受范圍內(nèi)就用變分推斷如果差異巨大檢查是不是原始模型的非光滑環(huán)節(jié)比如離散變量、硬截斷導(dǎo)致變分分布族嚴重失配。另外變分推斷對初始化也比較敏感多試幾組隨機種子如果結(jié)果的不確定性變動異常說明優(yōu)化過程本身可能有問題。從工程實踐角度看我的策略是變分推斷用在探索階段MCMC 用在發(fā)布階段。每次項目保質(zhì)期內(nèi)的最終結(jié)論都以 MCMC 驗證過的版本為準。5.4 業(yè)務(wù)落地時的不確定性傳播最后分享一個經(jīng)常被忽略但實際影響很大的問題——概率編程的輸出如何進入下游系統(tǒng)。一種常見做法是只把后驗均值傳到下游這等于又把不確定性丟掉了。較好的做法是把后驗分布的完整樣本傳遞過去。具體而言如果你在做庫存優(yōu)化不是傳一個“預(yù)測需求 1300 件”而是傳一組后驗樣本“需求可能分布在 [1180, 1420]”下游的庫存優(yōu)化代碼可以對每一個樣本跑一遍最優(yōu)庫存決策得到一整套決策的分布最終選擇覆蓋某個分位的庫存量。如果你在做 A/B 實驗的后續(xù)部署判斷不是簡單用“贏了的版本全量上”而是用后驗差異的分布來確定上線節(jié)奏和回滾閾值。在做預(yù)算分配時可以把每個渠道的收益參數(shù)后驗樣本輸入優(yōu)化器這樣得到的分配方案本身就是帶不確定性的決策者可以看到最壞情況和最好情況之間的差異。很多業(yè)務(wù)場景本身并不缺少優(yōu)化方法缺少的是給優(yōu)化方法提供不確定性輸入的意識。概率編程給出了這個輸入它是從純粹統(tǒng)計推斷到實際業(yè)務(wù)決策的橋梁。6. 實戰(zhàn)進階從簡單模型走向生產(chǎn)級應(yīng)用6.1 從單模型到貝葉斯工作流把概率編程真正用起來之后你會發(fā)現(xiàn)自己需要的不只是一個采樣函數(shù)而是一整套工作流?;谖业膶崙?zhàn)經(jīng)驗一個可復(fù)用的貝葉斯建模流程應(yīng)該是這樣的問題形式化把業(yè)務(wù)問題翻譯成一個可以計算概率的可量化問題。明確哪些是觀測變量、哪些是待推斷參數(shù)、哪些是派生業(yè)務(wù)指標。數(shù)據(jù)審查做描述性統(tǒng)計和可視化識別數(shù)據(jù)缺失、異常值和量綱問題。先驗設(shè)立基于業(yè)務(wù)常識設(shè)定合理的弱信息先驗。向業(yè)務(wù)方了解“過去的歷史水平”讓它成為先驗均值的重要參考。模型搭建與迭代從簡單模型開始驗證核心結(jié)構(gòu)正確后再加復(fù)雜度。每次增加結(jié)構(gòu)后都用后驗預(yù)測檢查來檢驗改進效果。推斷與診斷執(zhí)行 MCMC或變分推斷檢查r_hat、ess_bulk、發(fā)散情況。必要時做參數(shù)化調(diào)整。敏感性分析改變先驗設(shè)定、修改數(shù)據(jù)子集做重跑確認結(jié)論的穩(wěn)定性。產(chǎn)出與溝通把后驗分布翻譯成業(yè)務(wù)語言輸出“差異為 X 的概率是多少”、“落在哪個區(qū)間的把握有多大”這類可操作的結(jié)論。部署與監(jiān)控把模型固化為可重復(fù)執(zhí)行的推斷管道對新數(shù)據(jù)周期性地更新后驗并監(jiān)控分布漂移。這八步聽起來多但每一步的時間消耗并不大。熟練之后一個中等復(fù)雜度的模型從開始到出結(jié)論基本在一兩天內(nèi)就能完成。6.2 讓模型可復(fù)現(xiàn)的幾個工程習(xí)慣概率編程的模型天然具有隨機性因此可復(fù)現(xiàn)性是工程化繞不開的問題。沒有固定隨機種子的模型運行每次結(jié)果都有細微差異這在調(diào)試時會讓人分不清改動模型前的效果差異是真實的還是隨機波動造成的。我現(xiàn)在在項目里嚴格遵守以下三條規(guī)則固定隨機種子所有采樣、變分推斷調(diào)用都顯式傳入random_seed并使它在每次運行中保持一致。記錄模型版本和依賴環(huán)境把模型代碼、依賴庫版本號、數(shù)據(jù)文件哈希寫進自動生成的模型元數(shù)據(jù)文件。后期排查問題這一條能救命。用緩存機制保存推斷結(jié)果同一個模型、同一份數(shù)據(jù)跑過一次后把 trace 序列化為文件。后續(xù)所有人的分析都基于這份緩存結(jié)果保證溝通討論的前后一致。這些習(xí)慣沒有任何高深之處但它們直接決定了概率編程能否從小項目走向正式工作流。7. 一些關(guān)于概率編程的實話做概率編程這些年最大的感受是這項技術(shù)真正改變?nèi)说牡胤讲辉诓蓸铀惴ǖ男识谒季S方式的轉(zhuǎn)換。當你開始用分布去思考所有問題你會意識到以前那些點估計的結(jié)論有多少是被“平均”掩埋掉的。但也要說句公道話概率編程不是萬能的。它的計算成本遠高于傳統(tǒng)統(tǒng)計方法在超大樣本場景下需要投入相當?shù)臅r間去做工程優(yōu)化。如果業(yè)務(wù)場景的數(shù)據(jù)量足夠大、模型相對簡單、只需要一個穩(wěn)定輸出傳統(tǒng)方法可能更實惠。概率編程適合的是對不確定性本身有決策價值、數(shù)據(jù)規(guī)??煽?、模型復(fù)雜度中上的場景——用一句直白的話說當你發(fā)現(xiàn)“上線后預(yù)測總偏差大、又不知道偏差從哪來”的時候概率編程大概能給你答案。另外概率編程也不該被你當成一個“萬能答案機”。它給出的后驗分布依然建立在模型假設(shè)之上。假設(shè)錯了后驗再漂亮也不可靠。這也是為什么做概率編程的人一半時間在建模另一半時間在驗證假設(shè)、檢查模型與數(shù)據(jù)的匹配度。這些約束和代價不是說概率編程不好而是說你有必要想清楚我到底需要的是不確定性信息還是僅僅需要結(jié)論的確定性去完成眼前這一步。想明白了再上手這個工具給你的回報會遠高于你的預(yù)期。