云之南

风声，雨声，读书声，声声入耳；家事，国事，天下事，事事关心

日志

关于我

云之南

专业背景：计算机科学研究方向与兴趣: JavaEE-Web软件开发, 生物信息学, 数据挖掘与机器学习, 智能信息系统目前工作: 基因组, 转录组, NGS高通量数据分析, 生物数据挖掘, 植物系统发育和比较进化基因组学

文章分类

BLAST+实战

2010-07-19 16:00:31| 分类：生物信息学 | 标签： |举报 |字号大中小订阅

下载LOFTER 我的照片书 |

BLAST+实战
2010-04-29 16:08

http://hi.baidu.com/lidaof/blog/item/9d29eade7f340a5895ee3724.html

项目要求，开始使用BLAST+

我们从格式化数据库开始

1，我先获得了一个植物rRNA的序列文件，为fasta格式，用makeblastdb格式化之，如下

makeblastdb -in plant26.rna.fna -dbtype nucl -parse_seqids -hash_index -out plant_rna

注意：在BLAST+2.2.24中不加参数 -parse_seqids ，-hash_index 不然成死循环

Building a new DB, current time: 04/29/2010 16:09:11
New DB name: plant_rna
New DB title: plant26.rna.fna
Sequence type: Nucleotide
Keep Linkouts: T
Keep MBits: T
Maximum file size: 1073741824B
Adding sequences from FASTA; added 6747 sequences in 0.546889 seconds.

此时目录下还是如往常一下生成如下文件

plant_rna.nhr plant_rna.nin plant_rna.nnd plant_rna.nni plant_rna.nsd plant_rna.nsi plant_rna.nsq

2，我们有如下一个fasta序列想做blast,文件如下

more test
>gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probab
le (VaoD) mRNA, complete cds
ATGTTCAACGCGAAGAACGGTTTTTCTGAGGCACACGTGAGGGGATGTCAGACCAAACGACTCACCAAACAGAACTACGC
CGAACTTTCTCGATGTGACACGTTGGAAGACATCAAGACGTACTTGCAAACGATGAGTGATTATTCAGAATATGTTCGTG
ATCTTCAAGCGCCAGTGAGACCGGTTGACATTATTGAATGCTGCAGAAAGAGACAGATCGCAGAGTTTAATATTTGCTGT
CAGCAGGCTTCTTCCCCTTTGTCCAATTTTTTGGAGTATTTGACGTACGGATACATGATCGATAATCTTGTGTTGGCTTT
AAATGGCATGCTTCGTGGACGTACCACAGAGGCAATACTTGAGAAGTGTAGCCCCATTGGTTTTTTCGATTCTTTATCCG
CGGTTGTCGTGTCGAGTAGTGTCCAAGAACTCTACAGACTAGCTCTCGTGGATACACCGCTTGCCTCTTATTTCAGTAGC
TCGATTAAGGCAGAAGATCTGGATGAGTTAAATATTGAGCTCATACGGAACGTCCTATACAAGGAATATTTGCAAGATTT
CATGGTTTTCTGCAACAAAATGGATCAAAACACACGTCAATTGATGGAGAAACTACTTAGCATGGAGGCCGATCGGCACG
CGATAAGAATCACACTGAACTCTTTCGGAACAGAGCTTTCCAAGGCTGATCGAAGAAATCTTTATACGAATTTTGGCACC
ATGTACCCCGATGGCTTCGCGCGTCTTGCGAATTGTGAAACGGTAGATGAAGTGAAACGCATACTAGTAGCTTATCCAGA
ATTCAGAGAGTTGACGAAAAGTGATGATCCCCACTACATTGACAGGGGACTACGCGTTCTCGAACTGGAAGCATGTGGAC
AAGCACTCGATGAGCAATTCAATTTCGCTATCTTTTATGCTTTCGTAAAGTTTCAGGAGAACGAAATAAACAACCTGATG
TGGCTCACTGAGTGTGTTGCTCAAAGGCAAAAAAGTAGTCTAGGCGAGGGCATTGTCTACATACAATAG

3，blast搜索命令格式如下：

blastn -db plant_rna -query test -out test.blastn

4，我们来看看输入文件test.blastn的内容

BLASTN 2.2.23+

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: plant26.rna.fna
6,747 sequences; 8,392,753 total letters

Query= gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901
vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
Length=1029
Score E
Sequences producing significant alignments: (Bits) Value

ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar AT... 1901 0.0

>ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit
D, probable (VaoD) mRNA, complete cds
Length=1029

Score = 1901 bits (1029), Expect = 0.0
Identities = 1029/1029 (100%), Gaps = 0/1029 (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   GATAATCTTGTGTTGGCTTTAAATGGCATGCTTCGTGGACGTACCACAGAGGCAATACTT 360
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 301   GATAATCTTGTGTTGGCTTTAAATGGCATGCTTCGTGGACGTACCACAGAGGCAATACTT 360

Query 361   GAGAAGTGTAGCCCCATTGGTTTTTTCGATTCTTTATCCGCGGTTGTCGTGTCGAGTAGT 420
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 361   GAGAAGTGTAGCCCCATTGGTTTTTTCGATTCTTTATCCGCGGTTGTCGTGTCGAGTAGT 420

Query 421   GTCCAAGAACTCTACAGACTAGCTCTCGTGGATACACCGCTTGCCTCTTATTTCAGTAGC 480
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 421   GTCCAAGAACTCTACAGACTAGCTCTCGTGGATACACCGCTTGCCTCTTATTTCAGTAGC 480

Query 481   TCGATTAAGGCAGAAGATCTGGATGAGTTAAATATTGAGCTCATACGGAACGTCCTATAC 540
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 481   TCGATTAAGGCAGAAGATCTGGATGAGTTAAATATTGAGCTCATACGGAACGTCCTATAC 540

Query 541   AAGGAATATTTGCAAGATTTCATGGTTTTCTGCAACAAAATGGATCAAAACACACGTCAA 600
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 541   AAGGAATATTTGCAAGATTTCATGGTTTTCTGCAACAAAATGGATCAAAACACACGTCAA 600

Query 601   TTGATGGAGAAACTACTTAGCATGGAGGCCGATCGGCACGCGATAAGAATCACACTGAAC 660
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 601   TTGATGGAGAAACTACTTAGCATGGAGGCCGATCGGCACGCGATAAGAATCACACTGAAC 660

Query 661   TCTTTCGGAACAGAGCTTTCCAAGGCTGATCGAAGAAATCTTTATACGAATTTTGGCACC 720
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 661   TCTTTCGGAACAGAGCTTTCCAAGGCTGATCGAAGAAATCTTTATACGAATTTTGGCACC 720

Query 721   ATGTACCCCGATGGCTTCGCGCGTCTTGCGAATTGTGAAACGGTAGATGAAGTGAAACGC 780
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 721   ATGTACCCCGATGGCTTCGCGCGTCTTGCGAATTGTGAAACGGTAGATGAAGTGAAACGC 780

Query 781   ATACTAGTAGCTTATCCAGAATTCAGAGAGTTGACGAAAAGTGATGATCCCCACTACATT 840
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 781   ATACTAGTAGCTTATCCAGAATTCAGAGAGTTGACGAAAAGTGATGATCCCCACTACATT 840

Query 841   GACAGGGGACTACGCGTTCTCGAACTGGAAGCATGTGGACAAGCACTCGATGAGCAATTC 900
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 841   GACAGGGGACTACGCGTTCTCGAACTGGAAGCATGTGGACAAGCACTCGATGAGCAATTC 900

Query 901   AATTTCGCTATCTTTTATGCTTTCGTAAAGTTTCAGGAGAACGAAATAAACAACCTGATG 960
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 901   AATTTCGCTATCTTTTATGCTTTCGTAAAGTTTCAGGAGAACGAAATAAACAACCTGATG 960

Query 961   TGGCTCACTGAGTGTGTTGCTCAAAGGCAAAAAAGTAGTCTAGGCGAGGGCATTGTCTAC 1020
             ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
Sbjct 961   TGGCTCACTGAGTGTGTTGCTCAAAGGCAAAAAAGTAGTCTAGGCGAGGGCATTGTCTAC 1020

Query 1021 ATACAATAG 1029
|||||||||
Sbjct 1021 ATACAATAG 1029

Lambda K H
1.33 0.621 1.12

Gapped
Lambda K H
1.28 0.460 0.850

Effective search space used: 8286997432

Database: plant26.rna.fna
Posted date: Apr 29, 2010 4:09 PM
Number of letters in database: 8,392,753
Number of sequences in database: 6,747

Matrix: blastn matrix 1 -2
Gap Penalties: Existence: 0, Extension: 0

5，我们可以看出，如果不加任何参数，和以前的输出差不多相同，大家还是可以使用以前写的bioperl或biopython parser进行结果的分析

6,下面我们看看增加输出参数呢，我们使用-outfmt 7作为控制输出结果，其中选择性输出部分内容，命令和输出如下：

blastn -db plant_rna -query test -outfmt "7 qacc sacc evalue length pident"
# BLASTN 2.2.23+
# Query: gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
# Database: plant_rna
# Fields: query acc., subject acc., evalue, alignment length, % identity
# 1 hits found
gi|145342129|ref|XM_001416109.1| XM_001416109 0.0 1029 100.00
# BLAST processed 1 queries

这个文件是自解释的，大家可以看得很清楚每项是什么意思，呵呵，其他输出的一些部分特性如下：

blastn -db plant_rna -query test -outfmt "7 qacc sacc evalue length pident"
# BLASTN 2.2.23+
# Query: gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
# Database: plant_rna
# Fields: query acc., subject acc., evalue, alignment length, % identity
# 1 hits found
gi|145342129|ref|XM_001416109.1| XM_001416109 0.0 1029 100.00
# BLAST processed 1 queries

blastn -db plant_rna -query test -outfmt "7 qgi sgi evalue length pident" # BLASTN 2.2.23+
# Query: gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
# Database: plant_rna
# Fields: query gi, subject gi, evalue, alignment length, % identity
# 1 hits found
0 145342129 0.0 1029 100.00
# BLAST processed 1 queries

blastn -db plant_rna -query test -outfmt "7 qid sid evalue length pident"
# BLASTN 2.2.23+
# Query: gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
# Database: plant_rna
# Fields: evalue, alignment length, % identity
# 1 hits found
0.0 1029 100.00
# BLAST processed 1 queries

blastn -db plant_rna -query test -outfmt "7 qseqid sseqid evalue length pident"
# BLASTN 2.2.23+
# Query: gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
# Database: plant_rna
# Fields: query id, subject id, evalue, alignment length, % identity
# 1 hits found
gi|145342129|ref|XM_001416109.1| gi|145342129|ref|XM_001416109.1| 0.0 1029 100.00
# BLAST processed 1 queries

blastn -db plant_rna -query test -outfmt "7 qframe sseqid evalue length pident"
# BLASTN 2.2.23+
# Query: gi|145342129|ref|XM_001416109.1| Ostreococcus lucimarinus CCE9901 vacuolar ATP synthase subunit D, probale (VaoD) mRNA, complete cds
# Database: plant_rna
# Fields: query frame, subject id, evalue, alignment length, % identity
# 1 hits found
1 gi|145342129|ref|XM_001416109.1| 0.0 1029 100.00

评论这张

转发至微博

阅读(2018)| 评论(0)

历史上的今天

this.p={  m:2,
              b:2,
              loftPermalink:'',
              id:'fks_081066092087083071087080085095085094087074093087083071',
              blogTitle:'BLAST+实战',
              blogAbstract:'BLAST+实战<BR\>2010-04-29 16:08<BR\>  <TABLE\>  <TBODY\>  <TR\>  <TD\><A rel=\"nofollow\" href=\"http://hi.baidu.com/lidaof/blog/item/9d29eade7f340a5895ee3724.html\"  \>http://hi.baidu.com/lidaof/blog/item/9d29eade7f340a5895ee3724.html</A\><BR\>  <P style=\"TEXT-INDENT: 2em;\"  \>项目要求，开始使用BLAST+</P\>  <P style=\"TEXT-INDENT: 2em;\"  \>我们从格式化数据库开始</P\>  <P style=\"TEXT-INDENT: 2em;\"  \>1，我先获得了一个植物rRNA的序列文件，为fasta格式，用makeblastdb格式化之，如下</P\>  <P style=\"TEXT-INDENT: 2em;\"  \>makeblastdb -in plant26.rna.fna -dbtype nucl -parse_seqids  -hash_index  -out plant_rna</P\>  <P style=\"TEXT-INDENT: 2em;\"  \><FONT color=\"#ff0000\" \>注意：在BLAST+2.2.24中 不加参数 -parse_seqids </FONT\></P\></TD\></TR\></TABLE\>',
              blogTag:'',
              blogUrl:'blog/static/1869915420106194031949',
              isPublished:1,
              istop:false,
              type:0,
              modifyTime:1317309607478,
              publishTime:1279526431949,
              permalink:'blog/static/1869915420106194031949',
              commentCount:0,
              mainCommentCount:0,
              recommendCount:0,
              bsrk:-100,
              publisherId:0,
              recomBlogHome:false,
              currentRecomBlog:false,
              attachmentsFileIds:[],
              vote:{},
              groupInfo:{},
              friendstatus:'none',
              followstatus:'unFollow',
              pubSucc:'',
              visitorProvince:'',
              visitorCity:'',
              visitorNewUser:false,
              postAddInfo:{},
              mset:'000',
              mcon:'',
              srk:-100,
              remindgoodnightblog:false,
              isBlackVisitor:false,
              isShowYodaoAd:true,
              hostIntro:'专业背景：计算机科学                           \n\n研究方向与兴趣: JavaEE-Web软件开发,\n生物信息学, 数据挖掘与机器学习, 智能信息系统                                                                \n目前工作: 基因组, 转录组, NGS高通量数据分析, 生物数据挖掘, 植物系统发育和比较进化基因组学',
              hmcon:'1',
              selfRecomBlogCount:'0',
              lofter_single:'<iframe width="140" height="560" style="overflow:hidden;" src="http://www.lofter.com/mailEntry.do?blogad=1&blog" frameBorder="0"></iframe>'
            }

{list a as x}
    {if !!x}
    <div class="iblock nbw-fce nbw-f40">
      <a class="fc03 noul" target="_blank" hidefocus="true" href="http://blog.163.com/${x.visitorName}/">
      {if x.visitorName==visitor.userName}
      <img alt="${x.visitorNickname|escape}" onerror="this.src=location.f40" class="cwd bdwa bdc0" src="${fn1(x.visitorName)}&r=${visitor.imageUpdateTime}"/>
      {else}
      <img alt="${x.visitorNickname|escape}" onerror="this.src=location.f40" class="cwd bdwa bdc0" src="${fn1(x.visitorName)}"/>
      {/if}
      </a>
      <div class="cwd vname thide">
        {if x.moveFrom=='wap'}
          <a class="noul pnt" target="_blank" href="http://blog.163.com/services/wapblog.html?frompersonalbloghome"><span title="来自网易手机博客" class="iblock wapIcon"> </span></a>
        {elseif x.moveFrom=='iphone'}
          <a class="noul pnt" target="_blank"><span title="来自iPhone客户端" class="iblock iphoneIcon"> </span></a>
        {elseif x.moveFrom=='android'}
          <a class="noul pnt" target="_blank"><span title="来自Android客户端" class="iblock androidIcon"> </span></a>
        {elseif x.moveFrom=='mobile'}
          <a class="noul pnt" target="_blank" href="http://blog.163.com/services/emsblog.html?frompersonalbloghome"><span title="来自网易短信写博" class="iblock wapIcon"> </span></a>
        {/if}
        <a class="fc03 m2a"  target="_blank" hidefocus="true" href="http://blog.163.com/${x.visitorName}/">
          ${fn(x.visitorNickname,8)|escape}
        </a>
      </div>
    </div>
    {/if}
    {/list}

<#--最新日志，群博日志--> <#--推荐日志-->

<p class="fc06">推荐过这篇日志的人：</p>
    <div>
      {list a as x}
      {if !!x}
      <div class="iblock nbw-fce nbw-f40">
        <a class="fc03 noul" target="_blank" hidefocus="true" href="http://blog.163.com/${x.recommenderName}/">
        <img alt="${x.recommenderNickname|escape}" onerror="this.src=location.f40" class="cwd bdwa bdc0" src="${fn1(x.recommenderName)}"/>
        </a>
        <div class="cwd thide">
          <a class="fc03 m2a" target="_blank" hidefocus="true" href="http://blog.163.com/${x.recommenderName}/">
            ${fn(x.recommenderNickname,6)|escape}
          </a>
        </div>
      </div>
      {/if}
      {/list}
    </div>
    {if !!b&&b.length>0}
    <p  class="fc06">他们还推荐了：</p>
    <ul>
    {list b as y}
      {if !!y}
        <li class="rrb"><span class="iblock">·</span><a class="fc03 m2a" target="_blank" href="http://blog.163.com/${y.recommendBlogPermalink}/?from=blog/static/1869915420106194031949">${y.recommendBlogTitle|escape}</a></li>
      {/if}
    {/list}
    </ul>
    {/if}

<#--引用记录--> <#--博主推荐--> <#--随机阅读--> <#--首页推荐--> <#--历史上的今天--> <#--被推荐日志--> <#--上一篇，下一篇--> <#-- 热度 -->

{list a as x}
    {if !!x}
    <div class="hotItem iblock nbw-fce nbw-f40">
      <a class="fc03 noul" target="_blank" hidefocus="true" href="http://blog.163.com/${x.publisherUsername}/">
      {if x.publisherUsername==visitor.userName}
      <img alt="${x.publisherNickname|escape}" onerror="this.src=location.f40" class="cwd bdwa bdc0" src="${fn1(x.publisherUsername)}&r=${visitor.imageUpdateTime}"/>
      {else}
      <img alt="${x.publisherNickname|escape}" onerror="this.src=location.f40" class="cwd bdwa bdc0" src="${fn1(x.publisherUsername)}"/>
      {/if}
      </a>
      <div class="cwd vname thide">
        <a class="fc03 m2a"  target="_blank" hidefocus="true" href="http://blog.163.com/${x.publisherUsername}/">
          ${fn(x.publisherNickname,8)|escape}
        </a>
      </div>
      <a class="f-myLikeIcons hottype {if x.type==1} js-liketype{elseif x.type==2} js-reblogtype{elseif x.type==3} js-sharetype{else}{/if}" target="_blank" hidefocus="true" href="http://blog.163.com/${x.publisherUsername}/"> </a>
    </div>
    {/if}
    {/list}

<#-- 网易新闻广告 --> <#--右边模块结构--> <#--评论模块结构--> <#--引用模块结构--> <#--博主发起的投票-->

页脚

我的照片书 - 手机博客 - 下载LOFTER APP - 订阅此博客

云之南

导航

日志

BLAST+实战

历史上的今天

最近读者

热度

评论

页脚