英语原文共 16 页，剩余内容已隐藏，支付完成后下载完整资料

用P-样条和混合效应模型来聚类纵向截面以及应用于时间过程基因表达数据

N. Coffey， J. Hinde，E. Holian

摘要

由于纵向数据变得越来越普遍，因而提出了各种方法来分析这类型的数据。对时间过程基因表达谱进行研究，其中聚类分析在鉴别随时间共表达基因组中起着重要作用。许多程序已经被用来聚类时间过程基因表达数据，然而先前描述的技术存在许多限制。因而提出了一种替代方法，目的是减少其中的一些限制。该方法利用线性混合效应模型和P样条平滑之间的连接，来同时平滑基因表达数据，目的是消除任何测量误差/噪声，并使用混合效应模型的有限混合来聚类表达谱。这种方法具有许多优点，包括减少计算时间以及在标准软件包中更容易实现。

关键词：纵向面 时间过程基因表达 聚类 混合效应模型 有限混合模型

1.介绍

多元聚类方法如k均值聚类（Hartigan和Wong，1978），层次聚类（Eisen等，1998；Spellman等，1998），自组织图上的聚类（Kohonen，1997；Tamayo等，1999），有限混合模型（FraleyandRaftery，2002；McLachlan等2002，2003，2006），模糊平均聚类（FutschikandCarlisle，2005）和紧密聚类（Tseng和Wong，2005）等都有利于降低基因表达数据的维数，并鉴定共表达基因和共调节基因组（由于具有相似表达谱的基因更倾向于受相同的调节机制控制）。参考Thalamuthu等（2006）和Hand 以及 Heard（2005）对众多聚类技术的讨论与比较。时间过程基因表达研究涉及随时间重复测量数千个基因的表达水平，因而导致极高维的数据。上述方法将每个基因的测量序列作为不同点的载体处理，因此序列元件的任意排列不会影响聚类结果。然而，在时间过程基因表达研究中，需要着重考虑数据的时间排序以及获得的相应聚类。此外，时间过程基因表达数据表现出诸如缺失值，不等的取样时间或较大的测量误差等问题。上述提到的许多技术在处理缺失值方面存在困难，并且对所有基因需要相同的取样时间，又或者不能解释随着时间的推移在相同基因上进行的测量之间的相关性。这促进了诸如CAGED（Ramoni等，2002），隐马尔科夫模型（Schliep等，2005），贝叶斯混合模型（Wakefield等，2003），线性混合效应模型 （Celeux等，2005；Ng等，2006；Qin和Self，2006；Nueda等，2007），基于形状相似性的聚类（Hestilow和Huang，2009）和使用SOM的时间数据聚类（Chen，2009）等技术的发展。然而，这些技术不便于从测量数据中去除噪声，因而忽略在基因表达谱中任何明显的平滑性。因此，基于曲线的聚类方法已经被用于聚类时间过程基因表达数据。这样的方法假定随时间的基因表达是连续过程，因此可以由连续平滑曲线或函数表示。将表达谱作为连续函数可以确保缺失值和不规则采样的数据被适当地处理，同时可以使用平滑技术去除测量误差，并且还可以考虑在相同基因上随时间推移的测量之间的相关性。

在时间过程基因表达分析中基于曲线聚类的一些早期的例子，如Bar-Joseph 等（2003），Luan和Li（2003）和James以及Sugar（2003）使用基函数（例如三次样条，B样条等）的线性组合对每个簇中的平均表达谱进行建模，并对估计的基函数系数进行聚类。其他作者，如Kim等（2008）使用傅立叶基函数的线性组合来表示聚类的表达谱，Kim等（2008）基于傅里叶级数的导数系数聚类。Song等（2007）使用基函数扩展的函数型主成分（FPC）确定了基于FPC评分的聚类，而Leng和Muuml;ller（2006）使用FPC的线性组合表示表达谱，并进行函数型逻辑回归分类，从而将表达谱分类成组。然而，为了估计群均值曲线，这些论文中使用的方法需要选择基函数的数量：K和这些函数的连接点：结。选择基函数（或等效的节点）的数量的最佳值是复杂的问题，并且难以对应用于数据的平滑程度进行控制。使用太多的基函数导致估计的表达谱中过拟合和缺乏平滑性，而使用太少导致欠拟合并且过度平滑数据。一个解决方案是使用平滑样条回归，其中在每个唯一的时间点放置结，并且通过添加一个惩罚项来优化控制的标准。在统计文献中，常见的方法是惩罚表达式曲线估计的曲率（即积分平方二阶导数）（参见Wahba，1990；Green和Silverman，1994；Ramsay和Silverman，2005，的全部细节）。通过平滑参数lambda;对数据的拟合和平滑度之间的权衡进行控制，并且可以使用诸如交叉验证（CV）或广义交叉验证（GCV）的方法来选择lambda;的最优值。在时间过程基因表达研究中，Ma等（2006）使用平滑样条回归聚类基因表达谱，称为他们的方法为SSClust。Ma（2008）将其扩展到贝叶斯算法，Ma和Zhong（2008）在聚类算法中添加了额外的协变量。Deacute;jeanet 等（2007）在进行基于离散化导函数的主分量分数的聚类之前，使用平滑样条回归来估计基因表达谱的导数。Tarpey（2007）和Kayano等（2010）通过将非正交基函数转换为正交基函数（分别使用SVD和Clesles分解），并分别与k均值和SOM相结合，以确保基函数系数的聚类与原始数据的聚类相当。使用平滑样条回归的聚类的主要缺点是与这些方法相关的高计算量。由于平滑样条通常使用与唯一时间点相同数目的基函数，并且优化准则需要数值评估与惩罚项相关联的积分，所以随着样本大小或唯一观察时间的数目增加，计算可能变得麻烦。此外，平滑样条聚类需要为每个聚类选择lambda;的最优值，这也增加了拟合模型的计算负担。

相比之下，本文使用如Eilers和Marx（1996）以及Ruppert等（2003）讨论的惩罚样条曲线（P样条）平滑，以模拟每个簇中的基因表达谱。P样条平滑是低秩平滑法，其需要使用基函数的数量K相对较大，但仍小于遇到的唯一时间点数。如前所述，当基函数的数量较大时，平均表达谱的无约束估计将导致数据的过拟合和过量噪声的拟合。P样条平滑保留所有的基函数，但是使用对估计系数的离散惩罚来限制它们的影响。同样，对数据拟合和平滑度之间的权衡是由平滑参数lambda;控制。使用P样条表示平滑问题来减小数据的维度，从而减少了计算负担。根据Ruppert（2002），一旦使用足够的基函数，P-样条也对所选择的基函数的数量相对不敏感。第2.1节给出了选择节的数量和位置的方法。因为惩罚项是离散的，P样条也很容易计算，而不是像平滑样条的连续积分。此外，将受惩罚的平滑问题表示为线性混合效应模型具有许多附加的优点。在聚类背景下，将平滑问题写为混合效应模型提供了一个框架，用于同时确定每个聚类中平均表达谱的平滑估计，并通过使用额外随机的方法确定群内基因特异性表达谱的估计效应（例如每个基因的随机截距），包括额外的协变量，以及使用混合效应模型混合的聚类表达谱。平滑参数lambda;的最佳值可以通过受限最大似然（REML）自动选择，从而进一步减少计算量，并且可以使用标准统计软件包来实现模型拟合。本文使用线性混合效应模型表示的P-样条平滑来聚类时间过程基因表达谱。虽然该方法是在时间过程的基因表达数据的背景下提出的，但它可以应用于需要聚类分析的任何纵向截面数据集。

本文的其余部分概述如下。第2节描述了如何将特定基因的原始时间过程基因表达数据表示为平滑曲线，并讨论如何将P样条平滑问题实现为线性混合效应模型。第3节讨论使用有限混合模型和混合效应模型的模型混合来同时平滑和聚类时间过程基因表达谱。第4节介绍了进行模拟研究的结果，以检验所提方法的性能。第5节将该方法应用于大型时间过程基因表达数据集。可以根据相应作者的要求提供R代码。

2.方法

2.1平滑

随时间变化的单个基因的基因表达可以被认为是来源于一个平滑的基函数。然而，基因表达数据通常基于一个离散的时间点测量，且常受到大量的测量误差(噪声)影响。因此，假定观测数据可以建模为

(1)

其中是在时间点评估的平滑表达谱的值，并且是测量误差。非参数平滑技术提供了从原始数据去除测量误差和估计平滑表达谱的灵活方式。平滑包括将未知谱表示为一个由K个基函数组成的线性组合

其中是要估计的基函数系数。有许多类型的基函数可以使用，然而，许多应用程序使用样条基函数。样条基包括在时间间隔中的指定点处平滑连接的分段多项式的集合。这些连接点被称为节点，并由表示。基函数K的总数取决于所使用的结数。为了演示目的，描述了第p个自由度的截断幂基，虽然该方法可以扩展到包含其他基函数，如B样条基函数。

第p个自由度的幂基函数由以下基函数组成

其中p是幂的自由度（即p = 1意味着一次截断幂基，p = 2意味着二次截断幂基等），t表示时间，表示第m个节，M是使用的结的总数，并且为截断线性基函数。增加节数M即增加截断线性基函数的数量，这导致更灵活的拟合。基函数的总数取决于节数。

基于第p个自由度的截断幂基的平滑基因表达曲线可以表示为

（2）

其中是第j个观察时间，并且对确定基函数系数的估计具有重要性。

估计基函数系数，(2)可以写成一个线性回归模型

（3）

y包含原始表达式的值，设计矩阵X的形式为

回归系数的矢量

并且假设。

系数可以以多种方式来估计。第一个需要选择基函数数量K以及最小化残差的平方和

（4）

然而，如第1节所述，选择K的最优值是一个复杂的问题，缺少选择K的合适方法。该问题的解决方案是使用平滑样条回归，其中基函数的数量被设置为等于唯一观察时间的数量（即K = n），并且通过将惩罚项添加到优化标准来控制过拟合。然后通过最小化惩罚项的残差平方和来获得beta;的估计

（5）

其中表示的二阶导数，lambda;是控制拟合数据与平滑度之间的权衡的平滑参数。可以通过诸如交叉验证（CV）或广义交叉验证验证（GCV）（Ramsay和Silverman，2005）的方法来选择lambda;的最佳值。由于基函数的数量通常等于唯一观察时间的数量，并且拟合过程需要积分项的数值估计，所以随着样本大小和或观察数量次数增加，基于公式（5）可能产生较大的计算量。为了解决这个问题，采用P样条平滑来替代，其中选择相对大的基函数的数量K（其中K lt;n），并且通过使用岭惩罚如常数来限制截断线性基函数的影响从而考虑过拟合。一旦所使用的数量足够多，P样条平滑对于所选择的基函数（或等效的节点）的数量将相对不敏感（Ruppert，2002）。Ruppert（2002）还提出了一个选择结的数量的经验法则：，其中结被放置在数据的分位数上，因此

(6)

其中D是惩罚基函数系数的惩罚矩阵。其中对于第p个自由度的截断幂基

(7)

这意味着（2）中的多项式的系数是未被惩罚的的，而截断线性基函数的剩余系数被惩罚。

是所得到的P样条估计的粗糙度的度量，并且平滑参数lambda;可以使用如CV或GCV方法来选择。 下面部分探讨P样条平滑和线性混合效应模型之间的密切联系。

2.2平滑和混合效应模型

线性混合效应模型的一般形式

(8)

其中X和Z分别是固定效应和随机效应设计矩阵，beta;是固定效应的向量，u是随机效应的向量，其中和误差项。为了使拟合标准（6）与混合效应模型相关，平滑基因表达谱如（2）中所示。设X是由前（p 1）个基函数组成的设计矩阵，Z是由剩余的M个基函数组成的设计矩阵，

，

其中和。

惩罚平滑问题的最小化标准可以表示为

和是混合模型的BLUPs

, (9)

其中，，且和被假定为相互独立。为了方便起见，假设误差项ε是独立的，尽管在实践中更通用的结构为AR（1），AR（2）等，这些也可以通过最小化标准的适当修改来使用。然后将beta;和u的估计确定为混合模型方程的解（Henderson，1950）以及方差分量和的估计，因此平滑参数通过REML（参见Ruppert等，2003，第4.9节的进一步细节）确定。拟合值由下式给出

（10）

其中C = [X Z]，D如（7）中所述，，称为平滑矩阵。上述指定的模型可以通过R中的amer包来拟合，并且可以推广到模拟基因表达聚类。

在聚类环境中使用P样条的混合模型来表示有几个优点。这些包括在统一框架中建模聚类平均表达曲线和基因特异性表达曲线，处理缺失数据的优势，通过REML自动选择平滑参数，使用标准软件如S-plus，SAS，R等，以及将在第3节讨论的用于模型拟合和使用混合效应模型

剩余内容已隐藏，支付完成后下载完整资料

资料编号：[26625]，资料为PDF文档或Word文档，PDF文档可免费转换为Word