View source for Read split.pl

[[Category:Bioinformatics ]]
<source lang=perl>
#! /usr/bin/perl -w
# *************************************************************
# This Programme does the following in sequence
# -> Parses a multi fasta
# -> Count characters and splits as (1) 0 (2) <=80 (3) >80
# -> Capture the results in files in 3 catogories
# -> Generates the statistics
# -> Generates Tab delimited file with ID and size
# Author : Rajkumar (itc@rajkumar.in)
# Release: MAR 2010
# *************************************************************
 
# Libraries
	use strict;
	use Bio::Perl;
	use Bio::SeqIO;
 
# Variables
	my $file = "set9"; # change this for each file	
	my $path = "/home/raj/bio/user/raja";
	my $multifasta = "$path/data/$file.txt";
	my $file1      = "$file.reads_0bp.txt";
	my $file2      = "$file.reads_1-80bp.txt";
	my $file3      = "$file.reads_80+bp.txt";
	my $statfile   = "$file.stat.txt";
	my $idfile     = "$file.id_size.txt";
 
# Getting the grip of files
 	if (-e "$file1") {`rm  $file1`;}
	if (-e "$file2") {`rm  $file2`;}
	if (-e "$file3") {`rm  $file3`;}
 	if (-e "$idfile") {`rm  $idfile`;}
 
# Let's Hack
	my $in  = Bio::SeqIO->new(-file => "$multifasta", -format => 'Fasta');
	my ($identifier, $sequence) = $_;
 	while ( my $seq = $in->next_seq() ) {
		$identifier = $seq->id;
  		$sequence = $seq->seq;
 
 
		my $count_seq = length($sequence);
		my $fasta = ">".$identifier."_Size:".$count_seq."bp\n".$sequence."\n";
		my $id_size = "$identifier\t$count_seq\n";
 
		if ($count_seq eq 0){ 
			open (FH, ">>".$file1);
			print FH $fasta ;
			close FH ;
 
			open (FH, ">>".$idfile);
			print FH $id_size ;
			close FH ;
 
		}
 
		if ( $count_seq < 81 && $count_seq > 0){ 
			open (FH, ">>".$file2);
			print FH $fasta ;
			close FH ;
 
			open (FH, ">>".$idfile);
			print FH $id_size ;
			close FH ;
		}
 
		if ( $count_seq > 80){ 
			open (FH, ">>".$file3);
			print FH $fasta ;
			close FH ;
 
 
			open (FH, ">>".$idfile);
			print FH $id_size ;
			close FH ;
		}		
	}
# statistics
	my $stat1 = `grep -c ">" $file1`; chomp $stat1;
	my $stat2 = `grep -c ">" $file2`; chomp $stat2;
	my $stat3 = `grep -c ">" $file3`; chomp $stat3;		
 
 
	my $stat = "Statistics\n==========\n$file1    $stat1 sequences\n$file2 $stat2 sequences\n$file3  $stat3 sequences\n";
	open (FH, ">".$statfile);
	print FH $stat ;
	close FH ;		
 
# END

</source>