problam in nesting loop

Discussion in 'Perl Misc' started by Rita, Nov 18, 2005.

  1. Rita

    Rita Guest

    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
     
    Rita, Nov 18, 2005
    #1
    1. Advertising

  2. Rita

    Rita Guest

    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
     
    Rita, Nov 21, 2005
    #2
    1. Advertising

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

It takes just 2 minutes to sign up (and it's free!). Just click the sign up button to choose a username and then you can ask your own questions on the forum.
Similar Threads
  1. Hardrock

    Dynamic loop nesting

    Hardrock, Sep 2, 2004, in forum: C Programming
    Replies:
    8
    Views:
    436
    Hardrock
    Sep 13, 2004
  2. Trans
    Replies:
    10
    Views:
    327
    Sean O'Halpin
    Sep 16, 2005
  3. Rita

    problam in subroutine?

    Rita, Dec 5, 2005, in forum: Perl Misc
    Replies:
    2
    Views:
    109
  4. Rita

    Problam in extract data

    Rita, Dec 14, 2005, in forum: Perl Misc
    Replies:
    18
    Views:
    182
    Tad McClellan
    Dec 16, 2005
  5. Isaac Won
    Replies:
    9
    Views:
    457
    Ulrich Eckhardt
    Mar 4, 2013
Loading...

Share This Page