获取物种分类信息的方法（TaxonKit/ete3/Biopython）

TaxonKit

中文介绍 - TaxonKit - NCBI Taxonomy Toolkit (shenwei.me)

TaxonKit是采用Go语言编写的命令行工具， 提供Linux, Windows, macOS操作系统不同架构（x86-64/arm64）的静态编译的可执行二进制文件。 发布的压缩包不足3Mb，除了Github托管外，还提供国内镜像供下载，同时还支持conda和homebrew安装。 用户只需要下载、解压，开箱即用，无需配置，仅需下载解压NCBI Taxonomy数据文件解压到指定目录即可。

• 源代码 https://github.com/shenwei356/taxonkit ，
• 文档 http://bioinf.shenwei.me/taxonkit （介绍、使用说明、例子、教程）

选择系统对应的版本下载最新版 https://github.com/shenwei356/taxonkit/releases ，解压后添加环境变量即可使用。或可选conda安装

conda install taxonkit -c bioconda -y
# 表格数据处理，推荐使用 csvtk 更高效
conda install csvtk -c bioconda -y

# 数据库准备
# 有时下载失败，可多试几次；或尝试浏览器下载此链接
wget -c https://ftp.ncbi.nih.gov/pub/taxonomy/taxdump.tar.gz 
tar -zxvf taxdump.tar.gz

# 解压文件存于家目录中.taxonkit/，程序默认数据库默认目录
mkdir -p $HOME/.taxonkit
cp names.dmp nodes.dmp delnodes.dmp merged.dmp $HOME/.taxonkit

TaxonKit为命令行工具，采用子命令的方式来执行不同功能， 大多数子命令支持标准输入/输出，便于使用命令行管道进行流水作业， 轻松整合进分析流程中。

子命令	功能
list	列出指定TaxId下所有子单元的的TaxID
lineage	根据TaxID获取完整谱系（lineage）
reformat	将完整谱系转化为“界门纲目科属种株"的自定义格式
name2taxid	将分类单元名称转化为TaxID
filter	按分类学水平范围过滤TaxIDs
lca	计算最低公共祖先(LCA)
taxid-changelog	追踪TaxID变更记录
version	显示版本信息、检测新版本
genautocomplete	生成shell自动补全配置脚本

• 输出：
  • 所有命令输出中包含输入数据内容，在此基础上增加列。
  • 所有命令默认输出到标准输出（stdout），可通过重定向（>）写入文件。
  • 或通过全局参数-o或--out-file指定输出文件，且可自动识别输出文件后缀（.gz）输出gzip格式。
• 输入：
  • 除了list与taxid-changelog之外，lineage, reformat, name2taxid, filter 与 lca 均可从标准输入（stdin）读取输入数据，也可通过位置参数（positional arguments）输入，即命令后面不带 任何flag的参数，如 taxonkit lineage taxids.txt
  • 输入格式为单列，或者制表符分隔的格式，输入数据所在列用-i或--taxid-field指定。

通过taxoid返回分类地位信息

head Test.txt
NZ_BMNY01000001.1	399579
NZ_BMNY01000002.1	399579
NZ_CP095390.1	1419050
NZ_CP095394.1	1419051
NZ_QPLR01000008.1	2282130
NZ_QPLR01000786.1	2282130


cut -f 2 Test.txt \
| taxonkit lineage \
| taxonkit reformat -F -P | cut -f 1,3 >Test_taxoid2list.tsv 

#用csvtk进行美化
cut -f 2 Test.txt \
| taxonkit lineage \
| taxonkit reformat -f "{k}\t{p}\t{c}\t{o}\t{f}\t{g}\t{s}" -F -P \
| csvtk cut -t -f -2 \
| csvtk add-header -t -n taxid,kindom,phylum,class,order,family,genus,species \
| csvtk pretty -t >Test_taxoid2list2.tsv 

选项：

-F/--fill-miss-rank：

用更高层级的分类单元信息来补齐缺失的层级

-P/--add-prefix：

给每个分类学水平添加前缀，比如s__species。

-t/--show-lineage-taxids：

输出分类学单元对应的TaxID。

-r/--miss-rank-repl: 替代没有对应rank的taxon名称

-S/--pseudo-strain: 对于低于species且rank既不是subspecies也不是stain的taxid，使用水平最低taxon名称做为菌株名称。
 

ete3

ete全称为Environment for Tree Exploration，直译就是树探索环境，此工具可以直接在终端输入
pip install ete3 进行安装即可。ete包主要功能与构建系统发生树有关，若是有相关需求可以查看其介绍文档，地址：The ETE tutorial。我主要使用到了其中的分类工具，即处理NCBI 的Taxonomy数据库的工具。此工具用于物种信息和分类号的转换十分简便，使用时是根据NCBI的最新分类文件来运作的，因此分类信息十分可靠。

下载升级数据库
ete3使用NCBI 的Taxonomy数据，因此需要首次使用时需要先下载数据。在应用中主要用到的是ete3包的NCBITaxa模块。首次使用NCBITaxa模块时会检测是否有分类数据存在，没有的话会自动下载。长期未更新时可以直接使用升级选项获取最新的分类数据。即：

from ete3 import NCBITaxa                       # 导入此模块
ncbi = NCBITaxa()
ncbi.update_taxonomy_database()                 # 升级

使用脚本

from ete3 import NCBITaxa 
ncbi = NCBITaxa()
#ncbi.update_taxonomy_database()

def taxoid2taxo(taxo):
    lineage = ncbi.get_lineage(taxo)
    names = ncbi.get_taxid_translator(lineage)
    list = []
    for taxid in lineage:
        list.append(names[taxid])
    return list
def name2taxo(name):
    try:
        taxoid = ncbi.get_name_translator([name])[name][0]
        taxo_list = taxoid2taxo(taxoid)
    except:
        taxo_list = "NA"
    return taxo_list

函数使用说明：

taxoid2taxo：通过taxoid返回taxo列表

name2taxo：通过科学名返回taxoid，再由taxoid返回taxo列表

通过ncbi号查找taxoid

from Bio import Entrez
from Bio import SeqIO
from ete3 import NCBITaxa 
ncbi = NCBITaxa()

Entrez.email = 'XXXXXX@XXXX.com'
ncid_list = [i.strip() for i in open('ncid.list')]
output = []
for ncid in ncid_list:
    handle = Entrez.efetch(db='nucleotide', id=ncid,rettype="gb", retmode="text")
    records = list(SeqIO.parse(handle, "genbank"))[0]
    name = records.annotations['organism']
    taxoid = ncbi.get_name_translator([name])[name][0]
    output.append("{}\t{}\n".format(ncid,taxoid))

with open('ncid_taxo.txt','w') as f:
    for line in output:
        f.write(line)