一个生物模型有多少个参数？少则三五个，多则几十上百。
捕食者怎么繁殖，猎物怎么被捕食；蛋白质怎么磷酸化，mRNA 怎么降解；感染者怎么传播疾病，康复者怎么获得免疫——每个环节都要用一个数字来描述。这些数字从哪里来？文献里翻的，实验里测的，或者干脆猜的。
问题来了：哪个参数最重要？
直觉上大概会说"都很重要啊"。但计算资源不是无限的，实验经费更不是。如果只有精力精确测定三个参数，你选哪三个？如果模型结果对某个参数的变化毫无反应——花一个月测的那个数字，是不是根本无所谓？
这就是 Sobol’ 敏感性分析要解决的问题。
不只是"敏感"
敏感性分析（Sensitivity Analysis, SA）这个概念不新鲜。最朴素的做法是局部敏感性分析（Local SA）：把每个参数调高 1%，看输出变化多大。听上去很合理，问题在于——你只在一个点上做了测试。换个初始值，排名可能完全颠倒。对于非线性系统，这种局部视角约等于管中窥豹。
Sobol’ 指数属于全局敏感性分析（Global SA）。它不问你"在这个点上谁敏感"，而是问——在整个参数空间里，谁在操控输出的方差？
它的做法很直接：把模型输出 $Y$ 的总方差按来源拆开。
$$Y = \mathcal{M}(X_1, X_2, \dots, X_m)$$$\mathcal{M}$ 是你的模型，$X_i$ 是输入参数。Sobol’ 分解告诉你：
- $S_i$：$X_i$ 自己贡献了多少方差（first-order index）
- $S_{ij}$：$X_i$ 和 $X_j$ 的交互贡献了多少（second-order index）
- $S_{Ti}$：$X_i$ 的总贡献——包括它自己，加上和所有其他参数的两两交互、三三交互……一直加到 $m$ 阶（total-order index）
如果 $S_{Ti} \approx 0$，恭喜你，这个参数可以当成常数，模型对它几乎不敏感。如果 $S_i$ 大但 $S_{Ti} - S_i$ 小，说明这个参数是独狼，自己就很能打。如果 $S_i$ 小但 $S_{Ti}$ 大——这就很有意思了：它自己不重要，但和别人搞在一起就很关键。
这种参数的微妙之处，局部方法永远抓不到。
计算上的麻烦
Sobol’ 指数好看，但不好算。
理论上，你可以用 Monte Carlo（MC）来估算所有方差。问题是，对于 $m$ 个参数，要算到第 $k$ 阶交互，需要的模型评估次数是指数级的。如果你的模型求解一次就要几十分钟——比如一个 stiff ODE 系统——那 MC 的时间成本会让人想关电脑。
Tosin 等人的 tutorial 介绍了一种取巧的方法：用多项式混沌展开（Polynomial Chaos Expansion, PCE）构建一个 surrogate model，然后在 surrogate 上算 Sobol’ 指数。
PCE 的思路是这样的：把模型输出 $Y = \mathcal{M}(X)$ 展开成一组正交多项式的加权和
$$Y \approx \mathcal{M}^{PC}(X) = \sum_{\alpha=0}^{P} y_\alpha \Psi_\alpha(X)$$这里的 $\Psi_\alpha(X)$ 是一组多元正交多项式，选择规则取决于输入参数的分布——Uniform 对应 Legendre 多项式，Gaussian 对应 Hermite，等等。系数 $y_\alpha$ 通过回归（Ordinary Least-Squares）拟合，只需要少量模型评估作为训练样本。
关键在于：一旦 PCE 建好了，Sobol’ 指数的计算变成了——系数的平方和除以总平方和。
$$\hat{S}_u = \frac{\sum_{\alpha \in u} y_\alpha^2}{\sum_{\alpha=1}^{P} y_\alpha^2}$$不需要再跑模型。而且因为只涉及正数的求和，不存在 MC 方法中常见的减法消去误差（cancellation error）——你不会看到负的 Sobol’ 指数这种没有物理意义的东西。
打个比方：MC 在整片森林里一棵一棵数树；PCE 是先画一张大致准确的地图，然后在地图上数。
Tutorial 里有一个很直观的对比：在 Lotka-Volterra 捕食者-猎物模型上，MC 用了 75000 个样本跑了两分钟，PCE 只用 1000 个样本跑了 45 秒——还顺带画了验证图。（当然，两分钟对 Lotka-Volterra 算快的。换个 stiff 的系统，75000 次可能是一个周末和一个周末的区别。）
三个例子，一层比一层深
Tutorial 用三个生物模型展示了整个框架。
Lotka-Volterra 捕食者-猎物模型有四个参数：$a$（猎物出生率）、$b$（搜索率）、$c$（捕食者死亡率）、$d$（能量转化效率）。Sobol’ 分析的结果有点反直觉——最重要的参数随时间变化。
分析捕食者种群时：第一年 $d$ 最重要，第六年变成 $b$，第十年又变成 $a$。
这不奇怪，只是我们平时不太这样想问题。一个有 time-dependency 的系统，不同阶段由不同的机制驱动。如果一个干预措施只能影响一个参数，那你就得想清楚——现在是什么阶段？该碰哪个参数？
NF-$\kappa$B 信号通路模型就复杂多了。Hoffmann 等人 2002 年在 Science 上发表了原始模型，26 个分子物种、64 个反应系数。Krishna 等人后来做了简化，提取出核心反馈回路：7 个物种，12 个参数。这个模型的 QoI 不是某个时刻的浓度，而是峰值持续时间（peak duration）——NF-$\kappa$B 浓度高于其均值的时间长度。
Sobol’ 分析的结果是：最重要的五个参数是 $k_t$（转录率）、$k_{tl}$（翻译率）、$k_{Iin}$（I$\kappa$B 核输入率）、$k_m$（mRNA 降解率）、$k_{Nin}$（NF-$\kappa$B 核输入率）。
但更值得注意的是 $k_m$ 和 $k_{Nin}$——这两个参数的一阶指数都很小，总阶指数却很大。
这意味着什么？它们自己不直接控制输出，而是通过与其他参数的高阶交互来发挥作用。只看一阶指数就会漏掉它们。而漏掉它们，可能意味着漏掉了药物干预的最佳靶点。
SIR 流行病模型是最后一个例子。$S$（易感者）、$I$（感染者）、$R$（移除者），三个仓室，两个参数：$b$（传播率）和 $a$（移除率）。这次他们把初始感染者数量 $I_0$ 也当作一个随机变量。
分析感染人数 $I(t)$ 时：$b$ 和 $I_0$ 主导了整个方差，$a$ 的贡献几乎为零。
等等——移除率不重要？
不是说 $a$ 不会改变输出。重点在于：在所有参数同时浮动的不确定情境下，$b$ 和 $I_0$ 的变化解释了绝大部分输出波动，$a$ 的变化相比之下几乎看不出来。
资源有限的时候，该先搞准哪个参数，不言自明。
我的观察
读了这篇 tutorial 之后，有几个感触。
第一，Sobol’ 分析的妙处在于它把"参数重要性"这个模糊的直觉变成了可计算的东西。没有它，大家对"重要参数"的判断基本靠吵架——做实验的说 $a$ 重要，做理论的说 $b$ 重要。有了数字，至少吵得有依据。
第二，PCE 这个 trick 真是好。MC 跑 75000 次模型在 tutorial 的例子中只要两分钟，但换成我见过的一些计算生物学模型——一次求解几小时的那种——75000 次等于永远。PCE 用 1000 次训练样本建一个近似，误差可控，还能顺手算出所有阶的 Sobol’ 指数。唯一的代价是多学一种数学工具。值。
第三，参数重要性随时间变化的发现，对搞模型校准的人来说应该是一种提醒。很多时候大家只在一个时间点校准参数，然后假设这个校准结果对所有时间都成立。但如果不同时间段的敏感参数不一样，你拿后半段的数据去校准一个只在前半段敏感的参数，等于浪费数据。
第四，也是我觉得最值得警惕的一点：一阶不重要不代表真的不重要。$k_m$ 在 NF-$\kappa$B 模型里一阶指数几乎为零，但不看总阶指数就会错过它——而它和五个其他参数都有显著的二阶交互。这让我想起统计里的"主效应不显著但交互效应显著"。在生物系统这种充满反馈回路和耦合的地方，这种"隐形参数"可能比我们想象的多得多。
SoBioS 与可复现性
Tutorial 的作者做了一个开源库叫 SoBioS——Sobol’ indices for Biological Systems——基于 MATLAB 的 UQLab 工具箱封装而成。所有三个例子的代码都在上面，读者可以自己跑一遍。
这年头发 paper，把代码藏着掖着已经越来越说不过去了。像这种 tutorial 类的文章更是——如果读者不能自己动手试，那还算什么 tutorial？
不过话说回来，我自己也没有真的去跑。（实在一点也不想用MATLAB）
模型的参数是"积木"，但积木和积木不一样——有些是承重墙，有些只是装饰。知道谁是谁，模型才不是一个黑箱。