泛基因组基因组组装比较基因组

JCVI(Python版MCScan)使用记录

2022-02-09  本文已影响0人  谁冒充我小叮当

关于JCVI的使用,唐老师在github上已经有很详细的介绍,网络上也有很多教程。
我写这篇文章的主要目的是记录自己在安装及使用JCVI的过程和问题,在这中间也参考了许多帖子。
主要参考如下:
https://www.jianshu.com/p/f7971dbf5f85
https://www.jianshu.com/p/d25fb05ada97
https://www.jianshu.com/p/6f1a191a0c3a
https://zhuanlan.zhihu.com/p/353490154
https://www.omicsclass.com/article/481
以及JCVI在github上的官方说明
https://github.com/tanghaibao/jcvi/wiki/MCscan-(Python-version)

正文部分

—————————————————————————————————————

JCVI的安装

JCVI依赖的东西非常多,不管是编译还是用conda都需要额外安装别的软件,建议还是用conda安装。
我自己是在ubuntu的虚拟机上安装成功了,因此后面的演示也都在虚拟机上进行。

conda create -y -c bioconda -n JCVI jcvi
conda activate JCVI
conda install -c bioconda bedtools
conda install seqkit emboss
conda install last

我在安装Latex过程中遇到了各种麻烦,后来看到https://www.jianshu.com/p/6f1a191a0c3a这篇帖子,别人已经有完整流程了<捂脸>
这里我就拾人牙慧,简单介绍一下Latex安装方法
在有root权限的情况下

# ubuntu
sudo apt-get install -y texlive texlive-latex-extra texlive-latex-recommended
# centos
sudo yum install -y  texlive texlive-latex texlive-xetex texlive-collection-latexrecommended

没有root权限的情况下,详见https://www.jianshu.com/p/a5e9ca5886a4
用yum安装Latex后,我在使用JCVI过程中发现字体缺少类型报错,如果只是运行JCVI查看共线性结果的话,这个报错是没什么影响的。但如果是要可视化的话,还需要额外下载字体文件。用非root安装的Latex没有报错。

JCVI的运行

从网上下载了Arabidopsis_halleri和Arabidopsis_lyrata两种拟南芥的cdna序列和gff3文件


图一.png

一、GFF3转bed文件

python -m jcvi.formats.gff bed --type=mRNA --key=transcript_id \
Arabidopsis_halleri.Ahal2.2.44.gff3 > halleri_all.bed
python -m jcvi.formats.gff bed --type=mRNA --key=transcript_id \
Arabidopsis_lyrata.v.1.0.44.gff3 > lyrata_all.bed
#根据cDNA文件的序列名以及类型调整type和key
#转换完成的bed文件和fasta文件中的序列名要保持一致

二、去重复

python -m jcvi.formats.bed uniq halleri_all.bed
python -m jcvi.formats.bed uniq lyrata_all.bed
#仅保留最长且唯一mRNA序列

三、提取CDS序列

seqkit grep -f <(cut -f 4 halleri_all.uniq.bed ) \
Arabidopsis_halleri.Ahal2.2.cdna.all.fa | seqkit seq -i > halleri.cds
seqkit grep -f <(cut -f 4 lyrata_all.uniq.bed ) \
Arabidopsis_lyrata.v.1.0.cdna.all.fa | seqkit seq -i > lyrata.cds

这样就获得了每个物种唯一的cds和对应的bed文件


图二.png

四、提取共线信息

这里首先要把之前获得的四个文件的名字改成以下形式,使JCVI能够识别


图三.png
python -m jcvi.compara.catalog ortholog --no_strip_names --cscore=.99 halleri lyrata
#cscore=.99:这个选项能够有效地筛选出best hit,建议加上
#no_strip_names:是必加的参数
#dbtype:序列类型nucl还是prot

主要产生四个结果文件以及一个PDF文件的点图


图四.png

.last文件是原始的 LAST 输出,.last.filtered是过滤的 LAST 输出,.anchors是种子同线性块(高质量),.lifted.anchors招募额外的锚来形成最终的同线性块。


图五 官方文档的配图.png

五、测试同线性模式是否确实是1:1

python -m jcvi.compara.synteny depth --histogram halleri.lyrata.anchors
图六.png

六、染色体水平的局部同线性图绘制

进行这一步需要准备三个文件
第一:生成一个.simple比文件更简洁的.anchors文件

python -m jcvi.compara.synteny screen --minspan=30 --simple \
halleri.lyrata.anchors \
halleri.lyrata.anchors.new
#minspan:是规定syntenic region的总长度(i.e. 基因个数)最低阈值,默认是0
#minsize:是规定syntenic region的anchors总个数(i.e.有共线性关系的基因对数目)最低阈值,默认是0

第二:seqids文件,记录需要展示的染色体,这里由于我的示例数据中有一个不是染色体水平,因此我只选取了几条contigs

FJVB01000001.1,FJVB01000002.1,FJVB01000003.1,FJVB01000004.1,FJVB01000005.1,FJVB01000006.1,FJVB01000007.1,FJVB01000008.1,FJVB01000009.1,FJVB01000010.1,FJVB01000015.1
1,2,3,4,5,6,7,8

第三:layout文件,记录绘图的各种信息
y, xstart, xend:两个物种染色体在x轴或y轴的位置如y轴的0.6,0.4,如果是0.5,0.5两个物种的染色体就会重叠
rotation:旋转角度,0代表水平

y, xstart, xend, rotation, color, label, va,  bed
 .6,     .1,    .8,       0,      , halleri, top, halleri.bed
 .4,     .1,    .8,       0,      , lyrata, top, lyrata.bed
edges
e, 0, 1, halleri.lyrata.anchors.simple

然后就是运行绘图模块

python -m jcvi.graphics.karyotype seqid.txt layout.txt
图九.png

基因水平的局部同线性图

一、计算基因级匹配的布局

python -m jcvi.compara.synteny mcscan \
halleri.bed \
halleri.lyrata.lifted.anchors \
--iter=1 -o halleri.lyrata.blocks
#iter=1表示我们正在提取与halleri每个区匹配的对应lyrata的唯一最佳block
#查看生成的.blocks文件,它的第一列为halleri的区块,第二为lyrata区块
#如果选项--iter设置为 2,则halleri每个区将匹配对应lyrata的2个block,依此类推

二、合并bed文件

cat halleri.bed lyrata.bed > halleri_lyrata.bed

三、选取可视化区域

head -150 halleri.lyrata.blocks > blocks

四、准备layout2文件

# x,   y, rotation,   ha,     va,   color, ratio,  label
0.5, 0.6,        0, leftalign, center,       m,      1,    halleri
0.5, 0.4,      0, rightalign, center, #fc8d62,  1,    lyrata
# edges
e, 0, 1

五、运行

python -m jcvi.graphics.synteny \
--glyphcolor=orientation \
--format=png --glyphstyle=box  --shadestyle=line \
--font=Helvetica --style=darkgrid --diverge=RdYlBu \
--figsize=8x7 \
blocks \
halleri_lyrata.bed \
layout2.txt
图十.png

先这样吧,中间一些参数的意义以及更多的可视化形式有空再来补充

上一篇下一篇

猜你喜欢

热点阅读