You can find the longest ORFs for a set of DNA sequences with HMMER2GO by using a simple command (more info on the wiki):

hmmer2go getorf -i genes.fasta -o genes_orfs.fasta

This would be the ideal use case if you have ESTs or gene sequences, and it should give you the same result as filtering the results of getorf in a separate script (though I haven't tested the other approaches). You'll want to rethink these approaches if your contigs are actually genomic contigs where you'd expect many genes per contig.

Hi,

As an exercise, I tried to do it with the GATB C++ library, and it can be done with the following code:

#include <gatb/gatb_core.hpp>

#include <map>
#include <list>
#include <set>
#include <string>

using namespace std;

/********************************************************************************/
struct Entry
{
    Entry (int index=0, int length=0) : index(index), length(length) {}
    int    index;
    int    length;
    static bool compare (const Entry& e1, const Entry& e2)  { return e1.length > e2.length; }
};

/********************************************************************************/

int main (int argc, char* argv[])
{
    try
    {
        OptionsParser parser;
        parser.push_front (new OptionOneParam (STR_URI_INPUT,  "input bank"));
        parser.push_front (new OptionOneParam (STR_URI_OUTPUT, "output bank"));

        /** We read command line options. */
        IProperties* props = parser.parse (argc, argv);

        /** We get a handle on the input bank. */
        BankFasta intputBank (props->getStr(STR_URI_INPUT));

        /** We get a handle on the output bank. */
        BankFasta outputBank (props->getStr(STR_URI_OUTPUT));

        map < string,list<Entry> > seqMap;
        set <int> indexes;

        /** We read the input bank sequences. */
        Iterator<Sequence>* it = intputBank.iterator();  LOCAL (it);
        for (it->first(); !it->isDone(); it->next())
        {
            /** Shortcut. */
            Sequence& s = it->item();

            TokenizerIterator token (s.getComment().c_str(), " .");
            token.first();  string k1 = *token;

            seqMap[k1].push_back (Entry (s.getIndex(), s.getDataSize()));
        }

        for (map < string,list<Entry> >::iterator it = seqMap.begin(); it != seqMap.end(); it++)
        {
            /** We sort the sequences by decreasing sequence lengths. */
            it->second.sort(Entry::compare);

            /** We memorize the index of the biggest sequence. */
            indexes.insert (it->second.front().index);
        }

        /** We iterate the input bank sequences again. */
        for (it->first(); !it->isDone(); it->next())
        {
            /** We keep sequences for only our set of computed indexes. */
            if (indexes.find(it->item().getIndex()) != indexes.end())  {  outputBank.insert (it->item());  }
        }

        outputBank.flush ();
    }
    catch (Exception& e)
    {
        std::cout << "EXCEPTION: " << e.getMessage() << std::endl;
        return EXIT_FAILURE;
    }

    return EXIT_SUCCESS;
}

You can get the whole material for building the binary here. Once unarchived, you have a readme file explaining how to build the binary (mkdir build; cd build; cmake ..; make)

You can then run the command like this: orf -in input.fa -out output.fa

I got the correct output bank with the tiny example you gave.

If you have issues to build it, I can give you my binary.

Hi, this perl script can be used for finding out the longest ORF. You can run the script inside the folder contating the EMBOSS generated ORF files (each having many ORFs generated from different frames). This requires Bioperl and Text::Wrap module for formatting, which you can replace by some REGEX based method. NOTE that this script does NOT output ALL longest ORFs for a given DNA sequence. There can be more than one longest sequences having same lengths.

use strict;
use warnings;
use Bio::SeqIO;
use Text::Wrap;
$Text::Wrap::columns = 61;
my @orf_files=<*.orf>;
foreach my $input_orf(@orf_files){
my $in  = Bio::SeqIO->new(-file =>  $input_orf , '-format' => 'Fasta');
my $longest_seq="";
my $longest_id="";
my $longest_desc="";
while ( my $seq = $in->next_seq ){
               my $id=$seq->id();
               my $seq_str= uc ($seq->seq() );
               if(length($seq_str) > length($longest_seq) ){
                $longest_id=$id;
                $longest_desc=$seq->desc();
                $longest_seq=$seq_str;
                }


        }
print  ">",$longest_id," ",$longest_desc,"\n",wrap('', '', $longest_seq),"\n";            
}

Transdecoder does what you're asking!