串行的blast+2.2.25本地化: blastdbcmd,makeblastdb,blastn.简介

串行的blast+2.2.25本地化系列: blastdbcmd,makeblastdb,blastn.使用简介

-----------------------------------------------------------------------

 

blast+中有许多命令.

blastdbcmd,makeblastdb是处理数据库较为常用的命令

blastn是核苷酸与核苷酸比对的命令.

以下结果来自 blast+ 2.2.25,

运行平台是 ubuntu 10.04

 GCC 4.4.3

------------------------------------------------------------------------

1 blastdbcmd 

 

• 作用: Retrieves sequences or other information from a BLAST database

 

它相当于以前的fastacmd.利用这个命令,可以从一个blast数据库中获得你想要的信息:

一般的使用如:

 例子0:

blastdbcmd -db refseq_rna -info  

可以查看数据库refseq_rna的信息

例子1:

blastdbcmd -db refseq_rna -entry 224071016 -out test.fa

 

可以从数据库中提取除gi号为224071016的序列,并且以fasta格式存入文件(当然也可以以其它格式获得序列)

 

• 注:

gi ID是许多用来标志序列的标识符中的一种.是数据库文件中普遍使用,通行有效的保持索引的形式.

所有来源于NCBI的序列都有一个gi号“gi|gi_identifier”.是绝对唯一的.

 

而自己利用makeblastdb命令 构建的数据库中,利用以下三种标识符

1.   gnl|database|identifier
2.   lcl|identifier
3.  identifier

这些标识符的作用是

1. 区别于gi号
2. 在本数据库中使得序列标识符唯一
3. 在查询和比对中分辨query序列与subject序列

• 更加详细的参数参照 blastdbcmd -help

------------------------------------------------------------------------------------------

2 makeblastdb

 

• 作用: Formats input FASTA file(s) into a BLAST database

 

顾名思义,它的作用是将各种格式的文件转换为一个可供blast算法使用的blast数据库

• 一般的使用如:

makeblastdb -in input_file -dbtype molecule_type -input_type type_string -title database_title -parse_seqids -out database_name -logfile File_Name

 

-in 后接输入文件，你要格式化的序列

-dbtype 后接序列类型，nucl为核酸，prot为蛋白

-title 给数据库的别名

-input_type 给出输入文件的序列格式.

-parse_seqids 将序列分列

-out 后接数据库名，自己起一个有意义的名字，以后blast+搜索时要用到的-db的参数

-logfile 日志文件，如果没有默认输出到屏幕

 

• 更加详细的参数参照 makeblastdb -help

----------------------------------------------------------

3 blastn简单介绍

• 作用: 核苷酸序列的比对以及结果分析

 

• 例子0: 最简单的使用格式

blastn -query test_query.fa -db refseq_rna -out test_results

 

• 例子1: 

blastdbcmd -db refseq_rna -entry nm_000249 -out  test_query.fa   blastn -query test_query.fa -db refseq_rna -task blastn -dust no -outfmt 7 -num_alignments 2 -num_descriptions 2  

 

-task 规定搜索采用的策略:可选为'blastn' 'blastn-short' 'dc-megablast'  'megablast' 'vecscreen'. 默认为`megablast'

• 例子2: 常用的格式:
• blastn -db database_name -query input_file -out output_file -evalue evalue -max_target_seqs num_sequences -num_threads int_value -outfmt "7 qacc sacc evalue length pident"
  
  

-evalue 预期的从hit筛选hsp的阈值

 

-outfmt 选择性输出的一些内容. 如 : 7 代表查询结果表格化并且带有注释行. 而 qacc 以及 sacc等则是自定义的所需的信息.可供选择的信息选项详见 -outfmt参数.-outfmt直接触发分析程序,而不用另外使用分析程序了.

 

 

更加详细的参数参照 blastn -help

------------------------------------------------------------------------

经过观察发现: 用于blastn的query序列使用 fasta格式即可. 而数据库则不能直接使用fasta格式.而且数据库实际上是一个文件夹,而非单个文件.

可以使用makeblastdb将fasta格式的文件转化为数据库

另一方面可以用blastdbcmd从数据库中获取信息

 

-------------------------------------------------------------------------

makeblastdb,blastdbcmd,blastn使用示例:

 

/**** 

fasta文件的格式,形如: 

> identifier1  

ATCG...ATCG

>identifier2

ATCG...ATCG

...

...

...

> identifier3

ATCG...ATCG

****/

 

现有fasta文件test.fasta如下:

>seq1

ATGTTCAACGCGAAGAACGGTTTTTCTGAGGCACACGTGAGGGGATGTCAGACCAAACGACTCACCAAACAGAACTACGC

CGAACTTTCTCGATGTGACACGTTGGAAGACATCAAGACGTACTTGCAAACGATGAGTGATTATTCAGAATATGTTCGTG

ATCTTCAAGCGCCAGTGAGACCGGTTGACATTATTGAATGCTGCAGAAAGAGACAGATCGCAGAGTTTAATATTTGCTGT

CAGCAGGCTTCTTCCCCTTTGTCCAATTTTTTGGAGTATTTGACGTACGGATACATGATCGATAATCTTGTGTTGGCTTT

>seq2

AAATGGCATGCTTCGTGGACGTACCACAGAGGCAATACTTGAGAAGTGTAGCCCCATTGGTTTTTTCGATTCTTTATCCG

CGGTTGTCGTGTCGAGTAGTGTCCAAGAACTCTACAGACTAGCTCTCGTGGATACACCGCTTGCCTCTTATTTCAGTAGC

TCGATTAAGGCAGAAGATCTGGATGAGTTAAATATTGAGCTCATACGGAACGTCCTATACAAGGAATATTTGCAAGATTT

CATGGTTTTCTGCAACAAAATGGATCAAAACACACGTCAATTGATGGAGAAACTACTTAGCATGGAGGCCGATCGGCACG

>seq3

CGATAAGAATCACACTGAACTCTTTCGGAACAGAGCTTTCCAAGGCTGATCGAAGAAATCTTTATACGAATTTTGGCACC

ATGTACCCCGATGGCTTCGCGCGTCTTGCGAATTGTGAAACGGTAGATGAAGTGAAACGCATACTAGTAGCTTATCCAGA

ATTCAGAGAGTTGACGAAAAGTGATGATCCCCACTACATTGACAGGGGACTACGCGTTCTCGAACTGGAAGCATGTGGAC

AAGCACTCGATGAGCAATTCAATTTCGCTATCTTTTATGCTTTCGTAAAGTTTCAGGAGAACGAAATAAACAACCTGATG

TGGCTCACTGAGTGTGTTGCTCAAAGGCAAAAAAGTAGTCTAGGCGAGGGCATTGTCTACATACAATAG

 

/***

注意:自己生成的数据库中序列命名有以下三种形式：

a） > gnl|database|identifier

b） > lcl|identifier

c)  > identifier

***/

 

1 makeblastdb使用示例:

1)用命令将test.fasta转化为数据库test_db:makeblastdb -in test.fasta -dbtype nucl -title lidachao -parse_seqids -out test_db

oboyo@oboyo-laptop:~/blast-lastest$ makeblastdb -in test.fasta -dbtype nucl -title lidachao -parse_seqids -out test_db

 

Building a new DB, current time: 05/19/2011 17:31:55

New DB name:   test_db

New DB title:  lidachao

Sequence type: Nucleotide

Keep Linkouts: T

Keep MBits: T

Maximum file size: 1073741824B

Adding sequences from FASTA; added 3 sequences in 0.00573689 seconds.

 

2)格式化数据库后，创建三个主要的文件——库索引(indices)，序列(sequences)和头(headers)文件。生成的文件的扩展名分别是：.pin、.psq、.phr（对蛋白质序列）或.nin、.nsq、.nhr（对核酸序列）。而其他的序列识别符和索引则包含在.psi和.psd（或.nsi和.nsd）中。

  该示例是核苷酸所以生成的数据库文件是: nin , nsq ,nhr , .nsi ,.nsd ,.nog

 

2 blastdbcmd使用示例

1)查看生成的test_db数据库的信息

1. oboyo@oboyo-laptop:~/blast-lastest$ blastdbcmd -db test_db -info  

 

Database: lidachao

3 sequences; 1,029 total bases

 

Date: May 19, 2011  5:56 PM Longest sequence: 389 bases

 

Volumes:

/home/oboyo/blast-lastest/test_db

 

 

2)从test_db中抽取一段 数据库的序列(seq1),并以fasta形式显示

oboyo@oboyo-laptop:~/blast-lastest$ blastdbcmd -db test_db -entry seq1  

 

>lcl|seq1 

ATGTTCAACGCGAAGAACGGTTTTTCTGAGGCACACGTGAGGGGATGTCAGACCAAACGACTCACCAAACAGAACTACGC

CGAACTTTCTCGATGTGACACGTTGGAAGACATCAAGACGTACTTGCAAACGATGAGTGATTATTCAGAATATGTTCGTG

ATCTTCAAGCGCCAGTGAGACCGGTTGACATTATTGAATGCTGCAGAAAGAGACAGATCGCAGAGTTTAATATTTGCTGT

CAGCAGGCTTCTTCCCCTTTGTCCAATTTTTTGGAGTATTTGACGTACGGATACATGATCGATAATCTTGTGTTGGCTTT

 

test_db是自定义的数据库.利用自己的标识符"seq1",可以成功抽取了在数据库中标识符为 lcl|seq1 的序列.

在使用下载的数据库时同样可以使用gi ID 抽取序列.

 

3)将抽取的序列以fasta格式保存到 test_query文件

1. oboyo@oboyo-laptop:~/blast-lastest$ blastdbcmd -db test_db -entry seq1 -out test_query  

 

 

 

3 blastn使用示例.

1)利用上面生成 fasta格式的 test_query与 自定义的数据库 test_db进行比对.

1. oboyo@oboyo-laptop:~/blast-lastest$ blastn -query test_query -db test_db  

 

BLASTN 2.2.25+

 

 

Reference: Zheng Zhang, Scott Schwartz, Lukas Wagner, and Webb

Miller (2000), "A greedy algorithm for aligning DNA sequences", J

Comput Biol 2000; 7(1-2):203-14.

 

 

 

Database: lidachao

           3 sequences; 1,029 total letters

 

 

 

Query= lcl|seq1

 

Length=320

                                                                      Score     E

Sequences producing significant alignments:                          (Bits)  Value

 

lcl|seq1                                                               592    2e-173

 

 

>lcl|seq1 

Length=320

 

 Score =  592 bits (320),  Expect = 2e-173

 Identities = 320/320 (100%), Gaps = 0/320 (0%)

 Strand=Plus/Plus

 

Query  1    ATGTTCAACGCGAAGAACGGTTTTTCTGAGGCACACGTGAGGGGATGTCAGACCAAACGA  60

            ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

Sbjct  1    ATGTTCAACGCGAAGAACGGTTTTTCTGAGGCACACGTGAGGGGATGTCAGACCAAACGA  60

 

Query  61   CTCACCAAACAGAACTACGCCGAACTTTCTCGATGTGACACGTTGGAAGACATCAAGACG  120

            ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

Sbjct  61   CTCACCAAACAGAACTACGCCGAACTTTCTCGATGTGACACGTTGGAAGACATCAAGACG  120

 

Query  121  TACTTGCAAACGATGAGTGATTATTCAGAATATGTTCGTGATCTTCAAGCGCCAGTGAGA  180

            ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

Sbjct  121  TACTTGCAAACGATGAGTGATTATTCAGAATATGTTCGTGATCTTCAAGCGCCAGTGAGA  180

 

Query  181  CCGGTTGACATTATTGAATGCTGCAGAAAGAGACAGATCGCAGAGTTTAATATTTGCTGT  240

            ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

Sbjct  181  CCGGTTGACATTATTGAATGCTGCAGAAAGAGACAGATCGCAGAGTTTAATATTTGCTGT  240

 

Query  241  CAGCAGGCTTCTTCCCCTTTGTCCAATTTTTTGGAGTATTTGACGTACGGATACATGATC  300

            ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||

Sbjct  241  CAGCAGGCTTCTTCCCCTTTGTCCAATTTTTTGGAGTATTTGACGTACGGATACATGATC  300

 

Query  301  GATAATCTTGTGTTGGCTTT  320

            ||||||||||||||||||||

Sbjct  301  GATAATCTTGTGTTGGCTTT  320

 

 

 

Lambda     K      H

    1.33    0.621     1.12 

 

Gapped

Lambda     K      H

    1.28    0.460    0.850 

 

Effective search space used: 307764

 

 

  Database: lidachao

    Posted date:  May 19, 2011  5:56 PM

  Number of letters in database: 1,029

  Number of sequences in database:  3

 

 

 

Matrix: blastn matrix 1 -2

Gap Penalties: Existence: 0, Extension: 2.5

 

2)选定模式进行输出

1. oboyo@oboyo-laptop:~/blast-lastest$ blastn -query test_query -db test_db -task blastn -dust no -outfmt 7 -num_alignments 2 -num_descriptions 2  

 

# BLASTN 2.2.25+

# Query: lcl|seq1

# Database: test_db

# Fields: query id, subject id, % identity, alignment length, mismatches, gap opens, q. start, q. end, s. start, s. end, evalue, bit score

# 1 hits found

lcl|seq1 seq1 100.00 320 0 0 1 320 1 320 2e-169 578

# BLAST processed 1 queries

转自：http://blog.csdn.net/lidachao1/article/details/6433053