偏相关、高斯图模型与收缩估计

偏相关去除的是什么

偏相关在控制其他变量后度量两变量的关联。在多元正态假设下，偏相关为 0 等价于条件独立，非零偏相关构成高斯图模型（GGM）的边。偏相关与精度矩阵的关系、它在链式结构与对撞结构中的表现，已在连接强度的估计和偏相关与高斯图模型中讲过，这里只强调基因表达数据的特点：$p\gg n$。

三基因算例。只控制一个变量时，一阶偏相关有闭式公式：

$$r_{13\cdot2}=\frac{r_{13}-r_{12}r_{23}}{\sqrt{(1-r_{12}^2)(1-r_{23}^2)}}.$$

设 1→2→3 为调控链，$r_{12}=0.8$，$r_{23}=0.7$。若 1 与 3 之间没有直接作用，期望 $r_{13}=0.8\times0.7=0.56$——相关网络会画出 1–3 这条边；但 $r_{13\cdot2}=(0.56-0.56)/\cdots=0$，偏相关正确地去掉了它。若观测到 $r_{13}=0.7$，则 $r_{13\cdot2}=0.14/\sqrt{0.36\times0.51}\approx0.33$，提示 1 与 3 之间在 2 之外还有联系。

$p>n$ 时的问题

全偏相关需要协方差矩阵的逆。样本协方差 $S$ 的秩至多为 $n-1$，当 $p\ge n$ 时不可逆；即使 $n$ 略大于 $p$，$S$ 的小特征值被严重低估，求逆后噪声被放大。常见对策有两类：收缩（把 $S$ 拉向一个简单目标，保证可逆且条件数好）与稀疏正则化（图形 LASSO、邻域选择，见正则化估计：EBICglasso）。基因网络中，Schäfer 与 Strimmer（2005）的收缩方法因为有解析的最优收缩强度、计算快，而被广泛使用（R 包 GeneNet、corpcor）。

Schäfer–Strimmer 收缩估计

对标准化数据（即相关矩阵 $R$），以单位阵为目标：

$$R^{*}=\lambda I+(1-\lambda)R,\qquad\text{即 } r^{*}_{ij}=(1-\lambda)r_{ij}\ (i\neq j).$$

使均方误差最小的收缩强度为

$$\lambda^{*}=\frac{\sum_{i\neq j}\widehat{\operatorname{Var}}(r_{ij})}{\sum_{i\neq j}r_{ij}^{2}},$$

其中 $\widehat{\operatorname{Var}}(r_{ij})=\dfrac{n}{(n-1)^3}\sum_{k=1}^{n}(w_{kij}-\bar w_{ij})^2$，$w_{kij}=z_{ki}z_{kj}$ 是第 $k$ 个样本中两基因标准化值的乘积。直观地说：相关系数整体越“不可靠”（分子大）、相关整体越弱（分母小），就收缩得越多。$\lambda^*$ 截断到 $[0,1]$。得到 $R^*$ 后求逆，按 $\rho_{ij}=-\omega_{ij}/\sqrt{\omega_{ii}\omega_{jj}}$ 换算偏相关。推导见数学推导与算法。

收缩只给出非稀疏的偏相关估计，还需要判断哪些边显著。Schäfer 与 Strimmer（2005, Bioinformatics）提出用经验贝叶斯混合模型拟合全部偏相关的分布（大部分来自“零边”的零分布，少部分来自真边），对每条边给出局部错误发现率（local fdr），以此选边。

模拟中的表现

代码示例的“收缩估计偏相关”单元在一个含 8 个枢纽的稀疏 GGM（80 个基因）上比较 $|r|$ 与收缩偏相关：

（随机预测的 AUPR 为 0.02。）结果说明两点：其一，同一枢纽的两个邻居之间没有直接边却有明显相关（$n=400$ 时平均 $|r|\approx0.21$），偏相关把它们压到约 0.04，这是偏相关的价值；其二，样本很少时偏相关并不总是优于相关排序——估计一个精度矩阵比估计相关矩阵需要更多信息，收缩又把信号压得很弱。$n$ 增大后偏相关继续接近真值，相关排序则停在间接边造成的上限附近。