#! /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