将基因组组装到染色体水平无非就是两种方式:
- 独立组装(de novo);
- 基于参考基因组的组装(reference-guided);
** 即用近缘或者同一物种的基因组染色体同源性进行组装。
独立组装中的contigs或scaffolds一般用遗传图谱、BioNano光学图谱或者HiC技术进行染色体水平锚定。而基于参考基因组的组装一般是将contigs或scaffolds与近缘或同一物种的基因组进行比对,根据同源相似性,将长度提升至染色体水平。有的时候,我们并没有相关图谱数据,只能采用后者将contig或scaffold提升至染色体水平,使用工具--RagTag。
官网:
malonge/RagTag: Tools for fast and flexible genome assembly scaffolding and improvement (github.com)
RagTag简介:
RagTag流程大体分为四个部分:correct--scaffold--patch--merge
1. correct:
矫正是使用参考基因组来鉴定和矫正contigs中的潜在的组装错误,该步骤不会擅自凭空增加或减少序列,仅仅是将序列在错误的组装位置进行打断。
ragtag.py correct reference.fa query.fasta -t 80
重要的结果文件:
- ragtag.correct.fasta:是矫正后的序列文件;
- ragtag.correct.agp:序列打断位点坐标,AGP文件格式
说明:若reference和待延长的scaffold是不同的基因型(或单倍型),则会引起错配(Reference-guided misassembly signatures),这是由真实的生物学结构变异引起的。实际操作中,最终由用户确定是否对错误组装进行更正。作者也建议:如果可行的话,应该用独立的数据(例如物理、管光学和遗传图谱)验证所有的RagTag结果。
2.scaffold:
将相邻的contigs序列用100个N连起来,序列的位置和方向需要根据参考基因组的比对来确定
ragtag.py scaffold reference.fa ragtag_output/query.corrected.fasta -t 80 -C
#-C会将用到的contig/scaffold连在一起放在chr0中(中间用100个N连接)
重要的结果文件:
- ragtag.scaffold.fasta
- ragtag.scaffold.agp:连接使用的contig/scaffold位置信息以及gap信息。
3.patch:
将上一步延长的scaffold序列进行gap填补,大致方法也是使用reference对上一步的结果进行fills(利用一开始的contig序列去填补上一步的assembly gaps)和jions(连接不同的contigs);实在是没办法补填的gap,最终用N。
ragtag.py patch ./ragtag_output/ragtag.scaffold.fasta query.fasta -t 80
重要的结果文件:
- ragtag.patch.fasta:The final FASTA file containing the patched assembly
- ragtag.patch.agp:The final AGP file defining how is builtragtag.patch.fasta
4.merge:
在scaffolding过程中,可能会根据不同的参数或者图谱数据产生多个不同版本的基因组组结果,该步骤将不同的结果根据权重进行最终结果的生成。
如果有HiC数据,还可以加入HiC数据生成比较好的组装结果;
#多参考基因组的情况:
ragtag.py scaffold -o out_1 ref1.fasta query.fasta
ragtag.py scaffold -o out_2 ref2.fasta query.fasta
ragtag.py merge query.fasta out_*/*.agp
#有Hi-C数据的情况:
ragtag.py merge -b hic.bam query.fasta out_*/*.agp
参考:
RagTag基于参考基因组的组装和GeMoMa注释 - 简书 (jianshu.com)
使用初探:
现有某单菌的全基因的二代测序数据,使用spades从头组装至contig和scaffold水平,后试想是否能用RagTag锚定NCBI上的同物种的基因组,从而延长contig/scaffold长度至全长基因组。
#spades.py组装:
spades.py --careful --pe1-1 Ery-1_1.clean.fq.gz --pe1-2 Ery-1_2.clean.fq.gz -o ./Ery
#使用Quast评估:
quast.py contigs.fasta -o quast_out
#使用prokka对组装contig进行注释:
prokka contig.fasta --outdir prokka_annotation --prefix ery --kingdom Bacteria --cpu 20
#从注释结果中找到16S基因序列,进行NCBI blast rRNA数据库,并下载近缘种的全基因组序列
#用作下一步RagTag延长的reference;这里我16S比对链球菌(Reference),覆盖率和相似性均99.999%。
#RagTag使用:
ragtag.py correct reference.fa contig.fasta
ragtag.py scaffold reference.fa ragtag_output/query.corrected.fasta
ragtag.py patch ./ragtag_output/ragtag.scaffold.fasta contig.fasta
#最终“拟定的组装基因组”为:ragtag.patch.fasta
对组装结果评估:
1:基本情况:
reference全长为1,778,134bp;ragtag.patch.fasta全长为1,721,030bp,中间存在部分N的gap;要知道组装的最长的contig和scaffold不过303,902bp和303,902bp。可见使用ragtag大幅延申序列长度。
2:reads回帖率:
- reference的回帖率:34.23% overall alignment rate
-
ragtag.patch.fasta回帖率:96.86% overall alignment rate
可见,reference和我们真实的样本差别其实是很大;ragtag.patch.fasta从原始组装的contig出发,更能贴近我们测序样本的真实情况。
3:contig回帖情况:
以ragtag.patch.fasta为参考的,scaffold基本上能够回帖上去。
4:注释情况:
现在分别对reference、contig.fasta、ragtag.patch.fasta进行注释。以下为注释情况:
#reference:
organism: Genus species strain
contigs: 1
bases: 1778134
CDS: 1694
rRNA: 15
tRNA: 55
tmRNA: 1
#contig.fasta:
organism: Genus species strain
contigs: 127
bases: 1743566
CDS: 1655
rRNA: 4
tRNA: 58
tmRNA: 1
#ragtag.patch.fasta:
organism: Genus species strain
contigs: 1
bases: 1721030
CDS: 1650
tRNA: 52
tmRNA: 1
可见,使用RagTag组装之后,相比于contig.fasta,丢失了部分信息,特别是rRNA,这点在ragtag.patch.fasta中得到证实,即本该16S序列出现的region出现gap,具体原因还不清楚。
总的来说,由于是对象是细菌,无法使用HiC技术辅助,而目前市面上得到完整的细菌基因组的途径通常是二代+三代(二代负责矫错,三代延长,最终得到细菌完成图,全长无gap)。
而RagTag补长序列的方法,说白了就是利用近缘或同种的基因组为reference,将二代或三代的scaffold的按正确的顺序和方向排列,中间无法补全的用gap,可以通过锚定多个参考基因组或Hic技术来增强最后总补全的序列的准确性;这种方法本身有一定的合理性,但在缺乏图谱数据(例HiC)的情况下,可能只是旁门左道。
同时就研究目的而言,组装至染色体甚至基因组级别在研究片段倒位、穿插等中会用的多一点,而对于单基因或菌种研究中,费力通过二代测序组装全基因组可能是见吃力不讨好的事,第一,最终得到的“全长”可能会丢失部分信息(如上);第二,NCBI上Complete Genome几乎清一色来自三代数据,你用二代拼成的“全长”是否能够被认可?因此,要保证二代测序信息不被丢失,直接使用二代组装的contig/scaffold应该是最好的选择,如果非要全长,而你二代初步组装的scaffold的远达不到基因组长度,那还是老老实实用三代吧!