python 版本的MCScanX分析物种之间分共线性

手把手教你多个基因组共线性作图。

所有内容可与  基因家族分析  课程中学习!

MCscan可以用于查看多个基因组(或单个基因组内)中共线性区段的关系。今天我们介绍python版的MCscan(JCVI 包),其中包装了很多的绘图功能。这里介绍具体使用方法(以葡萄和桃子两基因组分析为例)。

官方网站:https://github.com/tanghaibao/jcvi/wiki/MCscan-(Python-version)

1. 安装

这里可以分两步,分别是安装LASTAL和jcvi python包。

LASTAL下载地址: http://last.cbrc.jp/,编译后将 lastal 和 lastdb 添加到PATH中。

jcvi安装比较简单,在命令行直接运行下列命令:

pip install jcvi


2. 数据下载

运行MCscan需要序列以及坐标文件,分别是fasta格式以及BED格式的,例如,可以从Phytozome下载cds的序列文件以及gff文件,再由gff文件生成BED文件。当然JCVI有相应的方法,可以直接下载,省去我们自己下载整理数据的麻烦,具体操作如下:

$ python -m jcvi.apps.fetch phytozome

...

         Acoerulea               Alyrata             Athaliana

       Bdistachyon                 Brapa           Cclementina

           Cpapaya          Creinhardtii              Crubella

          Csativus             Csinensis Csubellipsoidea_C-169

          Egrandis                Fvesca                  Gmax

        Graimondii        Lusitatissimum            Mdomestica

        Mesculenta             Mguttatus     Mpusilla_CCMP1545

   Mpusilla_RCC299           Mtruncatula          Olucimarinus

           Osativa               Ppatens              Ppersica

      Ptrichocarpa             Pvirgatum             Pvulgaris

         Rcommunis              Sbicolor              Sitalica

     Slycopersicum       Smoellendorffii            Stuberosum

            Tcacao            Thalophila              Vcarteri

         Vvinifera                 Zmays         early_release

...

$ python -m jcvi.apps.fetch phytozome Vvinifera,Ppersica

...

$ ls

Ppersica_139_cds.fa.gz  Ppersica_139_gene.gff3.gz  Vvinifera_145_cds.fa.gz  Vvinifera_145_gene.gff3.gz


GFF转BED:

$ python -m jcvi.formats.gff bed --type=mRNA --key=Name Vvinifera_145_gene.gff3.gz -o grape.bed
$ python -m jcvi.formats.gff bed --type=mRNA --key=Name Ppersica_139_gene.gff3.gz -o peach.bed

如果参考基因组没有cds的序列文件,可下载染色体序列,根据gff文件提取cds序列。

3. 成对共线性搜索

准备好输入文件后,就可以进行共线性搜索了。首先要将工作目录切换到输入文件所在目录,然后运行如下命令:

$ python -m jcvi.compara.catalog ortholog grape peach
20:33:42 [base] lastdb peach peach.cds
20:34:13 [base] lastal -u 0 -P 64 -i3G -f BlastTab peach grape.cds >grape.peach.last
20:34:30 [synteny] Assuming --qbed=grape.bed --sbed=peach.bed
20:34:31 [blastfilter] Load BLAST file `grape.peach.last` (total 403868 lines)
20:34:31 [base] Load file `grape.peach.last`
20:34:38 [blastfilter] running the cscore filter (cscore>=0.70) ..
20:34:39 [blastfilter] after filter (294217->31229) ..
20:34:39 [blastfilter] running the local dups filter (tandem_Nmax=10) ..
20:34:39 [blastfilter] after filter (31229->21089) ..
...
A total of 30654 (NR:18661) anchors found in 678 clusters.
Stats: Min=4 Max=683 N=678 Mean=45.2123893805 SD=80.8407828815 Median=16.0 Sum=30654
NR stats: Min=4 Max=460 N=678 Mean=27.5235988201 SD=48.0461479616 Median=11.0 Sum=18661


运行结束后会得到以下几个文件

$ ls grape.peach.*
grape.peach.lifted.anchors  grape.peach.anchors  grape.peach.last.filtered  grape.peach.last

4. 成对同线性可视化

观察成对的共线性的最佳可视化方法是使用dotplot,这里只需要一条命令。

$ python -m jcvi.graphics.dotplot grape.peach.anchors


生成一个 grape.peach.pdf 点图文件,如图:

attachments-2018-07-fiD8MsHS5b472174f2f24.jpg

5. 线性可视化

我们还可以使用同样的共线性输出文件做不同的可视化图片。这里需要在准备两个输入文件,seqids 文件和 layout 文件。

前者为两个基因组都包括哪些的染色体,这里最好移除较短的scaffolds,文件格式如下所示:

chr1,chr2,chr3,chr4,chr5,chr6,chr7,chr8,chr9,chr10,chr11,chr12,chr13,chr14,chr15,chr16,chr17,chr18,chr19
scaffold_1,scaffold_2,scaffold_3,scaffold_4,scaffold_5,scaffold_6,scaffold_7,scaffold_8

其中第一行为葡萄所包含的19个染色体,第二行为桃子包含的8个染色体。

后者为画图所需的配置文件,格式如下:

# y, xstart, xend, rotation, color, label, va,  bed
 .6,     .1,    .8,       0,      , Grape, top, grape.bed
 .4,     .1,    .8,       0,      , Peach, top, peach.bed
# edges
e, 0, 1, grape.peach.anchors.simple

其中前三列指定了轨迹的位置,然后是rotation, color, label, vertical alignment (va), 和 基因组 BED 文件.最后是连线的绘制,信息来自grape.peach.anchors.simple文件。

grape.peach.anchors.simple由以下命令生成:

$ python -m jcvi.compara.synteny screen --minspan=30 --simple grape.peach.anchors grape.peach.anchors.new 


然后,我们就可以绘图了!

$ python -m jcvi.graphics.karyotype seqids layout

结果图看这里

attachments-2018-07-ZhJR5Ggu5b4721b26273f.jpg

个性化设置显示颜色,最近做的一个基因组的共线性分析:

attachments-2019-03-OMH88Fv95c93448ed646a.jpg


好了,到这里我们就结束了,想了解MCscan python版本的更多内容,可以查看MCscan 官网:https://github.com/tanghaibao/jcvi/wiki/MCscan-(Python-version)


6. 共线性结果说明


1.基因组共线性详细结果文件,过滤掉不可靠的共线性区块之后的文件:*.anchors.new,第一列为一个 基因组基因ID,第二列为另一个基因组的基因ID文件,即是两个基因组共线性区域内基因对应关系;

不同的共线性区块用“###” 隔开

attachments-2018-10-l5dIj8zk5bd7c81815ebb.jpg

2.基因组共线性简化结果文件:* anchors.simple;一行为一个共线性区块,前两列表示一个基因组中两个基因组之间的区域,与后两列(3,4列)另一个基因组中的两个基因之间的区域存在共线性关系;

第5列物区域跨度,最后一列为: +为正向,-为反向

attachments-2018-10-CEbz59gu5bd7cab39775a.jpg


3.如果要设置颜色,可修改2.中的文件,也就是在需要修改颜色的行头部加上颜色值,具体格式如下:

注意颜色值为16进制,*隔开; 没有添加颜色值的行默认灰色:

attachments-2018-10-QuMv3jGD5bd7cc4645753.jpg

添加一行后:

attachments-2018-10-BeAttuUL5bd7cc7976508.jpg



此外,我们在网易云课堂上有各种教学视频,有兴趣可以了解一下:

1. 文章越来越难发?是你没发现新思路,基因家族分析发2-4分文章简单快速,学习链接:基因家族分析实操课程

2. 转录组数据理解不深入?图表看不懂?点击链接学习深入解读数据结果文件,学习链接:转录组(有参)结果解读转录组(无参)结果解读

3. 转录组数据深入挖掘技能-WGCNA,提升你的文章档次,学习链接:WGCNA-加权基因共表达网络分析

4. 转录组数据怎么挖掘?学习链接:转录组标准分析后的数据挖掘

5. 微生物16S/ITS/18S分析原理及结果解读

6. 更多学习内容:linux、perl、R语言画图,更多免费课程请点击以下链接:

http://m.study.163.com/provider/400000000234009/index.htm?share=1&shareId=1031480318

  • 发表于 2018-07-12 17:14
  • 阅读 ( 19062 )
  • 分类:软件工具

7 条评论

请先 登录 后评论
安生水
安生水

252 篇文章

作家榜 »

  1. omicsgene 541 文章
  2. 安生水 252 文章
  3. Daitoue 167 文章
  4. 生物女学霸 120 文章
  5. 红橙子 78 文章
  6. CORNERSTONE 72 文章
  7. 生信老顽童 48 文章
  8. landy 37 文章