Question

求生物信息学相关的 应用Perl 语言的 简单案例

Answer 1

这是我刚开始学习perl的时候写的一个应用bioperl做本地blast的脚本，希望有点帮助

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

my $filein=shift;

open OUT,">test.out" or die "$!";
my @para=(program=>'blastn',database=>"hg18/chr1.fa",expect=>1e-10);

my $blast=Bio::Tools::Run::StandAloneBlast->new(@para);

my $str = Bio::SeqIO->new(-file=>$filein,-format =>'fasta');
while (my $input = $str->next_seq())
{
my $id=$input->display_id();
my $len=$input->length();
my $report=$blast->blastall($input);
#print "$report\n";
while(my $result=$report->next_result)
{
while(my $hit = $result->next_hit)
{
while(my $hsp = $hit->next_hsp)
{
next if($hsp->percent_identity<85);
next if($hsp->length('total')<300);
#next if($hsp->hit->strand eq $hsp->query->strand);
print OUT $id,"\t",$hit->name,"\t",$hsp->length('total'),"\t",$hsp->hit->strand,"/",$hsp->query->strand,"\t",$hsp->hit->start,"\t",$hsp->hit->end,"\t",$hsp->query->start,"\t",$hsp->query->end,"\t",$hsp->percent_identity,"\t",$hsp->evalue,"\n";
}
}
}
}

close OUT;