Make data file ssr.pl

From RAJ INFO
Jump to navigation Jump to search
#! /usr/bin/perl -w


# ****************************************************
# This Programme does the following in sequence
# -> Open a multi fasta files
# -> combine camal and grandis sequeces
# -> Count the number of bases in each fasta
# -> Filter off the seq that have less than 350bp | more than 10% of N
# -> Output as fasta
# Author: Rajkumar (itc@rajkumar.in)

# ****************************************************

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

# Variables
	my $path         = "/home/raj/bio";
	my $datadir      = "$path/euca/data/camal";
	my $infile       = "$datadir/camal_all.fa";
	my $outdir       = "$path/euca/ssr/glo-gra-ssr";
	my $grand_infile = "$path/euca/data/grandis/scaffold/ssr/scaf_ssr.fa";
	my $temp         = "$path/euca/ssr/temp";
	my $outfile      = "camal_grand.fa";
	

# Making the res dir clear
	
	if (-e "$outdir/$outfile") {`rm $outdir/$outfile`; }
	
	
# Making a common inputfile for SSR mining
	
	`cat $infile $grand_infile >$temp`;
	
	
# Further Triming (Removing the N component)

	# Clear the variables
	my ($in, $identifier, $sequence) = ();
	
	# Opening
	$in  = Bio::SeqIO->new(-file => "$temp", -format => 'Fasta');

	
	while ( my $seq = $in->next_seq() ) {
		$identifier = $seq->id;
  		$sequence = $seq->seq;
  		
  		my $length = ($sequence =~ tr/[A-Z]|[a-z]//);
		my $count_N = ($sequence =~ tr/N//);
  		my $math = ($count_N * 100) / $length;

		# conditions to remove seq with <350 and > 10% of N
		if ( $length < 350 ) { next; }
		if ( $math > 10 ) { next; }
		
		# Removing N at beging and end
		$sequence =~ s/^N+//g;
		$sequence =~ s/N+$//g;
		
		# Trimming and naming the identifier
		$identifier =~ s/\|/_/g;
		if ($identifier !~ m/^GRA/) {$identifier =~ s/(.+)/CAMAL_$1/g;} 
		$identifier =~ s/(.+)_$/$1/g;

		
		
		#print "$identifier\t$length\t$count_N\n$sequence\n";
		
		# Write the results
		my $fasta  = ">$identifier\n$sequence\n";
		open (FH, ">>".$outdir."/".$outfile);
		print FH $fasta;
		close FH ;

  	}	
	
	`rm $temp`;
# END