No error in haplotype with incompatible variations

[ljournot]

Expected behavior

Hello,
I need to check the compatibility of multiple variations of the same sequence. I expect two incompatible variations to error whether there is another variation 5' or not.

Current behavior

I got some surprising behavior as illustrated by the MWE hereafter. I do not understand why some haplotypes with 3 variations do not error. The ones with only the 2 conflicting variations, or with a third variation after, do error, which is what I am expecting. It looks like the presence of a variation before, not after, the incompatible variations prevents the detection of the incompatibility.
I surely misunderstood something and will appreciate any help.
Laurent

Possible solution

No response

Steps to reproduce

julia> using BioSequences, SequenceVariation

julia> s = dna"AAGGCCTTCCGGTTAA"
16nt DNA Sequence:
AAGGCCTTCCGGTTAA

julia> v1 = Variation(s, "Δ2-3")
Δ2-3

julia> v2 = Variation(s, "Δ8-9")
Δ8-9

julia> v3 = Variation(s, "8AAA")
8AAA

julia> Haplotype(s, [v1, v2, v3])
Haplotype{LongSequence{DNAAlphabet{4}}, DNA} with 3 edits:
  Δ2-3
  Δ8-9
  8AAA

julia> Haplotype(s, [v2, v3])
ERROR: Multiple insertions at same position
Stacktrace:
 [1] error(s::String)
   @ Base .\error.jl:35
 [2] Haplotype{LongSequence{…}, DNA}(ref::LongSequence{DNAAlphabet{…}}, edits::Vector{SequenceVariation.Edit{…}})
   @ SequenceVariation C:\Users\Laurent\.julia\packages\SequenceVariation\y7qGk\src\Haplotype.jl:33
 [3] Haplotype(ref::LongSequence{DNAAlphabet{4}}, vars::Vector{Variation{LongSequence{DNAAlphabet{4}}, DNA}})
   @ SequenceVariation C:\Users\Laurent\.julia\packages\SequenceVariation\y7qGk\src\Variation.jl:50
 [4] top-level scope
   @ REPL[7]:1
Some type information was truncated. Use `show(err)` to see complete types.

julia> v4 = Variation(s, "G4C")
G4C

julia> Haplotype(s, [v4, v2, v3])
Haplotype{LongSequence{DNAAlphabet{4}}, DNA} with 3 edits:
  G4C
  Δ8-9
  8AAA

julia> v5 = Variation(s, "T14G")
T14G

julia> Haplotype(s, [v2, v3, v5])
ERROR: Multiple insertions at same position
Stacktrace:
 [1] error(s::String)
   @ Base .\error.jl:35
 [2] Haplotype{LongSequence{…}, DNA}(ref::LongSequence{DNAAlphabet{…}}, edits::Vector{SequenceVariation.Edit{…}})
   @ SequenceVariation C:\Users\Laurent\.julia\packages\SequenceVariation\y7qGk\src\Haplotype.jl:33
 [3] Haplotype(ref::LongSequence{DNAAlphabet{4}}, vars::Vector{Variation{LongSequence{DNAAlphabet{4}}, DNA}})
   @ SequenceVariation C:\Users\Laurent\.julia\packages\SequenceVariation\y7qGk\src\Variation.jl:50
 [4] top-level scope
   @ REPL[13]:1
Some type information was truncated. Use `show(err)` to see complete types.

Context

No response

Package version used

0.2.2

Julia version used

1.11.0

Operating system and version

Windows 10

Link to your project

No response

Please list installed packages here

(@v1.11) pkg> st
Status `C:\Users\Laurent\.julia\environments\v1.11\Project.toml`
  [6e4b80f9] BenchmarkTools v1.6.0
  [7e6ae17a] BioSequences v3.4.1
  [336ed68f] CSV v0.10.15
  [13f3f980] CairoMakie v0.13.1
  [944b1d66] CodecZlib v0.7.6
  [a93c6f00] DataFrames v1.7.0
  [31c24e10] Distributions v0.25.117
  [c2308a5c] FASTX v2.1.7
  [899a7d2d] GenomicFeatures v3.0.0
  [7073ff75] IJulia v1.26.0
  [d8418881] Intervals v1.10.0
  [033835bb] JLD2 v0.5.11
  [4345ca2d] Loess v0.6.4
  [f8716d33] MultipleTesting v0.6.0
  [eef6e190] SequenceVariation v0.2.2
  [2913bbd2] StatsBase v0.34.4
  [37e2e46d] LinearAlgebra v1.11.0

[ljournot]

A follow up on this issue. I realized that the order of variations matters when testing their compatibility, which is not what I was expecting.

julia> s = dna"ACGTAA"
6nt DNA Sequence:
ACGTAA

julia> var1 = Variation(s, "2TTT")
2TTT

julia> var2 = Variation(s, "Δ2-3")
Δ2-3

julia> Haplotype(s, [var1, var2])
Haplotype{LongSequence{DNAAlphabet{4}}, DNA} with 2 edits:
  2TTT
  Δ2-3

julia> Haplotype(s, [var2, var1])
ERROR: Multiple insertions at same position
Stacktrace:
 [1] error(s::String)
   @ Base .\error.jl:35
 [2] Haplotype{LongSequence{…}, DNA}(ref::LongSequence{DNAAlphabet{…}}, edits::Vector{SequenceVariation.Edit{…}})
   @ SequenceVariation C:\Users\Laurent\.julia\packages\SequenceVariation\y7qGk\src\Haplotype.jl:33
 [3] Haplotype(ref::LongSequence{DNAAlphabet{4}}, vars::Vector{Variation{LongSequence{DNAAlphabet{4}}, DNA}})
   @ SequenceVariation C:\Users\Laurent\.julia\packages\SequenceVariation\y7qGk\src\Variation.jl:50
 [4] top-level scope
   @ REPL[6]:1
Some type information was truncated. Use `show(err)` to see complete types.

I guess this behavior is linked to the _is_valid function. Because the behavior occurs only with Insertions, I wonder if this is linked to a problem I noted in the definition of the positions of Insertions. In the documentation, the constructor for Variations mentions that Insertions are coded the following way:

Insertions	<POS><ALT>	<ALT> is inserted between positions <POS> and <POS>+1	"11T". 

However, <ALT> is inserted between positions <POS-1> and <POS> as shown below.

julia> s
6nt DNA Sequence:
ACGTAA

julia> reconstruct(Haplotype(s, [var1]))
9nt DNA Sequence:
ATTTCGTAA

I hope this will help.

[ljournot]

A final MWE to pinpoint the problem.

julia> using BioSequences, SequenceVariation

julia> s = dna"ACGTGATC"
8nt DNA Sequence:
ACGTGATC

julia> v1 = Variation(s, "G3C")
G3C

julia> v2 = Variation(s, "3AAA")
3AAA

julia> reconstruct(Haplotype(s, [v1]))
8nt DNA Sequence:
ACCTGATC

julia> reconstruct(Haplotype(s, [v2]))
11nt DNA Sequence:
ACAAAGTGATC

julia> reconstruct(Haplotype(s, [v1, v2]))
11nt DNA Sequence:
ACCAAATGATC

julia> reconstruct(Haplotype(s, [v2, v1]))
11nt DNA Sequence:
ACAAACTGATC

The sequences reconstructed from [v1], [v2], and [v2, v1] are what I expect; the one from [v1, v2] look wrong to me despite the haplotype being valid (no error message). I understand this happens because the variations are reconstructed serially. Still, I suppose we should consider this is a problem. For now, I will check the order of the variations before reconstructing the sequences. Is there a better way to manage such a strange behavior?