problam in nesting loop

R

Rita

Hi Members
I am making a program and I ran into a bit of difficulty with it. my
program is Bioinformatics related.and i have two files of squences.i
took one seq from first file and second seq from second file.and i did
Blast to both sequence then i want Gap Sequence in result.and my
program is working for 2 sequences but i want to put in loop to whole
program .so it will take for 1 inputf irst seq from first file then 2
input first sequence from second file then blast it and find a gap
sequence then again it will take second sequence from first file and
second seq from second file.........
so the problam is that the gap sequence is always different some
sequence has 1 gap some has 2 gap and some doesn't has gap .so i want
to put this thing is in loop and store all output sequences in a file.i
am trying it but i cannot figur it out please can you give a look to my
programm-
use warnings;
use strict;
use Bio::SeqIO;
use Bio::SearchIO;

# Get the name of the text file from the user that has the
# seqid, source seq :
print "Please type the filename of the seqid and source seq : ";
my $filename = <STDIN>;

# Remove the newline from the filename
chomp $filename;

# either open the file, or exit from the program
open FILE,$filename ||die "Cannot open file \"$filename\"\n\n $!";
<FILE>;
chomp(my $line = <FILE>);

my ($seqid, $source_seq) = (split /\t/, $line)[1,4];
print $seqid,"\t";
print $source_seq,"\n";

#write First sequence in Fasta Format and store in file name tmp:
my $seqio = Bio::SeqIO ->new(-file =>'>tmp',-format =>'fasta');
$source_seq = Bio::Seq-> new(-seq => $source_seq, -disaply_id
=>'test1');
$seqio ->write_seq($source_seq);

my$filename1="$seqid";
open(FILE1,"C:/Ritu/FastaSeqs/$filename1.CONSENS.fasta")||die "Cannot
open file

\"$filename1\"\n\n $!";
<FILE1>;
chomp(my $sts_seq=<FILE1>);


#Blast to both sequences and store output in file name temp.blast:
system("c:/downloads/bin/bl2seq -i tmp -j
C:/Ritu/FastaSeqs/$filename1.CONSENS.fasta
-p blastn -o temp.blast ");

my $in = Bio::SearchIO ->new(-format => 'blast',
-file => 'temp.blast');
while( my $result = $in->next_result ) {
while( my $hit = $result->next_hit ) {

#Print Starting and Ending Point of 1 Hsps-

my $hsp1 =$hit->next_hsp;
my$start1=$hsp1->hit->start;
my $end1=$hsp1->hit->end;
print "start First Hsp ",$start1,"\n";
print "End First Hsp ",$end1,"\n\n";

#Print Starting and Ending Point of 2 Hsps-
my $hsp2 =$hit->next_hsp;
my$start2=$hsp2->query->start;
my $end2=$hsp2->query->end;
print "start Second Hsp ",$start2,"\n";
print "End Second Hsp ",$end2,"\n\n";

#Get Gap Sequence :
print "Gap Sequence- ",$sts_seq->subseq($end1+1,$start2-1);
}
}


#Close the file
close FILE;
close FILE1;
thanks
 
R

Rita

Hi
Your Programm is so good and it's working for two files. But my problam
is I don't want to give two input filename. I just want to give one
filename Data.txt which has 5 fields like
Submitter Name Sequence ID GB Accession # Comments
Source sequence
And my Data is
Ik-Young 12945 AW348161 good sequence
TGCAHJKLSURIOEPWIFJDIJDOIOKDSCMKCVCNCNCJLJDOFUOKFDLJLXVJNCLVJDLFDFJOFJUDOPPWPLDMGCCGCAACJLDJKLKLKJAOJODLOFOOQPQKDLSFJXNVNZMLPIRPEJIOFJFKLFVJODIPIEKKLCMMVLJOFPKAKLVMMLCVJFOLJOFPFKLDFKOEWIEOOGHHJGLFWPPEPFFFJPJGPGITIPEOPQFKKPROPEGOJVBJMKLNLCMLKASLKOPIIPRKVKMMBOPOEPLCDCVMLFBJTOJGTORIROEIRPEOIPWEIKEPRIPRPIRPWIRPEIPEIRWPIEPRIPE
..
..
..
..
..
..50 Data like this


I want to extract 2 fields in it. First one is Sequence ID and Second
one is Source Sequence.So i Use this code in programm but this code
reading only First Sequence

# Get the name of the text file from the user that has the seqid,
source seq :
print "Please type the filename of the seqid and source seq : ";
my $filename = <STDIN>;

# Remove the newline from the filename
chomp $filename;

# either open the file, or exit from the program
open FILE,$filename ||die "Cannot open file \"$filename\"\n\n $!";
<FILE>;
chomp(my $line = <FILE>);
my ($seqid, $source_seq) = (split /\t/, $line)[1,4];
print $seqid,"\t";
print $source_seq,"\n";

Now i want to write source sequence in fasta format in a temprary file
So I use

#write First sequence in Fasta Format and store in file name tmp:
my $seqio = Bio::SeqIO ->new(-file =>'>tmp',-format =>'fasta');
$source_seq = Bio::Seq-> new(-seq => $source_seq, -disaply_id
=>'test1');
$seqio ->write_seq($source_seq);

Then I want to give Second file name which is same as Sequence ID :
So I give code

my$filename1="$seqid";
open(FILE1,"C:/Ritu/FastaSeqs/$filename1.CONSENS.fasta")||die "Cannot
open file
\"$filename1\"\n\n $!";

This file has only one sequence . I have 1 folder with 50 files which
name is same as Sequence ID.

Then I Blast it
#Blast to both sequences and store output in file name temp.blast:
system("c:/downloads/bin/bl2seq -i tmp -j
C:/Ritu/FastaSeqs/$filename1.CONSENS.fasta -p blastn -o temp.blast ");

After Blast i Want toResult which has Nomatch Sequence

my $in = Bio::SearchIO ->new(-format => 'blast',
-file => 'temp.blast');
while( my $result = $in->next_result ) {
while( my $hit = $result->next_hit ) {

#Print Starting and Ending Point of 1 Hsps-

my $hsp1 =$hit->next_hsp;
my$start1=$hsp1->hit->start;
my $end1=$hsp1->hit->end;
print "start First Hsp ",$start1,"\n";
print "End First Hsp ",$end1,"\n\n";

#Print Starting and Ending Point of 2 Hsps-
my $hsp2 =$hit->next_hsp;
my$start2=$hsp2->query->start;
my $end2=$hsp2->query->end;
print "start Second Hsp ",$start2,"\n";
print "End Second Hsp ",$end2,"\n\n";

#Get Gap Sequence :
print "Gap Sequence- ",$sts_seq->subseq($end1+1,$start2-1);
}
}


#Close the file
close(FILE) or warn("Error closing $file1: $!");
close(FILE1) or warn("Error closing $file2: $!");
This Programm is giving me Result what I want But It is Good only for
two Sequence and I want to make for all 50 Sequences .
and for that first i have to use loop in whole program and second thing
is i have to put loop for get Gap Sequence too .Because this Sequence
has only one nonmatch or Gap Sequence but may be another 49 sequence
have 0 or 2,3.... Nonmatch or Gap Sequence.
I am trying but i couldn't figureout that.
still can i use your script in this programm?
Rita
 

Ask a Question

Want to reply to this thread or ask your own question?

You'll need to choose a username for the site, which only take a couple of moments. After that, you can post your question and our members will help you out.

Ask a Question

Members online

Forum statistics

Threads
473,969
Messages
2,570,161
Members
46,705
Latest member
Stefkari24

Latest Threads

Top