基因差异表达分析是生物标记物研究中最重要的统计分析之一，其结果的可靠性直接影响下一步实验研究的结果。然而，现有的蛋白质组数据中的基因差异表达分析，并没有考虑到蛋白组数据中特有的数据结构。同时对于是否保留单一多肽支持的蛋白进行统计分析并没有统一的规范。

针对这一难题，近日，瑞典斯德哥尔摩卡罗林斯卡研究所Janne Lehtio团队分享了一种新的蛋白质组学基因差异表达分析的统计方法：DEqMS，使得统计检验结果更加准确。相关研究结果发表于蛋白组学领域Top期刊《Molecular & Cellular Proteomics》上。

我们今天特别邀请到本文的第一作者，朱亚锋博士（现在哈佛医学院进行博士后研究）分享该研究的主要内容。

撰文 | 朱亚锋 博士

基因差异表达分析是生物标记物研究中最重要的统计分析之一，其结果的可靠性直接影响下一步实验研究的结果。蛋白质组数据中的基因差异表达分析一直以来都在沿用经典统计方法如t-test, ANOVA 和线性混合模型，或是之前为转录组数据开发的方法—Limma。然而这些方法并没有考虑到蛋白组数据中特有的数据结构。

另外，对于是否保留单一多肽支持的蛋白进行统计分析并没有统一的规范。有些研究为了降低统计检验结果的假阳性，会提前去掉这些单一多肽支持的蛋白，主要是由于其定量结果波动性很大。这些蛋白往往占到所有蛋白的10%-20%，且大部分是一些低表达的蛋白，很有可能参与着重要的生物调控过程，因此很有必要保留这些蛋白。

针对这一问题，我们开发了一种可用于蛋白质组学基因差异表达分析的统计方法——DEqMS。此方法基于在蛋白质组数据中观察到的基因方差与用于基因定量的多肽数量的相关性（如图1所示），应用贝叶斯定理，为不同数量多肽支持的基因估算一个更准确的先验方差（prior variance），来校正实际方差，使得统计检验结果更加准确。

图1: 蛋白质组数据中观察到的基因方差与用于基因定量的多肽数量的相关性。A) TMT标记蛋白组数据中，同组内生物学重复样本间的标准差随着蛋白PSM数量增加而降低。B) 红色和灰色曲线是分别使用DEqMS和Limma计算的先验方差(prior variance)。Limma为所有蛋白计算一个共同的先验方差，DEqMS计算的先验方差随着蛋白PSM数量增加而降低。图中散点代表蛋白的合并样本方差(pooled variance). A和B图中纵坐标经过自然对数ln变换。

DEqMS 方法的灵敏度和假阳率测试

首先我们测试了此方法对于统计检验后p-value 的影响。使用DEqMS 后得到的 p-values比t-test和Limma这两种传统方法更小(如下图所示)，可见DEqMS检测差异表达基因的灵敏度更高。

为了测试DEqMS在增加灵敏度的同时是否会增加假阳率，我们生成了一个由9个样本，6000个蛋白构成的模拟数据。并按照真实数据中的比例，将6000个蛋白分30个组，分别有1-30个PSMs定量。每组里的蛋白平均方差均借鉴真实数据中的拟合值。

在已知此模拟数据中没有任何差异表达蛋白的情况下，随机选择6个样本分为两组，使用t-test, Limma 和DEqMS进行差异表达分析 (3 vs 3)。在所有84种随机组合中，计算每次得到的FPR(false positive rate), 如上图boxplot所示，DEqMS的FPR中位值比Limma低，如预期中一样，t-test的FPR最低。

基准测试（Benchmarking）

我们进一步使用label-free 和 labelled 的定量添加实验数据, 来比较不同统计检验方法的准确性。label-free数据采用的是Matthias Mann实验室发表的定量添加实验[2]，在两种条件下（各有3个重复样本），10 μg 和30 μg E.coli 蛋白提取物被添加到等量 50 μg Hela细胞蛋白提取物中。在TMT10 标记的数据中，7.5 μg , 15 μg 和45 μg E.coli 蛋白提取物被添加到等量 70 μg Hela细胞蛋白提取物中，其重复样本数量分别是3，4，3。在label-free 和 labelled 的定量添加实验数据, 各有6566个 (1902 E. coli, 4664 human) 和11188个(2474 E. coli, 8764 human) 蛋白被检测和定量。随后比较不同统计检验方法的差异表达分析结果。使用两种方式对比，一个是ROC （receiver operating characteristic）曲线，和在 1% FDR条件下得到的true positive 和 false positive 的数量（如下图所示）。

图2: 基准测试结果。不同方法括号里显示的是partial area under curve (pAUC, in the range specificity > 95%) 除以理论最大值得到的百分比。

如图所示，在label-free和labelled 数据中，DEqMS的结果更好，或者至少与其他方法一样好。传统方法t-test虽然FPR最低，但是会同时丢失很多true positive。 比较在1%FDR条件下得到的true positive 和 false positive 的数量，两组数据中，DEqMS结果略微优于排在第二的Limma方法。

真实数据测试

由于以上定量添加实验中，所有蛋白的fold change 比较大（2-3 倍）且都相同，并且组内样本间只存在技术性差异(technical variation)，并无生物学差异(biological variation)。为了进一步测试不同方法，采用了一个已发表的TMT10蛋白组数据(包含8625蛋白的定量数据)【3】。此实验中，U1810 cells 经过了三种不同的microRNA mimic 处理，且同时产生了RNA-seq 数据。通过比较经miR-372处理的实验组与对照组(各有3个生物重复样本)，DEqMS ， Limma 和 ROTS三种方法分别产生了 201，120 和21个 显著差异表达的蛋白。在DEqMS 与 Limma 共同发现的109个差异表达的蛋白中，有30%的基因在RNA-seq数据中显示相同的调控方向且被预测为miR-372的靶基因。DEqMS单独发现的92个基因中，有同样的数据支持。然而11个仅在Limma发现的基因中，并没有此数据支持。

为了进一步解释DEqMS结果优越性的原因。我们画出了两种方法在估算方差时的残差平方和(residual sum of squares)，以及得到的后验方差(posterior variance)。Limma方法估算方差时的残差平方和明显大于DEqMS（下图B和D），且残差呈现非随机的分布，而DEqMS估算的残差则是随机分布在0左右。比较后验方差（下图C和E），DEqMS估算的后验方差向中心“收缩”，避免了数据中个别蛋白偶然出现的极低方差（例如C图中的红色圆圈标记出11个仅在Limma发现的差异表达蛋白）。

最后此方法还在临床蛋白组数据和磷酸化蛋白组数据中测试，发现基因方差与用于基因定量的多肽数量的相关性依然存在，说明DEqMS方法也同样适用于此类数据。值得一提的是，我们发现，当同组样本间的生物差异较大（如临床样本），DEqMS的优越性会逐渐减少。此方法以R package的形式，在Bioconductor共享。

参考文献

1. Yafeng Zhu, et al., 2020, DEqMS: A Method for Accurate Variance Estimation in Differential Protein Expression Analysis. Mol.Cell. Proteomics..

2. Cox, J., et al., (2014) Accurate proteome-wide label-free quantification by delayed normalization and maximal peptide ratio extraction, termed MaxLFQ. Mol.Cell. Proteomics.

3. Zhou, Y., et al., (2017) microRNAs with AAGUGC seed motif constitute an integral part of an oncogenic signaling network. Oncogene.