SNP数据相关性和主成分分析的高性能计算工具集

摘要

总结：全基因组关联研究被广泛用于研究疾病和性状的遗传基础，但它们带来了许多计算挑战。我们开发了gdsfmt和SNPRelate（用于多核对称多处理计算机体系结构的R包）来加速SNP数据的两个关键计算：主成分分析（PCA）和使用逐降一致性度量的相关性分析。我们的算法核心是用C/C++编写的，并且经过了高度优化。基准测试表明，PCA和identity-by-descent的单处理器实现分别比流行的EIGENSTRAT（v3.0）和PLINK（v1.07）程序中提供的实现快约8–50倍，并且通过使用八个核可以将速度提高30–300倍。SNPRelate可以分析成千上万个具有数百万SNP的样本。例如，我们的软件包用于对“基因-环境协会研究”联盟研究中的55 324名受试者进行PCA。

可用性和实施：可从R CRAN获得gdsfmt和SNPRelate(http://cran.r-project.org)包括一个小插曲。可以在以下位置找到教程：https://www.genevastudy.org/Complishments/software.

联系人： zhengx@u.washington.edu

1简介

全基因组关联研究（GWAS）广泛用于研究许多复杂疾病和性状的遗传基础，但基于芯片和序列的GWAS从数千个研究样本和数百万SNP中生成的大量数据带来了巨大的分析和计算挑战。一个重要的挑战是，当存在人口结构和隐性关联时，GWAS结果中出现的虚假阳性关联被夸大了(Cardon和Palmer，2003年；崔等。, 2009). 这些挑战可以通过使用主成分分析（PCA）来检测和纠正人口结构来解决(价格等。, 2006)以及确定每对研究样本之间关联程度的同一性逐次下降（IBD）方法。对于这两种方法，建议使用一组剪枝后的SNP，这些SNP彼此之间近似处于连锁平衡，以避免SNP簇的强烈影响(劳里等。, 2010). 然而，与这些方法相关的计算负担在大样本和SNP大小的情况下尤其明显，并且需要有效的数字实现和内存管理，特别是当芯片阵列大小增加且序列数据用于调用变量时。例如，1000基因组项目最近从全基因组测序技术中确定了约1500万个SNP位点(1000基因组项目联合会，2010年).

R是最流行的统计编程环境之一，但它通常没有针对高性能或并行计算进行优化，这将减轻大规模基于SNP的GWAS计算的负担。为了克服这些限制，我们启动了一个名为CoreArray的项目(http://corearray.sourceforge.net/)它包括两个R包：gdsfmt，用于为全基因组数字数据提供高效、平台依赖的内存和文件管理；SNPRelate，用于解决多核对称多处理计算机体系结构上的大规模、数值密集型GWAS计算（即PCA和IBD）。

2特点

为了支持全基因组数字数据的高效内存管理，gdsfmt为面向阵列的数据提供了基因组数据结构（GDS）文件格式。在这种格式中，每个字节最多编码四种SNP基因型，从而减少文件大小和访问时间。在扫描SNP图谱的过程中，可以同时对四种基因型进行操作。GDS格式支持数据阻塞，因此只有正在处理的数据子集才需要驻留在内存中，而且它还被设计用于对大型数据集的有效随机访问。虽然SNPRelate函数仅在GDS格式的数据文件上运行，但我们的包提供了从PLINK、排序变量调用格式、NetCDF和其他数据文件重新格式化数据的函数(达内塞克等。, 2011；劳里等。, 2010；Purcell公司等。, 2007).

SNPRelate为GDS基因型文件的PCA和IBD相关性分析提供了计算效率高的函数。遗传协方差矩阵和成对IBD系数的计算被划分为非重叠部分，并分配给多个核以加速性能，如所示图1PCA的SNP相关函数包括样本和SNP特征向量的基本计算，以及有用的辅助函数。样本特征向量和观察到的等位基因剂量之间的相关性可用于评估SNP对每个特征向量的影响的全基因组分布。SNP特征向量可用于计算新样本集的样本特征向量，这在具有实质相关性的研究中很有用(朱等。, 2008).

对于相关性分析，SNPRelate中的IBD估计可以通过矩量法（MoM）进行(Purcell公司等。, 2007)或最大似然估计(崔等。, 2009；米利根，2003年)通过国家身份。我们的经验表明，对于大规模数据分析，MLE比MoM计算量大得多，尽管MLE估计通常比MoM更可靠。此外，连锁不平衡剪枝函数生成了一个彼此近似连锁平衡的SNP剪枝子集，以避免SNP聚类在PCA和IBD分析中的强烈影响。任何一种方法都可以估计出个体的实际亲属关系矩阵，可用于下游关联分析(价格等。, 2010).

这两个R包都是用C/C++编写的，使用POSIX线程库在类Unix系统上进行共享内存并行计算，并且有一个R接口，在该接口中，内核通过阻塞计算来利用高速缓存进行了高度优化。对算法进行了优化，以逐块加载基因型，而不限制SNP的数量。这些算法仅受并行线程访问的主内存大小的限制，并且存储着遗传协方差矩阵或IBD系数矩阵。

GDS也被R/Bioconductor包GWASTools用作其数据存储格式之一(戈加滕等。, 2012). GWASTools为GWAS的质量控制和分析提供了许多功能，包括SNP或扫描统计、批量质量、染色体异常、关联测试等。

3性能

我们使用小型、中型和大型测试数据集来说明SNPRelate的性能。小样本集和中样本集是根据模拟数据构建的，分别包含500和5000个样本和100个K SNP标记。这一大组由55 324个课题组成，这些课题选自“基因-环境协会研究”（日内瓦）联合会的16个项目(科内利斯等。, 2010). 我们分别比较了SNPRelate与EIGENSTRAT（v3.0）和PLINK（v1.07）的PCA和IBD估计运行时间。这些实现是在一个系统上进行测试的，该系统具有两个四核英特尔处理器，运行2.27 GHz和32 GB RAM，并运行Linux Fedora 10。

如所示表1SNPRelate中PCA和IBD的单处理器实现分别比EIGENSTRAT和PLINK中提供的实现快约8–50倍。当使用八个核运行SNPRelate算法时，性能改善范围为～30到300。对大型数据集进行SNPRelate PCA(n个=55 324名受试者，具有～310 K的选择性SNP标记）。当使用八个核心时，计算遗传协方差矩阵（55 K×55 K）需要～64小时，使用单处理器版本的LAPACK in R计算特征值和特征向量需要～9天。对中小规模数据集的分析需要<1 GB的内存，对～55 K受试者的主成分分析需要～32 GB，因为遗传协方差矩阵存储在线程共享的主内存中。PCA运行时间的一个改进是使用多线程版本的BLAS来执行特征值和特征向量的计算，而不是默认的单处理器版本。尽管SNPRelate比EIGENSTRAT（用于PCA）或PLINK（用于使用MoM估计IBD）快得多，但结果在数值上是相同的（即相同的精度）。

致谢

作者感谢日内瓦财团成员(http://www.genevastudy.org)用于访问用于测试gdsfmt和SNPRelate包的数据。

基金：美国国立卫生研究院，基因，环境和健康倡议。遗传协调中心（U01 HG 004446）。

利益冲突：未声明。

工具书类

作者注释