
均值骗人PyMC 贝叶斯分位数回归抓需求上界【免费下载链接】pymcBayesian Modeling and Probabilistic Programming in Python项目地址: https://gitcode.com/GitHub_Trending/py/pymc大促前主管拍板备货按九成九不缺货来定。你算了日均 850 件就照它下单结果第三天断货——平均数压不住促销期抬起来的尾部。要的是一条九成概率不超它的线也就是条件分位数。用 PyMC 贝叶斯分位数回归就是直接把这条上界建模出来。读完你能写出 AsymmetricLaplace 模型、看懂收敛诊断再把后验区间翻译成能下单的备货量。一个帮你建立直觉的比喻先看个熟悉的东西体检报告上的某项指标医生不会只盯着人群平均而是看你的值落在第几分位——低于 P10 要警惕高于 P90 也要警惕。均值只告诉你人群中心在哪分位数告诉你整条分布长什么样、尾巴有多长。天气预报更直白。它说明天降雨 5 毫米对你没用你要的是会不会大到淋湿这种阈值判断而阈值判断靠的就是分位。回归也一样。传统线性回归给的是条件均值也就是 Y 的中心点。分位数回归给的是整条条件分布的轮廓你可以同时描出 P10、P50、P90 三条线把典型值和极端值一起框住。说白了非正态分布回归想做的事就是不逼数据服从正态而是承认分布有形状然后直接去描这个形状。你可以把它理解为均值是拍一张正面照分位数是绕着物体转一圈。30 秒跑通第一个模型下面这段能直接跑造一条带异方差噪声的线再用AsymmetricLaplace把 90% 分位当似然拟合出来。数据怎么造的一行带过重点看分布定义那几行。import numpy as np import pymc as pm rng np.random.default_rng(0) x rng.uniform(0, 10, 300) # 自变量比如促销力度 y 2.0 0.8 * x rng.normal(0, 0.5 * (x / 10 0.3), 300) # 模拟数据噪声随 x 变大 with pm.Model() as m: beta0 pm.Normal(beta0, 0, 10) # 截距先验 beta1 pm.Normal(beta1, 0, 10) # 斜率先验 sigma pm.HalfNormal(sigma, 2) # 尺度参数先验 mu beta0 beta1 * x # 线性预测条件分位数函数 pm.AsymmetricLaplace(y_obs, mumu, bsigma, q0.9, observedy) # q0.9 建模 90% 分位 idata pm.sample(1000, tune500, progressbarFalse) # MCMC 采样三个参数各管一件事mu就是你要预测的那条分位线本身业务上是九成不超它的需求量sigma是分布的离散程度越大说明同一个 x 下结果越散q决定你盯哪条线q0.9就是 90% 分位想看中位数就填 0.5。这里用 AsymmetricLaplace 当似然是因为它天然把 q 编码进损失里q 越大越重罚预测偏低正好匹配宁可多备也别缺货的心态。这个分布的实现就放在 continuous.py 里。PyMC 会把上面这些节点自动连成一张概率图模型随机变量是自由节点派生量挂在它们下游。怎么判断采样结果靠谱采样完别急着用先看三个信号。第一是 R-hat衡量各条链有没有混到同一处落在 1.01 以内算稳超过 1.1 就该重跑。第二是 ESS有效样本量它不是越大越好的装饰而是你这几千个样本里有多少个真正独立太低说明样本高度自相关后验区间不可信。第三是后验预测图让模型重新生成一批假数据叠到真实数据上看形状和尾部对不对得上对得上才说明模型真的描述了数据。一行命令把前两样调出来import arviz as az az.plot_trace(idata, var_names[beta0, beta1, sigma])上图左半是每个参数的可信区间粗横线是后验主体右半是 R-hat 点——都贴着 1说明这条链收敛了。读图比调参重要先看右半收没收敛再看左半区间宽不宽。从一条线到三条线多分位数建模业务里你很少只关心一条线。客服排班看的是中位数风险盯的是 P95清仓促销要看 P10。所以条件分位数建模常常一次画三条好把典型、偏高、偏低一起摆出来。技巧就一招用shape把参数拉成一维向量每条分位线配一套自己的beta和sigma再循环绑定似然。qs [0.1, 0.5, 0.9] # 同时看 10%、50%、90% 三条线 with pm.Model() as mq: beta0 pm.Normal(beta0, 0, 10, shapelen(qs)) # 每条线一个截距 beta1 pm.Normal(beta1, 0, 10, shapelen(qs)) # 每条线一个斜率 sigma pm.HalfNormal(sigma, 2, shapelen(qs)) mu beta0 beta1[None, :] * x[:, None] # 广播成 (n_obs, n_quantiles) for i, q in enumerate(qs): # 循环给每条线配似然 pm.AsymmetricLaplace(fy_{q}, mumu[:, i], bsigma[i], qq, observedy)这里有个反直觉的点三条线不是各画各的它们共享同一套数据MCMC 会把它们之间的相关性一并学出来。画出来你会看到三条线像扇面一样张开——x 小时彼此贴着x 大时上下分位离中位越来越远把越往右越不确定这件事直接描在了图上。一个落地场景电商补货的需求上界回到开头那个断货的坑。补货决策的核心诉求是需求上限不是均值。建模思路分三步。特征怎么选促销开关、折扣力度、星期几、天气这些直接影响当天需求但别急着塞上周销量进去它和当前需求高度相关容易把共线性搅进后验。为什么取q0.95而不是0.90因为缺货和滞销的成本不对称。大促断货会丢流量、招差评多压一点库存只是占用资金代价低得多所以上界要往更保守的那端取。怎么把后验翻译成备货量别拿点估计下单取后验的一个分位当安全备货量。X feats[[promo, discount, dow]] # 特征促销 / 折扣 / 星期 with pm.Model() as om: beta pm.Normal(beta, 0, 5, shapeX.shape[1]) # 每个特征一个系数 alpha pm.Normal(alpha, 0, 10) sigma pm.HalfNormal(sigma, 100) mu alpha X.dot(beta) # 线性预测95% 需求上界 pm.AsymmetricLaplace(demand, mumu, bsigma, q0.95, observeddemand) idata pm.sample()举个具体的数模型给出某天需求上界的后验95% 可信区间是 [800, 1200]中位数 900。如果按均值 900 备货等于赌需求不会超过 900而它超过 900 的概率不小粗算有 25% 左右会断货反过来直接按 1200 备又可能白压三成库存成本。正确做法是取后验的 95% 分位比如 1150当安全备货量把断货风险显式压到你愿意接受的 5%。这就是概率化预测比甩一个数值钱的地方。三个常见坑第一q 贴边。现象是q设成 0.99 甚至 1.0 时参数检查直接报错或后验剧烈漂移。原因是 q 必须严格落在 0 到 1 之间PyMC 内部把它换算成对称参数 kappaq 越接近 0 或 1kappa 越极端分布退成单边、梯度爆炸。避免办法是业务分位留点余地0.95、0.99 都行但别用 1.0真要够尾部就加大 tune。第二把 q 当置信区间。现象是有人直接拿q0.95的 mu 后验当95% 置信区间去汇报。原因是 q 描述的是数据分布的分位即需求落在它以下的概率而置信区间讲的是参数估计本身的不确定性两码事。避免办法是分位线只用来报预测需求的概率上界参数的不确定性另看 trace 图上的可信带别混着讲。第三shape 对齐出错。现象是多分位数一跑就报形状不匹配mu 广播对不上。原因是beta0、beta1用了shape拉成向量但 mu 里索引的维度没对齐或循环里取错了轴。避免办法是先在小数据集上把 shape 打印出来核对确认 mu 是 (n_obs, n_quantiles) 再喂给似然。接下来可以做什么三条能马上做的延伸。第一嵌进 A/B 测试把分位数回归当实验分析的一步直接比较新旧策略在不同分位上的差异而不只看均值提升。第二试非线性x 和分位线关系弯曲时把线性项换成样条或用分层模型同时给多个门店、区域各画一套分位线。第三接时间序列需求有周期和趋势时把 AsymmetricLaplace 似然挂到状态空间模型上让分位线随时间走。想继续深挖两个地方值得一读概率分布指南 里对不对称拉普拉斯的展开和 GLM 线性回归教程 用来对比传统回归和分位数回归的差别。PyMC 把分布、采样、诊断拆成了独立模块你以后想换分布、换采样器都是在这张图里挪一下的事。把这套流程跑通你就已经超过了大多数还只盯着 OLS 残差的人。【免费下载链接】pymcBayesian Modeling and Probabilistic Programming in Python项目地址: https://gitcode.com/GitHub_Trending/py/pymc创作声明:本文部分内容由AI辅助生成(AIGC),仅供参考