ALLMAPS 遗传图辅助组装
在挂载scaffold到染色体水平,scaffold的方向以及顺利是非常重要的。在该环节,可以依靠多种图谱,但是单一的图谱给的信息并不是很足,如有多个图谱信息,则会增加组装的准确性,ALLMAPS的功能就是利用多个图谱进行辅助挂载scaffold到染色体水平。
安装
Python依赖库(介意使用python3,作者将会很快淘汰python2):
上述大多数python模块可以使用pip安装
pip install --user biopython==1.70 numpy deap networkx matplotlib jcvi
简单流程
本次利用测试样本进行操作,具体数据可以点击这里, 包含3个文件:
- scaffold sequences (scaffolds.fasta)
- Map 1 (JMMale.csv)
- Map 2 (JMFemale.csv)
Step1 准备输入数据
## 格式如下,遗传图和scaffolds的对应关系
Scaffold ID, scaffold position, LG, genetic position
也可以通过excel进行获得,最后保存为csv格式
ALLMAPS支持多种图谱,具体可以看这里
Step2 合并两个图谱,得到权重文件(weights.txt)和bed文件
python -m jcvi.assembly.allmaps merge JMMale.csv JMFemale.csv -o JM-2.bed
Step3 修改权重文件,默认为1
cat weights.txt
JMFemale 1
JMMale 1
## 将重要的图谱放在第一行
Step4 scaffold的排序和方向
python -m jcvi.assembly.allmaps path JM-2.bed scaffolds.fasta
结果
得到的基因组相关文件
- JM-2.fasta 重新得到的染色体序列
- JM-2.agp AGP文件,用于上传Genbank
- JM-2.chain 用于更新scaffold坐标。比如,你用scaffold注释基因,可以利用liftOver去将其转变为染色体水平。具体戳ALLMAPS: How to lift over gene annotations
可视化比对
每条染色体对应得到一个PDF比对文件,左边为染色体和遗传图比对关系,可以非常容易看出有矛盾的marker。右边点图横坐标为染色体坐标,纵坐标为遗传图位置,右图可以非常清楚的看到连续性或者哪里出现断点。白色为可信区域。
数据汇总
JM-2.summary.txt 为一数据统计文件,该文件记录了比对的各种信息,传标记的密度,以及有多少序列被锚定到基因组,这些数据可以比较使用ALLMAPS前后的差异。
less JM-2.summary.txt
*** Summary for each individual map ***
======================================================================
o JMFemale JMMale
----------------------------------------------------------------------
Linkage Groups 1 1
Markers (unique) 61 59
Markers per Mb 4.2 4.3
N50 Scaffolds 5 5
Scaffolds 18 18
Scaffolds with 1 marker 7 7
Scaffolds with 2 markers 3 5
Scaffolds with 3 markers 1 0
Scaffolds with >=4 markers 7 6
Total bases 14,691,276 (96.2%) 13,781,640 (90.2%)
----------------------------------------------------------------------
*** Summary for consensus map ***
===================================================================================
o Anchored Oriented Unplaced
-----------------------------------------------------------------------------------
Markers (unique) 120 107 0
Markers per Mb 7.9 8.8 0
N50 Scaffolds 5 5 0
Scaffolds 22 13 0
Scaffolds with 1 marker 5 0 0
Scaffolds with 2 markers 7 3 0
Scaffolds with 3 markers 1 1 0
Scaffolds with >=4 markers 9 9 0
Total bases 15,270,543 (100.0%) 12,156,526 (79.6%) 0 (0.0%)
-----------------------------------------------------------------------------
## 可以看到对于consensus map比对率达到100%,比用任何一个单独图谱的组装效果都好
其它骚操作
拆分有矛盾的contigs
如果一些contigs存在组装错误,导致该contig比对到多个连锁群或染色体这些contig进行在运行ALLMAPS之前将其拆分。ALLMAPS可以对其进行拆分,详情请戳这里ALLMAPS: How to split chimeric contigs.
估计gap长度
在scaffold中,一般用100bpN表示gap。ALLMAPS可以根据遗传距离和物理距离对gap进行长度的估计,并得到一个新的AGP文件JM-2.estimategaps.agp
python -m jcvi.assembly.allmaps estimategaps JM-2.bed
详情戳这里ALLMAPS: How to estimate gap sizes.
利用不同类型的图谱
除了遗传图谱以外,ALLMAPS可以利用其它类型图谱。可以将多个BED文件的图谱merge为一个单独bed
python -m jcvi.assembly.allmaps mergebed opticalmap.bed geneticmap.bed synteny.bed -o merged.bed
详情戳这里 How to use different types of genomic maps.
参考
欢迎交流