1.安装

1.1 依赖软件安装

MCScan github链接

依照github上的说明需要以下软件:

需要的python包:

 

1.2本体安装

MCScan其实就是python的一个包(jcvi),安装只需要用pip就可以。

pip install jcvi

MCScan是支持python3的,且马上要不支持python2了,所以最好使用python3。

如果安装失败可以尝试使用conda安装python,再使用pip安装。

或下载源代码安装:

cd ~/code  # or any directory of your choice
git clone git://github.com/tanghaibao/jcvi.git
pip install -e .

 

  1. 示例文件运行

其实wiki上这部分写的很详细了,可以直接去上面看。

https://github.com/tanghaibao/jcvi/wiki/MCscan-(Python-version)

2.1 准备数据

MCScan需要的源文件是 cds序列文件(fasta格式)和gff3文件。且最好在phytozome上下载,用这个包自带的命令只能提phytozome数据库的数据。

示例使用的是葡萄,桃子吧(Vvinifera,Ppersica)。在phytozome上下载。需要账户,注册一个就好了。

python -m jcvi.apps.fetch phytozome
Phytozome Login: xxxxxxxx
Phytozome Password:
python -m jcvi.apps.fetch phytozome Vvinifera,Ppersica

这样就有了两个cds.fa和gff3文件(压缩)。

利用gff3文件生成.bed文件:

python -m jcvi.formats.gff bed --type=mRNA --key=Name Vvinifera_145_Genoscope.12X.gene.gff3.gz -o grape.bed
python -m jcvi.formats.gff bed --type=mRNA --key=Name Ppersica_298_v2.1.gene.gff3.gz -o peach.bed

生成cds文件:

python -m jcvi.formats.fasta format Vvinifera_145_Genoscope.12X.cds.fa.gz grape.cds
python -m jcvi.formats.fasta format Ppersica_298_v2.1.cds.fa.gz peach.cds

 

ls *.???
grape.bed grape.cds peach.bed peach.cds

这样就有了.cds文件 和.bed文件,后续共线性分析主要用这四个文件。

2.2 共线性分析

python -m jcvi.compara.catalog ortholog grape peach 

这样就生成了如下几个文件。

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

grape.peach.pdf:

 

2.3 绘图

在得到grape.peach.anchors 后就可以进行之后的绘图。

 

需要准备3个文件

第1个是 .simple 文件,由.anchors 生成。

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

第2个是 seqids 文件。里面是染色体的名字:

chr1,chr2,chr3,chr4,chr5,chr6,chr7,chr8,chr9,chr10,chr11,chr12,chr13,chr14,chr15,chr16,chr17,chr18,chr19
Pp01,Pp02,Pp03,Pp04,Pp05,Pp06,Pp07,Pp08

注意染色体名字必须在之前的bed文件中存在。

第3个是 layout 文件。里面是绘图的相关设定。

# 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

准备好这3个文件就可以绘图了:

python -m jcvi.graphics.karyotype seqids layout

生成了一个新的pdf文件:

这样简单地使用示例文件进行共线性分析及绘图就完成了。

 

scripts

  1. get_gene_bed.pl

use Getopt::Long;
my %opts;
use Data::Dumper;
GetOptions( \%opts, "in1=s","out=s","h");
if(!defined($opts{in1}) || !defined($opts{out})||defined($opts{h})){
	$USAGE;
}
open (IN1,"$opts{in1}")or die "cannot open the file\n";
open (OUT,">$opts{out}")or die "cannot open the file\n";
while(<IN1>){
	chomp;
	my @a = split /\t/,$_;
	if($a[2] eq "gene"){
		#if($a[2] eq "mRNA"){
			$a = ~m/ID=([^;]*)/;
			$id =$1;
			#$id=~s/\.TAIR10//;	
			print OUT "$a[0]\t$a[3]\t$a[4]\t$id\t$a[7]\t$a[6]\n";
	
	}
}



 

  1. get_fa_bed.pl

use Math::BigFloat;
use Bio::SeqIO;
use Bio::Seq;


#读入蛋白质序列
$in=Bio::SeqIO->new(
	-file => "$ARGV[1]",
	-format => 'Fasta'
);


#输出蛋白序列
$out = Bio::SeqIO->new(
	-file => ">$ARGV[2]",
	-format => 'Fasta'
);


##提取基因ID
my %keep=();
open IN,"$ARGV[0]"or die "cannot open the file";
while (<IN>){
	chomp;
	next if /^#/;
	my @a=split /\t/;
    ###下一行适当更改
	$keep{"$a[3].1"}=1;
	
}
close IN;
while(my $seq=$in->next_seq()){
	my ($id,$sequence,$desc)=($seq->id, $seq->seq, $seq->desc);
	if(exists $keep{$id}){
		$out->write_seq($seq);
	}
}
$in->close();
$out->close();