Blast+ xml格式解读

本地BLAST比对后,如果使用outfmt 5参数的话,会产生一个xml格式的文件,里面的比对信息不像tabular(outfmt 6)那样简显,但是对比信息却很完整。简单列举一些常用的信息。

  1. 每条序列的所有对比信息是以标记<Iteration>开始,</Iteration>结束,之间框在一起的则是一条序列的所有比对结果信息
    1. 后面接的数字表示序列的编号,1代表第一条序列
    2. 序列的header
    3. 序列的长度
  2. 每条序列的比对结果是以标记<Iteration_hits>开始,</Iteration_hits>结束,框架下面还有<Hit>表示序列的每个比对结果(是指一条序列比对上对应的一条满足阈值的序列),而<Hsp>则表示每个比对结果中的某一块的比对结果(比如一条序列上有好几处跟目标序列比对上了)
    1. <Hsp_num>比对上的序列的编号
    2. 比对结果的打分
    3. <Hsp_evalue>比对结果的evalue值
    4. 目标序列比对上的开始位置
    5. 目标序列比对上的结束位置
    6. 比对上序列的开始位置
    7. 比对上序列的结束位置
    8. 目标序列翻译的起始位置
    9. 比对上序列翻译的起始位置
    10. 比对上序列那部分的长度
    11. <Hsp_qseq>目标序列比对上的那部分序列信息
    12. <Hsp_hseq>比对上序列的那部分序列信息
  3. 简单的归纳就这样了,具体可翻看文档
  4. 知道这些标记的含义,我们就可以用perl将自己需要的信息提取出来,然后输出成表格。
#!/usr/bin/perl -w
use strict;

open my $fh, "blast.xml" or die "Cannot open the blast file\n";
local $/ = "<\/Iteration>\n";

while(<$fh>){
    my $query = $1 if ($_ =~ /<Iteration_query-def>(\S+)<\/Iteration_query-def>/);
    my $query_len = $1 if ($_ =~ /<Iteration_query-len>(\S+)<\/Iteration_query-len>/);
    while ($_ =~ /<Hit>(.*?)<\/Hit>/sg){
        chomp(my $hit = $1);
my $id = $1 if ($hit =~ /<Hit_id>(\S+)<\/Hit_id>/);
my $def = $1 if ($hit =~ /<Hit_def>(.*?)<\/Hit_def>/);
my $num = 0;
while ($hit =~ /<Hsp>(.*?)<\/Hsp>/sg){
    $num++;
    my $hsp = $1;
    my $query_from = $1 if ($hsp =~ /<Hsp_query-from>(\S+)<\/Hsp_query-from>/);
    my $query_to = $1 if ($hsp =~ /<Hsp_query-to>(\S+)<\/Hsp_query-to>/);
    my $hit_frame = $1 if ($hsp =~ /<Hsp_query-frame>(\S+)<\/Hsp_query-frame>/);
    my $qseq = $1 if ($hsp =~ /<Hsp_qseq>(\S+)<\/Hsp_qseq>/);
    my $align_len = $1 if ($hsp =~ /<Hsp_align-len>(\S+)<\/Hsp_align-len>/);
    print ">$query;orf$num len=$align_len frame:$hit_frame start:$query_from end:$query_to $id $def\n";
        }
    }
}
close $fh;

最后的输出打印出来的结果是临时的,可以随意修改

本文出自于http://www.bioinfo-scrounger.com转载请注明出处

文章目录
|