博客来源:https://www.vsni.co.uk/blogs/faster_GBLUP_ASReml-R我与ASReml 软件的故事我使用 ASReml 软件已经有10 年多了。从前基因组学时代开始,我在博士期间拟合了一个动物模型,其中系谱文件中有超过 50 万条记录,使用ASReml-R运行没有任何问题,在当时只花费了几分钟,就得到了一个包含所有育种值估计的结果文件。如今,计算机拥有更多的内存和计算能力,会更加缩短它的运行时间。随着基因组学时代的到来,我的研究方向也从基于谱系的动物模型转向基于基因组的动物模型。理论很简单:只需将您的分子关系矩阵A替换为基于基因组的关系矩阵G,您现在正在执行所谓的 GBLUP。然而,这个过程虽然看起来简单明了,但实际却并非如此。
我的第一个意外发现是生成G矩阵相关的复杂程度。首先,为了计算我的G矩阵,我必须预先处理我的分子数据(通常是 SNP)的技巧,即便是获取结果之后,仍然需要使用更多技巧(例如blending或bending)来调整我的G矩阵使其可逆且条件良好。然而,当我在使用 ASReml-R v. 4.1 的 GBLUP 分析中使用这个“调整过的” G矩阵时,最大的问题来了,通常,我的分析只考虑了两千到四千个基因型个体(一个比我旧的基于谱系的A矩阵小得多的矩阵,它有 50 万个体)。但令我惊讶的是,这个分析往往需要几个小时,或者因为 RAM 内存不足而无法工作。问题出在哪里?为了弄清楚这一点,让我们首先看一下为何旧的A矩阵运行的这么好。ASReml-R(以及 ASReml-SA)以稀疏形式生成A矩阵的逆矩阵。这意味着它被存储(和使用)为一个半对角矩阵,并且只保留非零值并用于计算。实际上,ASReml 使用稀疏矩阵的方法进行大多数计算。这不仅减少了存储所需的内存量,还节省了计算时间。需要注意的是,我过去拟合的A矩阵仅从三代生成,其中许多个体不相关,因此,我有一个非常空(或稀疏)的加性关系矩阵。这意味着即使我的矩阵非常大,但实际上它充满了零。估计只有 5-10% 的值是非零的。有限数量的非零数据使得分析非常有效。现在使用G矩阵,即便是只考虑 2,000 x 2,000 矩阵(并存储为半对角线),它也会充满非零值(可能全部为 400 万。
)。通常,它的倒数也充满了非零值。这是从基因组矩阵中预期的,因为任何和所有基因型都与所有其他基因型相关(即使只有非常低的相关性)。正如所料,这个矩阵比我旧的A矩阵密集的多,因此它的分析需要更多的内存,并且它无法使用 ASReml 中的稀疏矩阵方法,才导致我上面遇到的问题。如果我使用的是 R 环境,那么我会遇到其他问题。ASReml-R 可调用数据框对象(表型和矩阵)并生成潜在的大对象(例如解决方案)作为输出。这意味着 R 需要将所有分析前和分析后的对象保存在内存中。它将会占用我笔记本电脑的所有 RAM 内存,让我几乎没有空间进行统计分析。相比之下,ASReml-SA 读取表型数据并生成设计矩阵,但不保留原始数据。ASReml-SA 还会在生成文件时将输出写入文件(例如 BLUE 和 BLUP);因此,可以以更有效的方式使用内存。使 GBLUP 分析运行更快的 8 个技巧下面我将分享我使用或发现的进行模型拟合并获得结果的一些方法,它们使我的 GBLUP 分析得以成功运行。这不是一个详尽的技巧清单,但希望它能带给您一些启发。1) 使用替代软件确实是有一些软件包(和 R 包)对处理特定分析可能会更高效。但是,您必须确保考虑了复杂的线性混合模型。例如,某些例程不允许多个固定或随机效应,或者它们在误差结构规范中受到限制(例如,可能无法拟合异质或双变量误差)。就我而言,我很少选择它,因为我通常使用的数据类型需要高级线性混合建模。2) 选择GBLUP外的其他方法如果您的目标是拟合基因组预测模型,那么 GBLUP 是一个不错的选择。或者,您可以通过使用贝叶斯方法(例如,贝叶斯B)或多元回归算法进行拟合岭回归或 LASSO 模型。总的来讲,你需要依赖一组不同的统计假设。这个方法曾经对我有用,但用途有限,因为它只适用于相对简单的模型。3) 不使用R环境您可以在独立应用程序 ASReml-SA 中利用 ASReml 的强大功能和灵活性等优势功能。如上所述,该软件在处理问题的内存和空间方面更高效。与 ASReml-R 相比,符号和编码存在一些差异,但正如预期的那样,逻辑非常相似。我之前使用过此方法且效果不错,您可能还会注意到最新的 ASReml-SA 4.2 版比旧的 4.1 版本有几处增强的功能。4)让R更高效就我个人而言,我个人更喜欢 R 环境提供的灵活性和控制力,因此我经常不想脱离R环境。我发现了一些技巧。主要是从内存中删除不再需要的对象。例如,一旦您使用 ASReml-R 命令 ainverse()获得G逆矩阵,就不需要G矩阵。就可以消除任何中间矩阵或数据帧。此外,使用专门设计的库(例如 data.table)可以更轻松、更高效地处理大型数据集。然而,在过去几年中,R 的新版本改进了内存管理,我怀疑随着大型数据集成为常态,这种情况将继续下去。5) 转向云计算有时,我要分析的数据对于我的计算机来说仍然太大,或者由于时间紧迫而需要快速完成,这时,我就通过 AWS 使用了云计算(但 Microsoft Azure 是另一种选择),可以在惊人的机器配置上安装 RStudio 和 ASReml-R。如果您有这些服务的预算,这是一个不错的选择。但是您需要付出一些努力来建立这些系统。您可以点击下方链接查看获取有关如何实现最小化的一些好的建议。https://www.louisaslett.com/RStudio_AMI/ 6) 使用 ASReml-R 里的技巧ASReml-R(以及 ASReml-SA)中有几个功能选项,对拟合任何线性混合模型时都很有用。为了加速收敛,可能最重要的是提供一些好的初始方差分量——那些将接近最终估计值的分量。这在比较多个模型时特别有用。此外,您还可以将一些固定效应移到稀疏部分(通过使用 sparse),它可以减少计算负载。最后,如果您的模型太复杂(例如 MET 具有多个空间组件),那使用两阶段分析是一个不错的选择,第一步是分别拟合每个环境,第二步是对每个环境的预测(和权重)用于 GBLUP 分析,但此过程可能会丢失一些信息。7) 在分析中使用部分G矩阵在拟合线性混合模型时,唯一有助于估计方差分量或效应的基因型是那些具有表型响应的基因型。这意味着将G矩阵缩减为仅在表型数据中有记录的基因型可能会显著减小G矩阵的大小。尤其是在有很多历史基因分型与给定的目标分析分型无关的情况下,我会用这种方法。8) 对部分G矩阵进行条件预测当我仅使用具有表型观察结果的部分G矩阵来拟合线性混合模型(如上所述)时,我仍然需要对其他基因型进行基因组预测,这一直是我使用的主要“统计”策略。在这里,在您的模型拟合后,您可以使用估计的基因组育种值 (GEBV) 和完整的G矩阵来获得所谓的条件预测。本博客末尾针对“条件预测”提供了更多详细信息。每当您需要对未来基因分型个体进行预测时,这也是一个有用的策略,因为它不需要重新拟合原始基因组预测模型。如您所见,有多种方法可以提高您的分析速度,但“最佳”方法则取决于问题的类型。通常,我处理了少数具有基因型和表型的个体,使用部分G矩阵是一个很好的技巧;然而,我预计在接下来的几年里,我的G矩阵规模将大幅增长。尽管如此,在其他领域也取得了进展。当前VSNi 开发重点是使与 ASReml 核心算法相关的(ASReml 系列产品的基础计算引擎,包括 ASReml-R)这种类型的“密集”模型更容易、更快地适应。目前的工作包括:❆改进大型矩阵的内存管理;❆实现替代算法:❆混合模型方程的重新排序(减少填充,从而努力获取矩阵逆);❆几种矩阵运算(例如,用于反转某些特定矩阵类别的算法);❆在关键例程上增加并行化的使用(使用更多内核将加快进程);其中的一些功能已经在最新的 ASReml-SA 4.2 版中得以实现,我希望在接下来的ASReml-R 版本升级中,看到这些以功能及更多相关的计算改进。在不久的将来,我希望能在没有任何内存问题的情况下以合理的速度对数以万计的个体进行 GBLUP 分析。同样,在 R 环境下,我也期待在能这些方面取得进展。我的最后一个想法是有关未来的科学研究和发展。在 GBLUP 领域需要进行额外统计的、数值的甚至计算的研究。我确信,会有很好的“稀疏”近似算法来代替G矩阵,或者是一些巧妙的矩阵代数,可以更好地处理我们的混合模型方程。我们都应该关注最新发展,看看有没有新的可用的方法使这些分析在操作上更易于管理。关于作者Salvador Gezan 是一位统计学家/定量遗传学家,在育种、统计分析和遗传改良咨询领域拥有 20 多年的经验。他目前在英国 VSN International 担任统计顾问。Gezan 博士的职业生涯始于 Rothamsted Research(英国洛桑实验站),担任生物统计学家,在那里他使用 了Genstat 和 ASReml 统计软件。在过去的15年里,他为世界各地的公司和大学研究人员讲授 ASReml软件。希望本篇文章能给大家在GBLUP的运行速度方便带来一些启发,如果您需要使用ASReml R或ASReml SA,或者有其他的问题,可以在文末留言或私聊小编。
(图片来源网络,侵删)
0 评论