微生物文献解读：用问题形式串联文章主要思路及其结论，快速看懂正文主图。

MetaCHIP: community-level horizontal gene transfer identification through the combination of best-match and phylogenetic approaches

• 杂志：Microbiome [11.607]
• 发表时间：4 March 2019
• 第一单位：Centre for Marine Bio-Innovation, University of New South Wales
• 第一作者：Weizhi Song
• 通讯作者：Torsten Thomas
• 链接：https://microbiomejournal.biomedcentral.com/articles/10.1186/s40168-019-0649-y

简介

本文设计与开发算法MetaCHIP，结合序列比对和系统发育分析，可以实现无参考基因组的、在菌群层次上横向转移基因（HGT，horizontal gene transfer）的检测。

背景

• Q: 什么是横向转移基因（HGT）？
• A: 是指有机体之间遗传物质的转移，通常是微生物进化和适应性的重要原因（例如抗性基因，毒力基因）。
• Q: HGT和MGE（可移动遗传元件）的区别是什么？
• A: 虽然本文没有提及，但是我的理解是HGT可以是单个gene的转移。MGE是元件（element），指得是可以转移的（一般为）一簇（cluster）基因，他们一般包含行驶转移和整合功能的gene与所携带的gene（cargo gene），后者一般就是HGT关注的对象。
• Q: 为什么要研究bacterial community的HGT？
• A: 细菌本来就是以community发生作用，随着宏基因组学数据积累，可以让我们研究bacteria之间gene的转移，从而让我们从一个菌落的角度理解他们发挥的功能。
• Q: 以前都有什么方法？MetaCHIP的优势在哪里？
• A: 前人方法基本分为3种，以1.基因组成分特征（compositional features）为主的方法（如GIST，IslandViewer）；2.最佳比对（best-match）为主的方法（如DrakHorse，HGTector）；3. 显性系统发育分析为主的方法（如Ranger-DTL，AnGST），即对比gene tree和species tree之间的不一致性从而找到转移的gene。
  以往方法一是不适用于community，二是需要reference genome，这都限制了它的运用。MetaCHIP就是能解决以上两点。

方法

Figure.1 MetaCHIP的workflow.

算法主要分为两部分：best-match + phylogeny.

• best-match：

1. 将整个基因组数据中的ORF预测出来
2. 将ORFs两两用BLASTN比对，经过参数筛选后，得到两两的matches。
3. 将所有genomes分成group（按照分类学注释分组，例如A，B，C，这里的group可以是class/order/genus层面），各组的一条genome分别叫Ax，By，Cz; 每条genome上的一个gene叫Ax-N，By-N，Cz-N; 假设A组中有个gene叫A1-01，它的BLASTN结果在每组match为IAx，IBy，ICz。它与每组match的平均identity为IAA（不算自己），IAB，IAC。
4. 如果IAA不是所有结果中最大的，并且IAA不等于0（即除了自己在本组都没相似的gene，此时这个gene就可能是黑户混进来的），那么除了A之外的组就是这个gene为HGT的putative candidate group。然后这组（group）中BLASTN identity最高的gene就是putative HGT candidate。
5. 计算gene A1-01这组gene和HGT putative candidate group中gene的所有identities，利用百分数确定一个cutoff（例如我们把前10%分位数对应的identity确定为cutoff），当gene A1-01与putative HGT的identity > cutoff, 这个HGT才确定下来。

• phylogeny：

1. 通过best-match确定下来的gene pair，将这两个gene来自的两个groups中的所有orthologs拿出来建一棵gen/protein tree
2. 因为16S在genome bins时常常缺失，因此对于所有的genomes，先用CheckM识别出的43个通用的单拷贝gene（SCGs），，每个gene用HMMER进行多序列比对（MSA）后，串联构建所有genomes的species tree。然后，从这棵树中抽取只和HGT有关genes所在的genomes的subtree。
3. Ranger-DTL计算gene/protein tree与species subtree的一致性，从而判断是否发生HGT，并确定gene flow direction。

补充：

Ranger-DTL的大致算法：比较gene tree和species/organism tree的叶节点，看gene tree要怎么样经过duplication、transfer、loss事件的组合可以变成species tree，选择最简约的结果输出。

Q: 用SCG来建树真的可以反映species tree吗？
A: 为了看在多大completeness取值下，两种树可以相互替换，本文选择2个class下的genomes共 20条，将他们全genome划分成100个等长contigs，对于每个contigs依次随机抽出20，40，60，80个contigs，作为completeness取值范围从20%-80%的数据，用Mantel test（比较两个矩阵相似性的方法）比较SCG tree与species tree（基于16S rRNA）的相似性。

Figure.3 不同completeness下SCG tree与species tree的相似性。可以看到在completeness>40% 两种树一致性较高（>72.64%）。

结果

Q: 如何评价MetaCHIP模型的表现?
A: 本文使用了3种数据来进行模型测试（第3种数据更像是模型应用与数据挖掘：

1. 模拟数据 [genomes]：在genus level，选择10条来自Sphingobium (作为donor group) 与10条来自
   Sphingomonas (作为recipient group)的genomes；在class level，选择10条来自alphaproteobacterial
   (作为donor group) 与10条来自betaproteobacterial (作为recipient
   group) 的genomes（Table S1）。对每个donor genome选择10个gene随机转移到10个recipient genome去（donor genome的原gene自然要保留）。
2. 已验证数据 [genomes]：来自人体不同部位以及环境中2094完整bacterial genomes，这套数据前人分析过，在16S rRNA相似性>97%的genomes中把>500bp, identity>99%的DNA片段认为是HGTs。本文使用MetaCHIP在genus level进行分析。并将BLASTN identity与coverage=100%作为cutoff，将MetaCHIP结果和前人结果进行比较。
3. 真实未验证数据 [metagenomes]：一套human gut genome bins数据（去除contigs<2000bp）以及一套North Sea海水 genome bins（QC之后用metaSPAdes进行assemble）作为真实数据。

---

MetaCHIP在模拟数据上的表现

Q1: 有什么因素可能影响HGT的预测？在模拟数据上可以增加什么扰动？
A1: 目前本文考虑到了三个主要因素：recipient genome的突变程度（genetic diversity）、测序深度（depth）、测序覆盖度（coverage）以及genome assemble的软件。
因此 1）genetic dievrsity：本文对recipient genome进行0%，5%，…，30% genetic diversity的基因突变（使用HgtSIM完成），并且每个genetic diversity下bootstrap 10次。
2）depth & coverage：对于class level的数据，不同genetic divversity下从genome打撒生成不同丰度的reads，满足depth range: 3，6，9，12 milion reads,它们对应的average coverage range为 6，11，17，23X.
3）assemble software：IDBA_UD，metaSPAdes.

Q2: 上面说到将genome打撒成reads，但是可能在assemble中transfer的gene就拼不回去了，那什么情况下gene transfer才算重建了？
A2: 在assemble过程中，如果gene的任一边flanking region > 1 kbp 且与recipient genome能match，就算是在assemble中重建了这个gene。

Q3：怎么评价MetaCHIP在模拟数据集上的表现？
A3：从Fig4-5来详细看。
Figure.4 MetaCHIP在class level与genus level上（不同颜色）的模拟数据集上，在不同genetic diversity（横轴）情况下，对HGT复原比例（纵轴）。

从Fig.4 可以看到在genetic diversity >=10%时，recovery percentage在class和genus level都很好。class level在genetic diversity=30% 抢救无效。genus level在20%就无效了。总体上，genus level表现不如class表现好，也和前人研究结论一致。总而言之，当genetic diversity<=15%时，MetaCHIP可以找到>=40%的HGTs（我的问题：在这个数据集上为什么不探讨一下假阳性？）

Figure.5 MetaCHIP在不同测序深度（对应相应的coverage，横轴）下，不同genetic diversity以及assemble软件（不同颜色）的HGT复原比例（纵轴）

从Fig.5 可以看到，当没有突变时用metaSPAdes好，其他时候选用IDBA_UD（我的问题：当然两者在5%都很差，为什么？）。而且在无突变的6M reads时metaSPAdes效果最好。既然效果不是随着测序深度线性增加，所以根据复原比例看17X（9M）是最好的coverage（depth）。这就定下了什么genetic diversity用什么assembler，本文对17X的reads assemble
后的contigs构建genome bins。

Q4: contigs和bin有什么区别？
A4: 指从微生物群体序列中将不同个体的序列（reads或contigs等）分离开来的过程。简单来说就是把宏基因组数据中来自同一菌株的序列聚到一起，得到一个菌株的基因组。是的，可以达到菌株水平。这个过程就是binning，也就是说binning可以是contigs的分组。在本文是对genomes进行binning。

Q5: 那bin后MetaCHIP的效果如何？
A5: 看Fig.6-7。

Figure.6 assembly与binning后在不同assembled/binned genomes上不同genetic diversity（横轴）下的MetaCHIP HGT复原（纵轴）表现。

可以看到bin的结果中，genetic diversity=0%时，30%的gene在donor中找到，但是这30%都在recipient中找到了。genetic divesrity>5%时，不少于73.7% genes在donor和recipient都找到了。

Figure.7 MetaCHIP 对模拟数据（黑色）真实数据（灰色，即Fig4 class level那条线）在不同genetic diversity（横轴）上HGT复原（纵轴）表现。

可以看到，模拟数据在genetic diversity=0%时，26%的100个genes被best-match还原，其中9个（9%）被phylogeny证实，而在真实数据中这两个值是86.7%和30%；genetic diversity=5%时，模拟数据BM（best-match）和PG（phylogeney）为93.8%和81.3%；genetic diversity=10%时，MetaCHIP在模拟数据表现最好，BM和PG两个值分别是74%和69%，真实数据的指标为91.4%和85.2%。

Q6: 为什么看Fig.6复原率要分成donor和recipient？直接看recipient不行吗？
A6: 因为transfered gene是在donor genomes中的，assemble和bin可能使得HGT gene丢失，可能是通过对比看donor中的gene才能知道recipient gene中的gene有多少是可以assemble/bin出来的。

MetaCHIP在已验证数据上的表现

Q1: MetaCHI与前人验证结果相比如何？
A1: 见Fig8-9. Fig.8展示了在数量上HGT结果。Fig.9 在gene功能注释层面（COG）的HGT结果。
Figure.8 MetaCHIP从368 soil genomes上发现的HGTs（只做soil genomes是因为原数据metadata不完整）

因为原数据将genetic diversity<=1% 定义为HGT，即在本文中定义的recent HGT（定义为genetic diversity>1%的HGT），所以对应Fig.8中第一个bar，为433个recent HGT。原数据找到了368个HGTs，两者重合数目为248（占MetaCHIP的77.2%）【注意：这里的genetic diversity是在gene层面不是上部分的genomee层面】

Figure.9 MetaCHIP在368 soil genomes上找到的recent HGT和non-recent HGT的COG注释结果与所有genome COG比例的比较。

MetaCHI不仅找到了recent-HGT还找到了non-recent HGT。对所有2094条input genomes上的HGT进行COG注释，得到boxplot展示的分布，当HGT的COG比例大于背景COG的75%（也就是高于boxplot箱子最上横线–第三百分位数）那么就认为这个COG在HGT中富集。结果发现：non-recent HGT在COG的C（能量产生与转化）、E（氨基酸转运与代谢）、I（脂质转运与代谢）和L（复制、重组和修复）中富集，recent-HGT在K（转录）、P（无机离子转运与代谢）和U（细胞内交易、分泌和小泡转运）中富集。

---

MetaCHIP在真实数据上的表现

Human gut数据原数据即为1634 genome bins，其中completeness>40%，无污染的bin有138个。
North Sea数据为reads，assemble后得到contigs再binning后得到69 bins，其中completeness>40%，无污染的bin有37个。
Q1: MetaCHIP在真实数据中找到了多少HGT？结果如何？
A1: 在human gut数据上BM找到了560个HGT，PG找到了113个； 在North Sea数据上BM和PG找到了121和32个。见Fig10-11。

Figure.10 Human gut和North Sea数据中HGT的gene flow circos图。

genome bins数目越大的genome，HGT就越多。

Figure.11 对2套数据HGT做COG注释的boxplot结果。

human gut发现HGT主要富集在defence mechanisms
(V); energy production and conversion ©; translation, ribosomal structure and biogenesis (J); and nucleotide transport and metabolism (F)。其中防御系统V的富集主要是15个与ABC-type多药以及抗菌肽转移系统相关。说明human gut的AR gene是被转移的主要对象。
North Sea主要富集在energy production and conversion ©, carbohydrate transport and metabolism (G) and translation, ribosomal structure and biogenesis (J)。
C和G的注释和前人关于free-living微生物的研究结果一致。

Q2: MetaCHIP还有什么改进之处？
A2: 真实数据MetaCHIP找到的所有HGTs中40% genetic diversity>25%（见Table 2）, 但是在模拟数据中BM和PG都只检测到了<~20%（见Fig.7）的HGTs，说明HGTs数目还可以提高0.4x/0.2*1x=2倍。因此，对genetic diversity 较大的HGT检测还有待改进。

---

本文可改进处：

1. 对突变较大的数据低估了HGT的数目
2. 未与其他方法在同一数据集上做比较
3. 可不可以加一些参数让算法各自适用于human gut或者free-living的数据
4. 对于了解较为全面的human gut菌群是不是可以适当增加一些reference genome dataset增加结果的可信度
5. 未对算法时间和内存做出介绍，目测assemble将会是算法时间的较大瓶颈

---

总结

本文构建算法可以检测菌群中的HGT，摆脱了对reference genome的依赖，可以有效应用于快速积累的各种菌群数据。并在三种类型数据集上进行了算法表现的测试。算法思路虽然不复杂，但是对结果的探讨比较详细。最可惜的是并未与其他算法在性能上进行比较。