创建NR子库以及从NR库提取特定物种分类的序列

  1. 首先是软件的安装

    • 下载NCBI的TaxonKit软件,http://bioinf.shenwei.me/taxonkit/download/,linux系统直接解压,接着:

      将taxonkit放到环境变量中

      sudo cp taxonkit /usr/local/bin/
      

      将两个文件(names.dmp和nodes.dmp在ftp://ftp.ncbi.nih.gov/pub/taxonomy/taxdump.tar.gz下载)复制到用户目录下的隐藏文件夹.taxonkit中

      cp names.dmp ~/.taxonkit
      cp nodes.dmp ~/.taxonkit
      

      然后就能正常使用了

    • 下载NCBI的csvtk软件,http://bioinf.shenwei.me/csvtk/download/,linux系统也是直接解压,即可使用

    • 还有一个数据文件需要下载,里面有NCBI的accession与taxid的对应关系,prot.accession2taxid.gz

  2. 使用TaxonKit提取特定taxons下的所有taxid,以植物taxid 33090为例,–ids是指定taxid,–indent “”是将所列出的taxid左边的空格去除,以左对齐排列;然后也可以wc命令查看植物下有多少个taxid(有189153个taxids)

    taxonkit list --ids 33090 --indent "" > plant.taxid.txt
    wc -l plant.taxid.txt
    
  3. 使用csvtk在prot.accession2taxid.gz文件中提取plant.taxid所有的accession,参数可以在http://bioinf.shenwei.me/csvtk/usage/查询到

    zcat prot.accession2taxid.gz |csvtk -t grep -f taxid -P plant.taxid.txt |csvtk -t cut -f accession.version >plant.taxid.acc.txt
    

    简单的说,感觉csvtk就是一个脚本集合工具,因此这步通过自己写的脚本也能实现的。大概1分钟后就能获得plant下所有的accession号,wc查下总共有8652634条

    但是!!!如果只是想获得plant下的所有的accession号,其实上面3步都没必要做,因此NCBI网站上就可以直接下载,比如按照帖子http://www.biotrainee.com/thread-1818-1-1.html上的操作,在NCBI -> protein下检索txid33090[Organism] 即可显示出plant下所有的蛋白信息,然后单击Sent to -> File ->Accession List -> Creat File。所用时间可能还比上面的方法还要快。。。总共有8709559条,与上面的略有差异

    但是你是其他大类物种,比如细菌等,并且还要提取好几个物种的accession号,那么上述的方法将会节约你不少时间

  4. 如果是想建立NR子库,那么利用上面的得到的plant.taxid.acc.txt文件,使用blast+配套的程序即可,如:

    blastdb_aliastool -gilist plant.taxid.acc.txt -db nr -out nr_plant -title nr_plant       ##gi号
    
    blastdb_aliastool -seqidlist plant.taxid.acc.txt -db nr -out nr_plant -title nr_plant    ##accesstion号
    
  5. 如果是想提取特定物种(比如植物)下的所有NR序列,那么可以按照http://bioinf.shenwei.me/taxonkit/tutorial/的方法,主要的也是利用blast+的blastdbcmd工具:

    blastdbcmd -db nr -entry all -outfmt "%a\t%T" |csvtk -t grep -f 2 -P plant.taxid.acc.txt |csvtk -t cut -f 1 |blastdbcmd -db nr -entry_batch - -out nr.plant.fa