基于半参数层次混合模型的经验贝叶斯最优发现方法

摘要

多重检测已被广泛应用于全基因组研究，如微阵列实验。为了在这些全基因组研究中进行有效的基因选择，最佳发现程序(ODP)被Storey (英国皇家统计学会杂志，B辑，69卷,没有。3, pp. 347-368, 2007)。在本文中，我们开发了一种经验贝叶斯方法来实现基于半参数层次混合模型的ODP使用“通过粗糙平滑”的方法。在半参数层次混合模型下，先验分布可以灵活建模，ODP检验统计量和后验分布在分析上易于处理，计算容易实现。此外，我们在经验贝叶斯框架中提出了基于错误发现率(FDR)的显著性规则。介绍了两个临床研究的应用。

1.介绍

全面的基因表达微阵列是从大量候选基因中筛选不同表型(如临床亚型和疾病预后分类)之间差异表达基因的强大工具。基因筛选研究有可能有助于阐明疾病生物学和侵袭性，识别新的治疗靶点，并开发新的分子诊断，为个体患者优化药物[1.–4.]然而，微阵列数据的高维性在通过基因筛选分析从大量噪声变量中提取少量相关基因方面提出了特殊挑战。除了控制假阳性外[5.–7.，提高基因筛选分析的效率是很重要的。

为了有效筛选差异表达基因，Storey[8.]开发了最佳发现程序（ODP），可以将其解释为对单个假设检验的最强大检验的扩展[9到多次测试。Storey为多次检测定义了一个最优标准，最大限度地增加每个固定水平的预期假阳性(EFP)的预期真阳性(ETP)数量[8.］．ODP是作为一个测试程序开发的，实现了这种优化;它通过在多个测试中“借力”提高了多个测试的整体性能。然而，在应用ODP时，必须估计以下两个分量:(a)每个显著性检验(null或alternative)的真实状态和(b)每个检验对应的真实概率分布[8.］．为了解决这些估计问题，Storey等人[10构建广义似然比统计量[11，使用最大似然估计。Cao等[12]还提出了一种贝叶斯方法，Woo等人[13]开发了一种计算效率高的方法，称为“模块化ODP”

最近，Noma和Matsui [14]为ODP开发了一种经验贝叶斯方法，可以有效地规避(a)和(b)的估计问题。在一个分层的随机效应模型下，他们表明ODP是作为基于边际似然比统计量的测试规则推导出来的。在基于仿真的数值评价中，经验贝叶斯方法几乎达到了ETP的理论界限(即平均功率)，与现有方法相比，在大范围条件下表现良好[14］．

这种经验贝叶斯方法的一个关键方面是，在层次模型中必须指定一个完全参数形式。尽管之前在微阵列研究中讨论的大多数经验贝叶斯方法[15–21[假设先验分布的完全参数自然共轭模型，这一假设可能是限制性的，因为复杂的分子途径导致信息基因的效应大小可能存在很大的差异，而在微阵列研究中，这些差异的分布通常是未知的。在本文中，我们开发了一种基于半参数层次混合模型的ODP经验贝叶斯方法，其中效果大小的先验分布不是参数指定的。我们通过应用“通过粗糙平滑”的方法估计先验分布的非参数分量[22,23为ODP统计提供了一种有效的计算方法。

本文的组织结构如下2.，我们提供了一个关于Storey的ODP方法的理论结果的简要概述[8.]野马和松井[14］．节3.，我们描述了半参数层次混合模型和如何实现ODP方法。前列腺癌和淋巴瘤的临床研究4.．最后，我们将在本节中进行讨论5.．

2.多重显著性检验的最优发现程序

2.1。Storey的最佳发现程序

在本节中，我们简要回顾了ODP方法的理论结果。我们表示来自的数据集基因作为,在那里．在微阵列的研究中,对应的一组日志转换，规范化的基因表达样品的基因(见章节4.)此外，假设数据集是公共概率空间中定义的随机向量。识别差异表达基因的最佳目标是找到尽可能多的真阳性，而不会产生太多的假阳性[10］．层(8.]将多重检验的最优准则定义为使ETP在某一固定EFP水平上达到最大的统计显著性规则。

他导出了下列引理，给出了达到最优准则的多重检验过程。

引理1 (Storey [8.])。假设第一个测试有一个零密度和替代密度。此外，假设无效假设对th测试而另一种选择则是正确的，不失一般性。以下意义阈值统计数据达到ODP标准： 给定一个固定截止为达到可接受的EFP水平而选择的如果且仅当才被拒绝测试．

应该注意的是，ODP统计器由所有密度函数组成测试。换句话说，通过相互借鉴对单一假设检验而言，总体能力比传统最强大的检验能力提高更多[9,11］．并且，无论数据集之间的相关结构如何，都能得到最优结果．

如本节所述1.，供申请，有两个实际问题：（a）估计每个假设（无效或替代）的真实状态，对应于假设的分子和分母， (b)真实概率密度函数的估计和的,评估。Storey等[10]提供了一种实用的估计方法，其动机是单一显著性检验的广义似然比检验[11］．他们注意到重要性规则基于 等于在引理1.．为每个测试,让为零密度函数的估计，在零假设约束下，将所有未知参数替换为它们的最大似然估计，令使用无约束最大似然估计，对替代假设进行类似估计[10]建议使用规范插件估计: 在哪里是对现状的估计吗假设。估算， Storey等[10]根据单变量统计量(例如-用于两类比较的统计数据)，并估计为当第th次测试低于Storey估计的真零值的比例[6.]，及否则。其他估算程序Cao等人也提出了这些建议[12]和woo等。[13］．

2.2.Noma和Matsui的经验贝叶斯方法[14]

作为storey [8.]指出，ODP类似于众所周知的多点估计的收缩估计[24,25］．因为收缩估计被解释为经验贝叶斯估计[26–28]，层次模型将是跨测试共享信息的自然公式，无论是在推导ODP方面，还是在为实施它而开发估计方法方面[14］．

假设的概率密度函数具有相同的参数形式,在那里感兴趣的参数是和吗是的参数测试。我们假设参数的先验分布如下(）： 在哪里和是先验分布的超参数。

引理2 (Noma and Matsui [14])。在先验分布的经验贝叶斯框架下(4.），以下意义阈值函数实现了ODP标准： 对于固定截止线选择达到可接受的EFP水平，为零假设如果且仅当才被拒绝测试．

参见野马和松井[14解释;解释．

用于说明先验分布和，野玛和松井[14]采用经验贝叶斯方法[26,27］．使用经验贝叶斯实现来获得估计和从数据，然后可以代入ODP统计，, 在哪里和估算值是多少和分别地

与现有的ODP估计方法相比，基于经验贝叶斯方法有两个优点：第一，它只估计超参数和对于随机效果模型，而不是密度函数中的参数，必须按照Storey等人的程序估计大量参数（微阵列研究中的数千个或更多）[10］．其次，在经验贝叶斯方法中，对每个测试的真实状态(null或alternative)的估计是不需要的。

3.半参数层次混合模型与ODP

3．1.半参数层次混合模型

在本节中，我们考虑将ODP应用于具体的两类比较的微阵列分析(分为“0”和“1”两类)，例如，预后好和预后差，根据表达水平基因样本。这里考虑的基因表达数据包括来自双色互补DNA阵列的归一化对数比或来自寡核苷酸阵列(如Affymetrix genchip)的归一化对数信号。让和分别为0类和1类的样本量。让和分别为第0类和第1类的基因表达水平的平均值。我们认为差异表达基因的检测是在．也就是说，测试问题

我们假设基因的表达Th基因正态分布在每一类中， 在哪里是常见的类内方差()，如Storey等人[10］．让为两个类的平均值的差th基因(），，这是我们感兴趣的参数。让,分别表示类0和1的基因表达水平的样本方法。 使用此转换，将生成讨厌的参数的模型中消除了和．具体地说, 在哪里．请注意(7.)相当于测试 让和．

将ODP应用于测试问题在模型(10)，我们假设模型参数的随机效应分布如下(）： 在哪里在第一个组成部分为狄拉克函数，表示两个类之间的非微分表达式。第二部分表示微分表达式的一个组件，和和各成分的混合比例是否一致．基因特异方差的随机效应分布为便于解析处理，指定为参数共轭模型。的建模，而不是像Noma和Matsui那样采用全参数共轭正态分布模型[14]，我们采用Laird和Louis的“粗化平滑”方法[22]及沈和路易斯[23］．该方法用于获得经验贝叶斯分析中先验分布的非参数最大似然估计[29,30.]。该方法从先验分布的平滑估计开始，并使用EM算法将其“粗糙化”为非参数最大似然估计。用于建模，我们在一系列非零点处处于同样间隔的离散质量点的支持，以获得足够大的.我们假设 在哪里．我们将微分分量的分布函数的超参数表示为．沈和路易[23]给出了在数值计算的基础上实现粗化平滑方法的详细指导。沈和路易[23]建议网格点的数量至少为200。为了简单起见，他们建议网格点可以在分布的支持下等距分布。增加网格点的数量并根据先验密度将它们隔开将提供一个更好的近似因此会改进估计，特别是对尾部区域。然而，太多的数字不一定是有益的。有关分析和数值计算的更多细节，请参阅Shen和Louis的指南[23］．在本文中，我们将这种层次模型称为“半参数层次混合模型”。该模型的重要性质是后验分布和我们提供了截面中后验分布的解析表示3.4.．

超公数的最大可能性估计，由EM算法得到[31］．因为后面的分布和， EM算法的e步易于实现。EM算法的详细内容在附录中给出。

３．２．最优发现程序

ODP数据对于测试问题(11)为边际似然比。原假设的边际似然对应于非微分分量的边际似然: 同理，(11)对应于微分分量: 在半参数层次混合模型下。

因此，ODP统计量(边际似然比)为 通过插入的最大似然估计在此基础上实现了经验贝叶斯检验过程 注意，ODP统计量是在封闭形式下使用的最大似然估计．

3.3。评估意义

为了评估显著性，我们可以得到FDR的经验贝叶斯估计[5.–7.]的半参数层次混合模型。对于任何给定的阈值时，表示显著区为每个测试。让表示显著基因的一组指数。在贝叶斯选择规则框架下，FDR符合差分/非差分分类的误分类错误率[6.,7.,32]，以及罗斯福总统著名的贝叶斯估计[33,34] 是 在哪里 相应地，本文提出了一个经验贝叶斯估计对于重要区域可以通过插入最大似然估计来获得．

请注意，我们进行了小型模拟，以检查拟议的ODP程序的正确性，这些结果在可在线获得的补充材料中提供http://dx.doi.org/10.1155/2013/568480．

3．4．后验分布

在半参数层次混合模型下，后验分布和如上文所述，可通过分析获得；这为实施贝叶斯推断提供了计算优势。后验分布将用于推导估计值和置信区间的统计显著性指标之外，还包括-价值观。

让和指示指示变量，该变量Th基因组分属于．关节后部分布表示为差分和非差分成分的混合物： 在哪里．重量在第节中给出3.3.,．微分分量的后验分布表示为 在哪里 的边缘后验分布获取方式如下:

4.应用程序

我们使用来自两个微阵列临床研究的真实数据集来说明我们的方法。因为之前的模拟研究表明，Storey的方法[10比许多其他多种测试方法更强大[10,12,35，我们选择这个方法作为参考。

4．1.前列腺癌的例子

埃夫隆(36,37]分析了来自Singh等人的修饰过的基因表达数据集[38)与基因样品；前列腺癌患者健康对照组。对于提出的经验贝叶斯方法，我们放置网格点为等距的上* 0.01(除0;）.超参数估计为,,,，以及如图所示1.．估计的先验分布是偏态的和多模态的。这表明，在这种情况下，分析上易于处理的参数模型(例如，正态分布)可能是不充分的。

数字2.通过Storey等人的方法介绍多次测试的结果。[10]和本文发展的经验贝叶斯方法。这张图总结了重要基因的数量和- 值[6.,7.]的变化值的截止．有关公平比较，请参见相同的FDR控制方法3.3.应用于这两种方法。对于任何层次的-值时，所提方法识别的显著基因数量大于Storey等[10)方法。

4．2．淋巴瘤的例子

戴夫等．[39分析了Affymetrix HG-U133A和U133B (Affymetrix, Santa Clara, CA)芯片数据，来自未治疗的滤泡性淋巴瘤患者的191个活检标本．该数据可从BRB-ArrayTools人类癌症基因表达数据档案[40］．我们比较了5年内死亡的患者(预后差;)，生存5年以上(预后良好;）.

对于半参数先验分布的极大似然估计，超参数估计为和．的估计先验分布如图所示3.．(网格点被放置在上面除0外，等距0.01;。）与前列腺癌数据相比，大量估计都位于小区域吗因此，数据中包含的差异表达信号将比前列腺癌数据中的信号弱4.，重要基因的数量和它们之间的关系值。经验贝叶斯方法识别的显著基因数量大于Storey等人的方法[10］．

5.讨论

对于使用高维微阵列数据进行高效基因筛查，多种检测方法是控制假阳性结果的有效工具。尽管多重检测可以提供一个相关的框架，以确保对一组重要基因(如FDR [5.–7.]），研究人员想找到尽可能多的真正阳性基因，以进一步调查。由Storey引入的最佳发现程序（ODP）[8.，将为这一需求提供一个有效的解决方案。Storey等[10]在许多数值评估中提供理想的结果[10,12,35］．

本文通过放宽Noma和Matsui的经验Bayes ODP方法中的参数先验假设，提出了一种基于半参数层次模型的经验Bayes ODP方法[14］．层次混合模型和经验贝叶斯方法形成了基因间信息共享的基础，并为有效的多重测试提供了框架[15–21］．尽管经验贝叶斯方法比引出单个先验更健壮，但当真实先验不是来自假设的参数族时，参数经验贝叶斯方法会缺乏健壮性，正如Morris指出的[28，例如未能合并混合结构或未能指定效应-尺寸分布的形式．实际上，前列腺癌的估计先验分布(Section4．1）有多峰和倾斜的形状。在这样的实例中，天然缀合物家族可能不充分，以对未知的随机效应分布进行建模。“平滑逐渐粗糙化”方法[22,23是经验贝叶斯推理中灵活的建模方法之一[41–43]我们的方法的一个扩展是为随机效应分布指定一个非参数分布采用粗化平滑法，虽然估计大量参数需要大量计算，且参数估计不稳定。

此前，已有许多高效的基因选择方法被提出用于微阵列研究，而不仅仅是多种检测方法，如贝叶斯排序方法[20.,21,43］．为了评估这些方法的实际价值，全面的数值研究，包括大规模模拟和对许多真实数据集的应用，将是一个有价值的课题。

附录

A. EM算法

超参数的极大似然估计采用EM算法得到[31].关于年代,年代,年代,S为缺失变量，则日志完整数据似然表示为 总和中的第一个和第二个术语不依赖于普遍的参数。因此，每个目标函数-步为 在哪里,．上标表示迭代的次数。优化问题的求解-step表示为 同时，对和等效于反伽玛分布的形状和尺度参数的最大似然估计。它可以简单地用标准软件中提供的任何数值优化方法(例如，nlm或R中的optim)求解。

承认

本研究得到了日本科学促进会(24800081)的资助。

补充材料