UniProt
现在的生物信息分析中,对编码基因序列或蛋白产物序列进行功能注释是必不可少的。不管是基因组、转录组,还是代谢组、蛋白组等等,通过生物信息分析我们都能得到一大堆基因 ID 或序列,然而仅有这些“干巴巴”的序列是不够的,最终还是要回到生物学问题上,即这些序列都有什么功能。如果只有比较少的序列,可以通过与近缘物种或模式物种的已知功能的序列进行同源比对,得到序列可能的功能信息;但很多时候我们得到了大量的未知功能的序列,只与一两个物种的功能序列进行比对,这样注释很不全面,这个时候,可以通过与各类功能数据库进行大批量的序列比对,获取序列可能的注释信息。
UniProt 是目前资源最丰富、使用频率最高的蛋白序列数据库,今天,我们就介绍下 UniProt 数据库及其使用。
1、数据库介绍
UniProt (The Universal Protein Resource) 是信息最丰富、资源最广的蛋白质序列数据库,整合 Swiss-Prot、TrEMBL 和 PIR 三大数据库的数据而成。
UniProt 主要包含 3 个部分:
(1)UniProtKB(UniProt Knowledgebase)是蛋白质序列、功能、分类、交叉引用等信息存取中心;UniProtKB 主要由两部分组成:
UniProtKB/Swiss-Prot:高质量的、手工注释的、非冗余的数据集;主要来自文献中的研究成果和 E-value 校验过计算分析结果。有质量保证的数据才被加入该数据库;
UniProtKB/TrEMBL:该数据集包含高质量的计算分析结果,一般都在自动注释中富集,主要应对基因组项目获得的大量数据流以及人工校验在时间上和人力上的不足。注释所有可用的蛋白序列。在三大核酸数据库(EMBL-Bank/GenBank/DDBJ)中注释的编码序列都被自动翻译并加入该数据库中。它也有来自 PDB 数据库的序列,以及Ensembl、Refeq 和 CCDS 基因预测的序列;
(2)UniRef(UniProt Non-redundant Reference)将密切相关的蛋白质序列组合到一条记录中,以便提高搜索速度。目前,根据序列相似程度形成 3 个子库,即 UniRef100、UniRef90 和 UniRef50;
(3)UniParc(UniProt Archive)是一个综合性的非冗余数据库,包含了所有主要的、公开的数据库的蛋白质序列。由于蛋白质可能在不同的数据库中存在,并且可能在同一个数据库中有多个版本,为了去冗余,UniaraParc 对每条唯一的序列只存一次,无论是否为同一物种的序列,只要序列相同就被合并为一条,每条序列提供稳定的、唯一的编号 UPI。该数据库含有蛋白质的序列信息,而没有注释数据。
UniProt 数据库中,UniProtKB/Swiss-Prot 和 UniProtKB/TrEMBL 数据库运用得最多,今天我们主要介绍这两个数据库的使用。
2、数据库下载
2.1、完整数据库下载
完整数据库下载示例如下:(数据库较大,可以用其他下载工具或本地下载后再上传服务器)
wget -c ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/complete/uniprot_sprot.fasta.gz
wget -c ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/complete/uniprot_trembl.fasta.gz
2.2、分类数据库下载
完整数据库很大,直接用来比对耗时很长,UniProt 官网已经将完整数据库拆分成各个子库,可以直接下载,包括古菌、细菌、真菌、植物、人、哺乳动物、脊椎动物、无脊椎动物、啮齿动物等,这样非常方便进行比对,大大缩短比对时间。下载示例如下:
wget -nH -m -c --cut-dirs=5 ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/taxonomic_divisions/
2.3、其他资源下载
有时在比对时我们还想像 NR 数据库那样输出对应物种信息等,需要下载以下数据:
taxdump
wget -c ftp://ftp.ncbi.nih.gov/pub/taxonomy/taxdump.tar.gz
wget -c ftp://ftp.ncbi.nih.gov/pub/taxonomy/taxdump.tar.gz.md5
md5sum -c taxdump.tar.gz.md5
tar -zxvf taxdump.tar.gz.md5
prot.accession2taxid
wget -c ftp://ftp.ncbi.nih.gov/pub/taxonomy/accession2taxid/prot.accession2taxid.gz
wget -c ftp://ftp.ncbi.nih.gov/pub/taxonomy/accession2taxid/prot.accession2taxid.gz.md5
md5sum -c prot.accession2taxid.gz.md5
3、格式转换
等等,下载完数据库不是就可以直接比对吗?
其实是这样的,Swiss-Prot 和 TrEMBL 完整数据库官网提供了 FASTA 格式,下载后直接建立索引、然后比对:
uniprot_sprot.fasta.gz
uniprot_trembl.fasta.gz
但分类数据库官网只提供了 DAT 格式,下载后需要将其转换成 FASTA 格式,之后再建立索引、比对:
uniprot_sprot_archaea.dat.gz
uniprot_sprot_bacteria.dat.gz
uniprot_sprot_fungi.dat.gz
uniprot_sprot_human.dat.gz
uniprot_sprot_invertebrates.dat.gz
uniprot_sprot_mammals.dat.gz
uniprot_sprot_plants.dat.gz
uniprot_sprot_rodents.dat.gz
uniprot_sprot_vertebrates.dat.gz
uniprot_sprot_viruses.dat.gz
uniprot_trembl_archaea.dat.gz
uniprot_trembl_bacteria.dat.gz
uniprot_trembl_fungi.dat.gz
uniprot_trembl_human.dat.gz
uniprot_trembl_invertebrates.dat.gz
uniprot_trembl_mammals.dat.gz
uniprot_trembl_plants.dat.gz
uniprot_trembl_rodents.dat.gz
uniprot_trembl_unclassified.dat.gz
uniprot_trembl_vertebrates.dat.gz
uniprot_trembl_viruses.dat.gz
3.1、AWK
利用 awk 进行转换,示例如下:
pigz -dc uniprot_sprot_archaea.dat.gz|awk '{if (/^ /) {gsub(/ /, ""); print} else if (/^ID/) print ">" $2}' > uniprot_sprot_archaea.fasta
这种方法比较简单,只是得到了序列 ID 和蛋白序列,没有我们想要的功能描述信息。
>1A1D_PYRAB
MHPKVDALLSRFPRITLIPWETPIQYLPRISRELGVDVYVKRDDLTGLGIGGNKIRKLEF
LLGDALSRGCDTVITIGAVHSNHAFVTALAAKKLGLGAVLILRGEEVLKGNYLLDKLMGI
ETRIYEADNSWELMKVAEEVAEELKGEGKKPYIIPPGGASPVGTLGYIRGVGELYTQVKK
LGLRIDTVVDAVGSGGTYAGLLLGSAIVNAEWSVVGIDVSSATEKAKERVKNLVEKTKEL
LGINVKVQEPRIYDYGFGAYGKIVKEVAKLIKSVGTMEGLLLDPVYTGKAFYGLMDLAKK
GDLGESVLFIHTGGLPGIFHYGEEMLELLV
3.2、Swissknife PERL module
利用 Swissknife 模块进行转换,这种方法可以生成与完整数据库序列一致的 header,推荐使用该方法。
Swissknife 无法直接用 cpan 安装,需要下载源码包单独安装,安装命令示例如下:
wget ftp://ftp.ebi.ac.uk/pub/software/swissprot/Swissknife/swissknife_1.78.tar.gz
tar -zxvf swissknife_1.78.tar.gz
cd swissknife_1.78
perl Makefile.PL
make install
安装好后可以通过以下命令进行转换:
## 单独转换古菌数据库
pigz -dc uniprot_sprot_archaea.dat.gz|perl -e 'use strict; use SWISS::Entry; use SWISS::KW; my $inputfile = $ARGV[0]; open(IN,"$inputfile") or die "Cannot open input file $inputfile: $!"; local $/ = "\n//\n"; while() {s/\r//g; my $entry = SWISS::Entry->fromText($_); print $entry->toFasta();} close IN;' - | pig -c > uniprot_sprot_archaea.fa.gz
## 批量转换
ls *.dat.gz|while read line; do pigz -dc $line | perl -e 'use strict; use SWISS::Entry; use SWISS::KW; my $inputfile = $ARGV[0]; open(IN,"$inputfile") or die "Cannot open input file $inputfile: $!"; local $/ = "\n//\n"; while() {s/\r//g; my $entry = SWISS::Entry->fromText($_); print $entry->toFasta();} close IN;' - | pigz -c > ${line%%.*}.fa.gz; done
转换后的格式如下,包括了对应序列的功能描述:
>sp|Q9V2L2|1A1D_PYRAB Putative 1-aminocyclopropane-1-carboxylate deaminase OS=Pyrococcus abyssi (strain GE5 / Orsay) OX=272844 GN=PYRAB00630 PE=3 SV=1
MHPKVDALLSRFPRITLIPWETPIQYLPRISRELGVDVYVKRDDLTGLGIGGNKIRKLEF
LLGDALSRGCDTVITIGAVHSNHAFVTALAAKKLGLGAVLILRGEEVLKGNYLLDKLMGI
ETRIYEADNSWELMKVAEEVAEELKGEGKKPYIIPPGGASPVGTLGYIRGVGELYTQVKK
LGLRIDTVVDAVGSGGTYAGLLLGSAIVNAEWSVVGIDVSSATEKAKERVKNLVEKTKEL
LGINVKVQEPRIYDYGFGAYGKIVKEVAKLIKSVGTMEGLLLDPVYTGKAFYGLMDLAKK
GDLGESVLFIHTGGLPGIFHYGEEMLELLV
4、建立索引
DIAMOND 相比 BLAST 速度更快,官网测试快 500 - 20,000 倍,同时也支持 BLAST 的输出格式,获得和 BLAST 比较一致的结果,我们选用 DIAMOND 来进行数据库比对。DIAMOND 在 github 已经编译好了,可以直接下载来用:
wget -c https://github.com/bbuchfink/diamond/releases/download/v0.9.25/diamond-linux64.tar.gz
tar -zxvf diamond-linux64.tar.gz
首先要构建数据库索引,建立索引命令示例如下:
## Swiss-Prot
/data1/software/diamond/0.9.25/diamond makedb --in uniprot_sprot.fasta.gz --db sprot --taxonmap /data3/Databases/Taxonomy/accession2taxid/prot.accession2taxid.gz --taxonnodes /data3/Databases/Taxonomy/nodes.dmp --taxonnames /data3/Databases/Taxonomy/names.dmp --threads 20
## TrEMBL
/data1/software/diamond/0.9.25/diamond makedb --in uniprot_trembl.fasta.gz --db trembl --taxonmap /data3/Databases/Taxonomy/accession2taxid/prot.accession2taxid.gz --taxonnodes /data3/Databases/Taxonomy/nodes.dmp --taxonnames /data3/Databases/Taxonomy/names.dmp --threads 20
## 分类数据库
ls *.fa.gz|while read line; do name=${line%%.*}; /data1/software/diamond/0.9.25/diamond makedb --in $line --db ${name#*_} --taxonmap /data3/Databases/Taxonomy/accession2taxid/prot.accession2taxid.gz --taxonnodes /data3/Databases/Taxonomy/nodes.dmp --taxonnames /data3/Databases/Taxonomy/names.dmp --threads 20; done
5、序列比对
经过以上步骤,现在就可以直接进行序列比对了,比对命令示例如下:(详细参数请参照 DIAMOND 说明手册)
## bacteria of Swiss-Prot
/data1/software/diamond/0.9.25/diamond blastp --db /data3/Databases/UniProt/taxonomic_divisions/sprot_bacteria --out test_pep_sprot_diamond.out --outfmt 6 qseqid qlen sseqid slen pident length mismatch gapopen qstart qend sstart send evalue bitscore staxids sscinames salltitles --query test.pep.fa --max-target-seqs 1 --evalue 1e-5 --sensitive --max-hsps 1 --threads 10 --block-size 20 --index-chunks 1 --tmpdir ./tmp
## bacteria of TrEMBL
/data1/software/diamond/0.9.25/diamond blastp --db /data3/Databases/UniProt/taxonomic_divisions/trembl_bacteria --out test_pep_trembl_diamond.out --outfmt 6 qseqid qlen sseqid slen pident length mismatch gapopen qstart qend sstart send evalue bitscore staxids sscinames salltitles --query test.pep.fa --max-target-seqs 1 --evalue 1e-5 --sensitive --max-hsps 1 --threads 10 --block-size 20 --index-chunks 1 --tmpdir ./tmp
结语
关于 UniProtKB/Swiss-Prot 和 UniProtKB/TrEMBL 的介绍和使用示例今天就介绍到这里,当然,以上命令仅供参考;此外,UniProt 数据库每月更新一次,为获取更新更全面的注释信息,大家可以将下载、转换、建立索引等相关命令打包成一个脚本,可以定期更新(如果你的存储足够的话)。
扫一扫,了解更多资讯
武汉博越致和生物科技有限公司
个人微信 :hsm664422
电话:027-87705460
您身边的多组学科研助手!