SIGI-HMM crash with some datasets

[innovate-invent]

SIGI-HMM inconsistently crashes with some datasets.
#89 #134

[innovate-invent]

@klgray25
Colombo v3.8 did not fail for any of the data:
15584_genome.gbk 15597_genome.gbk 15598_genome.gbk 15599_genome.gbk 15600_genome.gbk 15602_genome.gbk

[innovate-invent]

It appears 15584_genome.gbk will crash SigiHMM. The catch is the tool being used to convert from gbk to embl. Bioperl generates slightly different (possibly older) embl output than Biopython. IslandCompare uses Biopython for the conversion while my test originally used Bioperl. I can now consistently reproduce this for SigiHMM version 3.8.
SigiHMM 4.0 also can't parse Biopython output but emits rather than crashing:

This line could not be parsed:                 seq:caa)

and discards a large amount of islands. This line is from 15584_genome.embl from Biopython:

FT                   /anticodon=(pos:complement(1123552..1123554),aa:Leu,
FT                   seq:caa)

The matching line from 15584_genome.embl from Bioperl:

FT                   /anticodon="(pos:complement(1123552..1123554),aa:Leu,seq:ca
FT                   a)"

The difference between biopython and bioperl is that bioperl quotes anticodons if they wrap.

Moving forward we have three options:

• Add a perl dependency specifically for SigiHMM (yuck)
• Discard the offending anticodon annotations and keep biopython
• Change biopython to line wrap the same as bioperl (this may actually be a bug in biopython)

I'll need @fionabrinkman and @cbertell to weigh in on this.

[fionabrinkman]

Yes, this crashing issue is known. Claire, can you confirm which SIGIhmm version you were wanting as there seems to be some uncertainty on this? I do think it’s confusing to have Islandviewer and islandcompare make different predictions but may have forgotten conversation about this and plans to change islandcompare too.

…

On Jun 5, 2019, at 4:45 PM, Nolan Woods ***@***.***> wrote: It appears 15584_genome.gbk will crash SigiHMM. The catch is the tool being used to convert from gbk to embl. Bioperl generates slightly different (possibly older) embl output than Biopython. IslandCompare uses Biopython for the conversion while my test originally used Bioperl. I can now consistently reproduce this for SigiHMM version 3.8. SigiHMM 4.0 also cant parse Biopython output but emits rather than crashing: This line could not be parsed: seq:caa) and discards a large amount of islands. This line is from 15584_genome.embl from Biopython: FT /anticodon=(pos:complement(1123552..1123554),aa:Leu, FT seq:caa) The matching line from 15584_genome.embl from Bioperl: FT /anticodon=(pos:complement(1123552..1123554),aa:Leu,seq:caa FT ) Moving forward we have two options: Add a perl dependency specifically for SigiHMM (yuck) Discard the offending anticodon annotations and keep biopython The difference between biopython and bioperl seem to be the criteria that they use when determining where to line wrap. Bioperl could potentially have this issue if the anticodon annotations are longer in other datasets. I will contact Genbank as I believe the anticodon annotations are actually not conforming to their own standard. I'll need @fionabrinkman and @cbertell to weigh in on this. — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub, or mute the thread.

[fionabrinkman]

Kristen, I should really have said “Claire and Kristen” below, not just Claire!

…

On Jun 6, 2019, at 12:51 AM, FionaBrinkman ***@***.***> wrote: Yes, this crashing issue is known. Claire, can you confirm which SIGIhmm version you were wanting as there seems to be some uncertainty on this? I do think it’s confusing to have Islandviewer and islandcompare make different predictions but may have forgotten conversation about this and plans to change islandcompare too. > On Jun 5, 2019, at 4:45 PM, Nolan Woods ***@***.***> wrote: > > It appears 15584_genome.gbk will crash SigiHMM. The catch is the tool being used to convert from gbk to embl. Bioperl generates slightly different (possibly older) embl output than Biopython. IslandCompare uses Biopython for the conversion while my test originally used Bioperl. I can now consistently reproduce this for SigiHMM version 3.8. > SigiHMM 4.0 also cant parse Biopython output but emits rather than crashing: > > This line could not be parsed: seq:caa) > and discards a large amount of islands. This line is from 15584_genome.embl from Biopython: > > FT /anticodon=(pos:complement(1123552..1123554),aa:Leu, > FT seq:caa) > The matching line from 15584_genome.embl from Bioperl: > > FT /anticodon=(pos:complement(1123552..1123554),aa:Leu,seq:caa > FT ) > Moving forward we have two options: > > Add a perl dependency specifically for SigiHMM (yuck) > Discard the offending anticodon annotations and keep biopython > The difference between biopython and bioperl seem to be the criteria that they use when determining where to line wrap. Bioperl could potentially have this issue if the anticodon annotations are longer in other datasets. I will contact Genbank as I believe the anticodon annotations are actually not conforming to their own standard. > > I'll need @fionabrinkman and @cbertell to weigh in on this. > > — > You are receiving this because you were mentioned. > Reply to this email directly, view it on GitHub, or mute the thread.

[innovate-invent]

Can we discuss the update in the relevant thread. It makes it difficult to follow the issue.
SIGI update discussion: #134

[innovate-invent]

How line wrapping is handled is not in the Genbank/EMBL format specification. I have contacted NIH requesting clarification. They have been good in the past at getting back quickly. If this is a bug in biopython I would prefer to fix it rather than the alternatives.

Edit: I corrected a mistake in my post above analysing the difference between biopython and bioperl output. I must have opened the wrong file when I inspected the bioperl output because when I reran it today all wrapping /anticodon annotations are now quoted.

[innovate-invent]

A change to the biopython converter has resolved the issue for now.
brinkmanlab/galaxy-tools@164d4c3

[innovate-invent]

SigiHMM 3.8 is now crashing for a different reason.
Something about 15596_genome.gbk is throwing:

Exception in thread "main" java.lang.StringIndexOutOfBoundsException: String index out of range: 10
	at java.lang.String.substring(String.java:1963)
	at io.RowEMBL.parseSingle(RowEMBL.java:409)
	at io.RowEMBL.parseLocation(RowEMBL.java:370)
	at io.RowEMBL.setLocation(RowEMBL.java:198)
	at io.FileEMBL.parse(FileEMBL.java:331)
	at io.BatchFileReader.readFile(BatchFileReader.java:57)
	at GenericBatch.execute(GenericBatch.java:140)
	at SigiHMM.main(SigiHMM.java:10)

This exception is thrown for data from both bioperl and biopython.
SigiHMM version 4.0 does not throw this exception.

[innovate-invent]

This issue is side stepped by filtering failed datasets from the workflow.
Issue #134 is now the best course of action.

Referencing relevant issues and closing:
biopython/biopython#2112
biojava/biojava#843