Compare commits

..

No commits in common. "2ff4181dd8a4f5fc173d2b694194cac9fa690d21" and "8a71715cd9fa3826da930f0caa965c109810232f" have entirely different histories.

7 changed files with 16 additions and 21 deletions

View file

@ -16,7 +16,6 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
- Code now follows [Blue style](https://github.com/invenia/BlueStyle) ([#28](https://github.com/BioJulia/SequenceVariation.jl/pull/28))
- :bomb: [BREAKING] Public and private API defined based on Blue style guidelines ([#28](https://github.com/BioJulia/SequenceVariation.jl/pull/28))
- :bomb: [BREAKING] Renamed type `Variant` to `Haplotype` ([#20](https://github.com/BioJulia/SequenceVariation.jl/issues/20)/[#29](https://github.com/BioJulia/SequenceVariation.jl/pull/29))
- :bomb: [BREAKING] Refactored `reconstruct` function to use `Haplotype`'s reference sequence ([#30](https://github.com/BioJulia/SequenceVariation.jl/pull/30))
### Removed

View file

@ -21,7 +21,7 @@ Insertion
Haplotype
reference(::Haplotype)
variations
reconstruct
reconstruct!
```
## Variations

View file

@ -31,10 +31,11 @@ bos_human_haplotype = Haplotype(bos_human_alignment)
If the alternate sequence of a haplotype is no longer available (as is often the
case when calling variants from alignment files), then the sequence can be
retrieved using the [`reconstruct`](@ref) function.
retrieved using the [`reconstruct!`](@ref) function.
```@repl call_variants
human2 = reconstruct(bos_human_haplotype)
human2 = copy(bovine);
reconstruct!(human2, bos_human_haplotype)
human2 == bovine
human2 == human
```

View file

@ -162,22 +162,21 @@ reference(h::Haplotype) = h.ref
Base.:(==)(x::Haplotype, y::Haplotype) = x.ref == y.ref && x.edits == y.edits
"""
reconstruct(h::Haplotype)
reconstruct!(seq::S, x::Haplotype{S}) where {S}
Apply the edits in `h` to the reference sequence of `h` and return the mutated sequence
Apply the edits in `x` to `seq` and return the mutated sequence
"""
function reconstruct(h::Haplotype)
len = length(reference(h)) + sum(edit -> _lendiff(edit), _edits(h))
seq = copy(reference(h))
function reconstruct!(seq::S, x::Haplotype{S}) where {S}
len = length(x.ref) + sum(edit -> _lendiff(edit), _edits(x))
resize!(seq, len % UInt)
refpos = seqpos = 1
for edit in _edits(h)
while refpos < leftposition(edit)
seq[seqpos] = reference(h)[refpos]
for edit in x.edits
while refpos < edit.pos
seq[seqpos] = x.ref[refpos]
refpos += 1
seqpos += 1
end
editx = _mutation(edit)
editx = edit.x
if editx isa Substitution
seq[seqpos] = editx.x
seqpos += 1
@ -185,14 +184,14 @@ function reconstruct(h::Haplotype)
elseif editx isa Deletion
refpos += editx.len
elseif editx isa Insertion
for i in editx.seq
for i in editx.x
seq[seqpos] = i
seqpos += 1
end
end
end
while seqpos length(seq)
seq[seqpos] = reference(h)[refpos]
seq[seqpos] = x.ref[refpos]
refpos += 1
seqpos += 1
end

View file

@ -32,7 +32,7 @@ export Substitution
export Variation
export altbases
export mutation
export reconstruct
export reconstruct!
export refbases
export reference
export translate

View file

@ -15,7 +15,7 @@ end
Deletion(x::Integer) = Deletion(convert(UInt, x))
Base.length(x::Deletion) = Int(x.len)
Base.:(==)(x::Deletion, y::Deletion) = length(x) == length(y)
Base.:(==)(x::Substitution, y::Substitution) = length(x) == length(y)
Base.hash(x::Deletion, h::UInt) = hash(Deletion, hash(x.len, h))
function _refbases(d::Deletion, reference::S, pos::UInt) where {S<:BioSequence}

View file

@ -43,10 +43,6 @@ var = Haplotype(align(seq1, seq2))
end
end
@testset "HaplotypeReconstruction" begin
@test reconstruct(var) == seq1
end
@testset "VariationPosition" begin
refseq = dna"ACAACTTTATCT"
mutseq = dna"ACATCTTTATCT"