ncbi批量下载基因家族，并按种族分类

以NBS-LRR基因家族为例
1、下载数据
从核苷酸数据库下载所有的NBS-LRR基因的概要文件
如果fasta数据下载，太慢了。如果网络好可以全部下载下来，然后本地写脚本

2、从概要中提取基因ac号，并分类

#! /usr/bin/perl
open IN1,"$ARGV[0]" or die "can`t open input file!\n ";
#第一步，提取所有种类和数量
while(<IN1>){
    chown;
    if(/\[(.+?)\]/){
        $species{$1}++;     
    }
}
#打印结果
close IN1;
map{printf "%-45s => %-s \n",$_,$species{$_};}sort keys %species;
#第二步，提取各物种的nbs-lrr的蛋白质的ac号
open IN2,"$ARGV[0]" or die "can`t open input file!\n ";
while(<IN2>){
    chown;
    if(/\[(.+?)\]/){
        $last=$1;
        $row=0;
        next;
    }
    $row++;
    if($row == 2){  
        if(/(.*?)\s/){
            $seq{$last}="$1\t$seq{$last}";
        }
    }
}
close IN2;
#打印结果并输出到文件
open OUT ,"> 2.out" or die "can`t not open output file!\n";
for(sort keys %seq){
    print OUT ">$_\n$seq{$_}\n";
}
close OUT;

图片.png

3、利用分类的基因ac号，找到感兴趣的物种，将其ac号复制到另外一个文件中去。整个为一行中间用\t隔开。

图片.png

4、用这个文件上传到ncbi下载序列
网址

https://www.ncbi.nlm.nih.gov/sites/batchentrez

图片.png
图片.png

5、如果你下载了全部序列，则可以自己写个脚本，本地提取更方便。

#! /usr/bin/perl -w

open IN,"$ARGV[0]" or die "$!";
open IN2,"$ARGV[1]" or die "$!";
open OUT,"> ac_seq.fasta";

map{$ac .= "$_"; }<IN>;
map{
    if(/>(\S+)/){
        $flag=0;
        $count+=$flag=1 if($ac=~/$1/)
    }
    print OUT "$_"if(1==$flag)

}<IN2>;
print"$count";

图片.png