Perl Print Nextline From Within While Loop?

Question

[User 3739]

I am trying to query a large FASTA file in which I start with a list of genes, search for the gene in the FASTA, and then I want to print out the header line (>gene) and the protein sequence from the subsequent line.

Due to the large size of the FASTA file, I am using a while loop. I can print the line that contains the string I am searching for, but I cannot figure out how to print the next line that contains the protein sequence.

Can someone please tell me the easiest way to do that?

Thanks.

[biobot 0.0.77.a.1099]

See the answers to this question. Do you realise that the next line after the header may not contain the entire sequence?

[Istvan Albert]

Hi Keith, I think your comment is appropriate as an answer as well. You might want to consider adding it.

[biobot 0.0.77.a.1099]

Duly added, with a pointer to another relevant answer too.

Answer 1

[biobot 0.0.77.a.1099]

See the answers to this this question on Fasta parsing in Perl and the answers to this question on extracting sequences from very large Fasta files.

Do you realise that the next line after the header may not contain the entire sequence? If you are unsure of the details involved in parsing a file format, it may be better to use a library such as Bioperl.

Answer 2

[Keeney]

This may be an inefficient way of doing it because I'm still fairly new to it, but one idea is to loop through it once, building an array by splitting up the gene inputs by newlines, then loop through your array looking for sequence matches. If it finds a match, just print the current position in the array +1.

So, something along the lines of this (the syntax is wrong):

@array = split(/\n/, ListOfGenes.fasta);

$i = 0;

# Now loop through the big FASTA file.
while ( <IN> ) {
  if ( $_ == $array[$i] ) {
    print OUT $array[$i+1]
  }
$i++;
}

It unnecessarily tries to match with everything in the array (which could be avoided by advancing $i by 2 or 3 or whatever's appropriate for your list, rather than by 1). Like I said, I'm still pretty new to it, so I don't know if that works, but it's quick to try!

That's my two cents, sorry if it's completely useless information!

[biobot 0.0.77.a.1099]

You realise that there are flaws in your proposal - good, that's the first step to improvement. Splitting on \n doesn't break the Fasta file on record boundaries. Assuming <IN> is reading in gene names, you're stepping both the array index and gene name at the same time when you want a nested loop (and what you really want is a hash lookup). Your $_ test should be eq for strings, not ==. The test will always fail because Fasta headers contain a leading >. The next element in the array will not always contain the entire sequence. The Fasta file is large and probably won't fit in memory.

[Keeney]

Hmm, so I may have misunderstood the input files; in this proposal, I was assuming the [?] was the large FASTA file with no newline characters, and ListOfGemes.fasta was an alternating list of genes and their protein sequences, separated by newlines (I assumed this because in the question he asks for the protein sequence from the "next line"). Maybe that's not what was meant...? I think you make a great point that a hash might be better, though. I'm still not too experienced with them so I didn't know how to propose that! And you're right, it should have been 'eq' rather than '=='.

[Keeney]

Maybe a pattern match in between > and whatever else (methionine, for example) would be better than splitting by newline, then?

Like (again, syntax is probably wrong):

if ( ListOfGenes.fasta =~ m/>M/g ) { ...

Answer 3

[Prakki Rama]

You can try this:

open FH,"filename.txt";

while(<FH>)
 {
    $line1=$_;
    $line2=&lt;FH&gt;;
    print "Line1: $line1";
    print "Line2: $line2";
}

close(FH);