Categorized | Biocompute

Tags | ,

Bioperl:解析blast结果

Posted on 27 五月 2009 by 柳城 ,阅读 563

下面是一个BioPerl小脚本, 用于解析Blast的结果.

#!/usr/bin/perl -w
use strict;
use Bio::SearchIO;

my $filename = 'MET30-vs-swissprot.3.WUBLASTP';
my $input = Bio::SearchIO->new(-format => 'blast',
          -file   => $filename);
#GOAL:
# Print out the query, hit names and percent identity, percent query aligned
while( my $r = $input->next_result ) {
  my $qname = $r->query_name;
  my $qlen  = $r->query_length;
  while( my $hit = $r->next_hit ) {
    my $hname = $hit->name;
    while( my $hsp = $hit->next_hsp ) {
      my $percent_id = $hsp->frac_identical * 100;
      my $query_aln_len = $hsp->length('query');
      # another way to do this
      # my $query_aln_len = $hsp->query->length;
      print join("\t", $qname, $hname, sprintf("%.1f",$percent_id),
               sprintf("%.2f",100 * ( $query_aln_len / $qlen))), "\n";
    }
  }
}

转载请注明 : 来源于 Bioperl:解析blast结果 | 柳城

赞助商

1条评论 于 “Bioperl:解析blast结果”

  1. Ying Ying Says:

    你好,我现在想通过blast结果判断基因的拷贝数多少。
    请问在score和e-value上应该怎样去cutoff值作为单拷贝和多拷贝的界限。谢谢~

    [回复]

Leave a Reply

广告招租

[强] [握手] [可爱] [ok] [呲牙] :) [偷笑] [流泪] [疑问] [亲亲] [擦汗] [得意] [衰] [可怜] [抱拳] [坏笑] more »

无觅相关文章插件,快速提升流量