兼顾原核生物和病毒的宏基因组物种组成分析软件KMCP
文章题目：KMCP: accurate metagenomic profiling of both prokaryotic and viral populations by pseudo-mapping
发表杂志：Bioinformatics
影响因子：6.931
发表时间：2022年12月29日
第一作者：沈伟
通讯作者：沈伟、胡鹏、任红
通讯单位：重庆医科大学附属第二医院病毒性肝炎研究所
原文链接：https://doi.org/10.1093/bioinformatics/btac845
软件链接：https://github.com/shenwei356/kmcp/
摘要
日益增加的微生物参考基因组促进了宏基因组物种组成分析准确性的提高，同时也对宏基因组物种组成分析软件的建库效率、数据库大小、分析时间带来了严峻的挑战。此外，大多数软件主要关注细菌、古菌和真菌，而缺少对病毒群落的关注。本研究中，我们推出KMCP，一个新的基于k-mer的宏基因组物种组成分析软件。它将参考基因组分割成多个相同大小的区块，并将每个区块的k-mer存储在修改和改进的COBS索引中用于序列的快速搜索。KMCP通过结合k-mer相似性和基因组覆盖信息来减少基于k-mer的宏基因组分类和物种组成分析软件的高假阳性。基于模拟和真实数据的基准测试显示，KMCP虽然比其它软件运行时间更长，但是KMCP不仅能在原核生物和病毒群落上准确地分析物种组成，还能在低深度临床样本中提供高可信度的病原微生物检测结果。
算法简介
KMCP先将每个微生物基因组分割成n（默认10）等份，每份称为基因组区块（genome chunk）。每个基因组区块的所有k-mer用一个布隆过滤器（Bloom filter）保存，所有Bloom filter保存在一个改进的COBS索引中。在reads搜索的时候，将查询序列的所有k-mer去和数据库中所有的Bloom filter比较，返回具有一定共有k-mer数量的基因组区块。最终，我们可以知道查询序列与哪些基因组的哪个基因组区块有较高相似度，这个过程我们称为pseudo-mapping。搜索完成后，匹配的参考基因组经过三轮过滤筛选（过滤条件包括k-mer相似性信息与基因组覆盖信息），最后再使用Expectation-Maximization（EM）算法估算匹配的参考基因组的相对丰度。
软件安装与使用
KMCP由Go语言实现，提供Linux、Windows、macOS操作系统的amd64和arm64架构的静态链接可执行程序，用户只需要下载解压，开箱即用。也可通过conda安装：
conda install -c bioconda kmcp
KMCP为命令行工具，采取“主命令-子命令”的架构，包含k-mer计算、数据库构建、Reads搜索、搜索结果合并、物种组成分析这5个主命令和一系列辅助命令。KMCP的5个主要子命令如下：
-
- Compute 从FASTA/Q文件计算k-mer，
-
- Index 从k-mer文件构建数据库，
-
- Search 从一个KMCP数据库中搜索序列，
-
- Merge 合并多个数据库的搜索结果，
-
- Profile 从搜索结果生成分类学组成谱。
我们在https://bioinf.shenwei.me/kmcp上提供了详细的文档，包括子命令用法、教程、常见问题等。
分类学数据
与其它方法不同，我们并没有将分类学数据信息完全集成到数据库中，而是采用相对独立的设计，仅仅需要NCBI格式的分类学数据文件（taxdump files）。这样做的好处是能在不重新建库的情况下使用更新的分类学数据版本，而且对于一些第三方或者自定义的分类学数据库比如GTDB Taxonomy或者ICTV Taxonomy，可以通过我们前期开发的Taxonkit工具中的create-taxdump命令从任意谱系数据创建分类学数据文件。甚至当参考基因组由来自GTDB的原核生物和来自NCBI的真菌和病毒组成时，还可以合并GTDB和NCBI分类学数据。
可扩展的搜索模式
基于Lowest Common Ancestor (LCA)的方法获得查询序列每个k-mer对应的分类学ID即TaxId，再计算LCA并给查询序列作标记。而在KMCP中，查询序列的所有k-mer去和数据库中所有的Bloom filter比较，返回具有一定共有k-mer数量的基因组区块。KMCP搜索过程不对查询序列赋予任何的分类学标签，因此，同样的查询序列在不同数据库中搜索的搜索结果可以合并在一起。搜索结果的可合并性（mergeability）允许用户从不同的参考基因组集合创建不同的数据库，并根据需要选择需要的数据库来搜索。另外，由于KMCP搜索速度与基因组区块的数量线性相关，而且所有基因组区块的索引数据（即Bloom filter）相互独立，用户可以将参考基因组分为多个较小的部分，每个分布分别创建数据库并用来搜索。这样，用户可以利用计算机集群进行加速搜索，每个计算节点负责一个小数据库的搜索。而内存容量有限的个人计算机也可以创建许多较小的数据库，从而利用大规模的参考基因组。
预建数据库
虽然KMCP建库效率高且简单方便，我们仍为用户提供了3个预建的数据库。对于细菌和古菌，我们采用GTDB r202的47,894条基因组。对于真菌，我们从Refseq r208下载了属于398个物种的403条基因组。对于病毒序列，我们从Genbank r246下载所有病毒序列，每个TaxId最多保留5个基因组，最终有属于23,632个物种的27,936条基因组用于构建病毒基因组数据库。分类学数据使用的是2021年12月6日下载的NCBI Taxonomy数据库的taxdump文件。
准确性评估
为了客观、综合、全面评估KMCP的宏基因组物种组成分析准确性，我们采用4个有金标准的公共数据集来进行基准测试，包括：CAMI挑战中的CAMI海洋数据集、模拟的原核数据集、真实的病毒组数据集、临床感染体液数据集。CAMI海洋数据集物种种类最多，且丰度跨度最大；模拟的原核数据集物种相对较少，丰度跨度适中；真实的病毒组数据物种数量较少，但物种序列较为相似；而临床感染体液数据主要特点是测序深度低。
① CAMI挑战中的CAMI海洋数据集
在属水平，KMCP在分类学识别（taxon identification）和相对丰度（relative abundance）估算上的准确性与mOTUs2和MetaPhlAn2类似。而在种水平，KMCP优于所有其它工具：KMCP的Completeness（Recall,召回率，0.915）比Bracken（0.944）略低，比mOTUs2（0.809）、MetaPhlAn2（0.850）和Centrifuge（0.887）高；而KMCP的Purity（Precision，精确度，0.830）比mOTUs2（0.887）和DUDes（0.871）低。在相对丰度上，KMCP的L1 norm error（越低越好）在属水平排在第二低，而在种水平最低，优于其它软件。另一个在所有分类学水平衡量软件预测丰度与金标准的差异的指标weighted UniFrac error，KMCP的值最低，优于其它所有软件。
图3. CAMI海洋数据集上的宏基因组物种组成分析准确性。（a）在属水平和种水平的分类学识别准确性指标：completeness/recall、purity/precision、and F1 score，以及丰度估计准确性指标：L1 norm error。（b）所有分类层级的相对丰度准确性指标。（c）以不同阈值（累计0到10%）过滤后的种水平分类学识别准确性指标。
② 模拟的原核数据集
这里我们用所有软件的最新版本，在另一个模拟的原核基因组数据进行测试。在属和种的水平，基于k-mer的方法（KMCP，Bracken，Centrifuge和Ganon）和基于标记基因的方法（MetaPhlAn3和mOTUs3）均展现出较高的Completeness，平均值均高于0.9；然而基于比对的方法包括DUDes和SLIMM效果稍弱。至于Purity，在属水平，KMCP、MetaPhlAn3和mOTUs3的值最高，均高于0.9，紧随其后的是Centrifuge（0.740）；在种水平，mOTUs3的Purity仍在0.9以上，而MetaPhlAn3和KMCP则分别降到了0.855和0.803。至于F1 score，在属水平KMCP比MetaPhlAn3和Centrifuge高，而在种水平比MetaPhlAn3稍低。至于L1 norm error，KMCP在属水平表现最好，紧接着是mOTUs和MetaPhlAn3；而在种水平，mOTUs3表现最好，紧接着是MetaPhlAn3和mOTUs。Weighted UniFrac error的趋势和属水平L1 norm error的趋势相同：KMCP在属水平表现最好，紧接着是mOTUs3和MetaPhlAn3。
图4. 模拟细菌数据集上的宏基因组物种组成分析准确性。（a）在属水平和种水平的分类学识别准确性指标：completeness/recall、purity/precision、and F1 score，以及丰度估计准确性指标：L1 norm error。（b）所有分类层级的相对丰度准确性指标。（c）以不同阈值（累计0到10%）过滤后的种水平分类学识别准确性指标。
③ 真实的病毒组数据集
由于该噬菌体群落是通过混合实验室培养的病毒颗粒而来，则测序数据中可能含有未消化完全的细菌DNA。然而图5a显示，在界（superkingdom）水平，各软件的组成结果差别较大。KMCP和MetaPhlAn3显示了较高的病毒丰度比例，分别为99.72%和99.02%；而其它软件的比例则在80%以下。然后，我们用TaxonKit移除非病毒预测后重新计算各基因组的相对丰度。如图5b所示，在属和种水平，使用了和KMCP相同的参考基因组的Bracken、Centrifuge、和Ganon与KMCP的Completeness均在0.9以上；使用官方数据库的Bracken和MetaPhlAn3则低很多。至于Purity，在属水平，KMCP和Centrifuge的值均在0.95以上,不过在种水平，它们分别降到0.876和0.691；而Bracken (KMCP)则从0.857降到了0.369。从种水平的F1 score来看，KMCP表现最好（0.891）,Centrifuge次之（0.814）。同时，KMCP的丰度估计准确性在属水平和种水平都表现最好，L1 norm error和weighted UniFrac error均最低（图5b和c）。
图5. 真实的病毒组数据集上的宏基因组物种组成分析准确性。（a）在属水平和种水平的分类学识别准确性指标：completeness/recall、purity/precision、and F1 score，以及丰度估计准确性指标：L1 norm error。（b）所有分类层级的相对丰度准确性指标。（c）以不同阈值（累计0到10%）过滤后的种水平分类学识别准确性指标。软件名称以（KMCP）结尾的表示该软件使用和KMCP同样的参考基因组创建数据库。
④ 临床感染体液数据集
因为采用了更新的数据库，Bracken的敏感性和特异性均得到较大提高。采用和KMCP相同参考基因组集合的Bracken，Centrifuge和Ganon的敏感性均为80.26%；KMCP，Bracken和Ganon的特异性均为100%。KMCP由于在S19样品上识别错误，遗漏了一个病原菌。因而敏感度稍低于Bracken和Centrifuge。而Centrifuge较低的原因是，它在3个样品（S33、S35、S82）中错误地检测出了一个真菌（Coccidioides immitis）。
另一方面，Bracken，Centrifuge和Ganon在分析结果中主要报告预测的基因组/物种匹配的reads数量。而KMCP则能报告更多的指标，包括基因组区块覆盖比例，相似度得分（所有唯一匹配reads的k-mer覆盖度的90百分位值），这两个指标的值也用于预测基因组的排序。这些额外的指标使得KMCP能够提供更可靠的结果：更少的预测结果（图6a），且病原微生物在预测结果中排名更靠前（图6b）。对于仅有少量reads匹配的预测结果，KMCP预测结果中的基因组区块覆盖比例也能提供有价值的证据来帮助用户判断病原菌存在与否。如图6c所示，4-10条reads也能达到0.4-0.5的基因组区块覆盖比例。
图6. （a）所有软件都准确检测的55个样品中的预测结果数量，样品按照Bracken (KMCP)的预测结果数量排序。（b）所有软件都准确检测的55个样品中，病原微生物在预测结果中的排名。（c）KMCP正确检测的60个样品中，reads数量和基因组区块覆盖度的关系。软件名称以（KMCP）结尾的表示该软件使用和KMCP同样的参考基因组创建数据库。
计算性能评估
综上所述，本文基于三代ONT序列提高了宏基因组装的质量、SV的发现能力，发现了大量包括插入突变和基因倒位在内的结构变异对于菌株水平上基因功能的影响，以及噬菌体、CRISPR-spacer等系统的的深度挖掘。这项研究是课题组利用三代Nanopore测序技术解析肠道病毒组（Cao et al., Medicine in Microecology, 2020, 4:100012; Cao et al. Gut Microbes, 2021, 13），近期发表的真菌组分析方法（Lu et al, Molecular Ecology, doi:10.1111/mec.16534），和靶向RNA检测病原微生物（Zhao et al., Advanced Science, 2021, 8, 2102593）之后，在利用三代Nanopore测序技术探索肠道微生物研究领域的新进展。
图 7. 宏基因组分析的用时和最大内存占用。KMCP的构建时间包括compute和index两个命令。软件名称以（KMCP）结尾的表示该软件使用和KMCP同样的参考基因组创建数据库。KMCP_HPC是KMCP在集群中运行的场景，每个计算节点搜索一部分参考基因组构建的数据库。
总结
我们开发了一款新的宏基因组物种组成分析软件KMCP，它使用改进的COBS数据结构来对微生物参考基因组进行索引，并结合k-mer相似性与基因组覆盖信息来提高宏基因组物种组成分析的准确性。KMCP数据库构建的高效性和序列搜索的可扩展性使得KMCP能有效的利用日益增多的参考基因组。在模拟和真实数据集上的基准测试表明：KMCP虽然比其它软件运行时间更长，但是KMCP不仅能在原核生物和病毒群落上准确地分析物种组成，还能在低深度临床样本中提供高可信度的病原微生物检测结果。
参考
