Hi Tomas, excellent question, more below... On Mon, 2010-04-05 at 21:36 -0300, tomas duffy wrote: > Hello, > I've aligned a few chlamydial genomes in GenBank format. When I use > the Export Orthologs function I get the two files I'm supposed to get, > the one with the positions of the orthologs in each genome, and the > one that also includes the sequences. > Unfortunately, the sequences do not correspond with the position > stated in the header. For example, this is the sequences of one of the > first group of orthologs: > >0:266-1267:TC_0001 + > tactgttctgaaggagatgagtttgctggagcttctttagccttttcacgacgataagctagatactaagctctttgatt
<snip> > acccaaaatcttaccgatcgatctctgcgagaaggaaaaattcctttttcaagaaact > > This gene (266-1267:TC_0001) when I look it up on EMBL is different > from the sequences above. > > What is it that I'm doing wrong? Am I miss interpreting something? I think you're using Mauve correctly, and although I can't say with certainty unless you tell me exactly which other genomes you were aligning, I believe the output was also generated correctly by Mauve. As you observed, the sequence coordinates listed do not correspond to the sequence itself. However, in some cases that is an expected result for the ortholog export functionality (but NOT for normal XMFA genome alignments). It's true that this behavior is counter-intuitive, however several design decisions led to the behavior. First, when Mauve generates an ortholog alignment file, it includes everything between the left-most and right-most alignment column covered by the orthologous features. The feature (e.g. a gene) may be annotated in one organism to contain a much larger range of sequence than it does in another genome. Rather than trim the alignment to the smallest range, Mauve includes the largest orthologous range possible. There are various reasons we chose to do this, but it may not be desirable behavior in all cases, so perhaps a good software improvement would be more flexibility here. But I digress. The second issue leading to the behavior you observed is how Mauve names the homologous features in the output file. Mauve defaults to using the locus_tag qualifier from the genbank file, however, for some feature classes and in some (degenerate) genbank format files, /locus_tag qualifiers may not be present. In fact, there may not be any reasonable name for the feature in any qualifier, but only feature type and the sequence range. I wish such GenBank files didn't exist, but they do. (and I wish more people would use XML, and I wish...). To handle such degenerate cases, the features are "named" in part by the coordinate range they were annotated with in the genbank file. Thus, the feature you pasted above is "named" 266-1267:TC_0001. In your case the genbank file also has a locus_tag, so TC_0001 is appended to the name. Admittedly, the formatting here looks like it should be specifying the sequence range given immediately after the header line, and this is confusing with the XMFA format. I think a good revision of the export orthologs function would be to give the full sequence range, then give the annotated feature range (and label it as such, in XML, hey!), then the feature name or locus tag. Anyhow, I think that's what happened in your output file, and that what you observed isn't a bug, it's just counter-intuitive behavior. Confusing FastA header lines are the bane of bioinformatics. -Aaron ------------------------------------------------------------------------------ Download Intel® Parallel Studio Eval Try the new software tools for yourself. Speed compiling, find bugs proactively, and fine-tune applications for parallel performance. See why Intel Parallel Studio got high marks during beta. http://p.sf.net/sfu/intel-sw-dev _______________________________________________ Mauve-users mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/mauve-users
