BES-NR blast.pl

From RAJ INFO
Revision as of 11:31, 16 July 2010 by Raj (talk | contribs)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Jump to navigation Jump to search
#! /usr/bin/perl -w
# *************************************************************
# This Programme does the following in sequence
# -> Parses a multi fasta
# -> Removes the sequence if they have <50bp and of 'N'
# -> Capured as individual files.
# Author : Rajkumar (itc@rajkumar.in)
# Release: JUL 2009
# *************************************************************

# Libraries
	use strict;
	use Bio::Perl;
	use Bio::SeqIO;

# Variables
	my $path        = "/home/raj/bio";
	my $path1       = "$path/user/raja";
	my $path2       = "$path1/blast/BES-NR";
	my $multifasta  = "$path2/BES.fa";
	my $resdir      = "$path2/multifasta";
	
# Getting the Grip of folders
	if (-e $resdir)   {`rm -fr $resdir`;}
	if (! -e $resdir) {`mkdir $resdir`;}


# making several Fasta files of GDNA
	my $in  = Bio::SeqIO->new(-file => "$multifasta", -format => 'Fasta');
	my $identifier;
	my $sequence;
	
	while ( my $seq = $in->next_seq() ) {
		$identifier = $seq->id;
  		$sequence = $seq->seq;
  		
  		#Trim the identifier
  		$identifier =~ s/_Size:[0-9]+bp//g; 
		
		my ($fasta) = $_;
		$fasta = ">$identifier\n$sequence\n";

		my $count_seq = length($sequence);
		my $count_N = ($sequence =~ tr/N//);
  		my $math = $count_seq * 0.2;
  		
		$sequence =~ s/-/N/;
		if ($count_seq < 50){ next; }
		if ($count_N > $math){ next; }
				 
		open (FASTA, ">".$resdir."/".$identifier.".fa");
		print FASTA $fasta;
		close FASTA;		
	
		# print "$identifier\t$count_seq\n";
	
	}
# END