Repository: cogent3/cogent3
Pull request #1412 (2023-05-05 06:07:32)

Annotation refactor -- consistent functionality and passing tests


Patch:
diff --git a/src/cogent3/core/alignment.py b/src/cogent3/core/alignment.py
index 7a2742ff1..6e4974cb2 100644
--- a/src/cogent3/core/alignment.py
+++ b/src/cogent3/core/alignment.py
@@ -53,7 +53,7 @@ from numpy.random import choice, permutation, randint
 
 import cogent3  # will use to get at cogent3.parse.fasta.MinimalFastaParser,
 
-from cogent3.core.annotation import Annotation, Map, _Annotatable
+from cogent3.core.annotation import Annotation, Map
 from cogent3.core.annotation_db import (
     FeatureDataType,
     GenbankAnnotationDb,
@@ -65,7 +65,6 @@ from cogent3.core.genetic_code import get_code
 from cogent3.core.info import Info as InfoClass
 from cogent3.core.profile import PSSM, MotifCountsArray
 from cogent3.core.sequence import ArraySequence, Sequence, frac_same
-
 # which is a circular import otherwise.
 from cogent3.format.alignment import save_to_filename
 from cogent3.format.fasta import alignment_to_fasta
@@ -644,8 +643,12 @@ class _SequenceCollectionBase:
             interval. This also causes dropping
             annotations.
         """
+        if isinstance(self, Alignment):
+            reversed = self.seqs[0].map.reverse
+        else:
+            reversed = self.seqs[0]._seq.reversed
         new_seqs = dict()
-        db = None if sliced else deepcopy(self.annotation_db)
+        db = None if reversed and sliced else deepcopy(self.annotation_db)
         for seq in self.seqs:
             try:
                 new_seq = seq.deepcopy(sliced=sliced, exclude_annotations=True)
@@ -926,11 +929,15 @@ class _SequenceCollectionBase:
             # copy only the specified seqs
             for r in seqs:
                 result[r] = get(r)
-        if result:
-            return self.__class__(result, names=seqs, info=self.info, **kwargs)
-        else:
+        if not result:
             return {}  # safe value; can't construct empty alignment
 
+        result = self.__class__(result, names=seqs, info=self.info, **kwargs)
+        if self.annotation_db:
+            result.annotation_db = self.annotation_db.__class__()
+            result.annotation_db.update(self.annotation_db, seqids=result.names)
+        return result
+
     def get_seq_indices(self, f, negate=False):
         """Returns list of keys of seqs where f(row) is True.
 
@@ -1160,13 +1167,8 @@ class _SequenceCollectionBase:
             version=__version__,
         )
 
-        try:
-            annotations = [a.to_rich_dict() for a in self.annotations]
-        except AttributeError:
-            annotations = []
-
-        if annotations:
-            data["annotations"] = annotations
+        if hasattr(self, "annotation_db") and self.annotation_db:
+            data["annotation_db"] = self.annotation_db.to_rich_dict()
 
         return data
 
@@ -1224,17 +1226,8 @@ class _SequenceCollectionBase:
             data=list(zip(self.names, concatenated)),
             info=self.info,
         )
-
-        if aligned:
-            left = [
-                a for a in self._shifted_annotations(new, 0) if a.map.end <= len(self)
-            ]
-            right = [
-                a
-                for a in other._shifted_annotations(new, len(self))
-                if a.map.start >= len(self)
-            ]
-            new.annotations = left + right
+        # cannot copy annotations on addition
+        new.annotation_db = None
         return new
 
     def add_seqs(self, other, before_name=None, after_name=None):
@@ -1684,39 +1677,24 @@ class _SequenceCollectionBase:
             raise ValueError(f"unknown moltype '{moltype}'")
 
         data = [s.to_moltype(moltype) for s in self.seqs]
-        return self.__class__(
+        result = self.__class__(
             data=data, moltype=moltype, name=self.name, info=self.info
         )
+        if self.annotation_db:
+            result.annotation_db = deepcopy(self.annotation_db)
+        return result
 
     def to_dna(self):
         """returns copy of self as an alignment of DNA moltype seqs"""
-        data = [s.to_moltype("dna") for s in self.seqs]
-        new = self.__class__(
-            data=data, moltype=cogent3.DNA, name=self.name, info=self.info
-        )
-        if isinstance(self, _Annotatable) and self.annotations:
-            new.annotations = self.annotations[:]
-        return new
+        return self.to_moltype("dna")
 
     def to_rna(self):
         """returns copy of self as an alignment of RNA moltype seqs"""
-        data = [s.to_moltype("rna") for s in self.seqs]
-        new = self.__class__(
-            data=data, moltype=cogent3.RNA, name=self.name, info=self.info
-        )
-        if isinstance(self, _Annotatable) and self.annotations:
-            new.annotations = self.annotations[:]
-        return new
+        return self.to_moltype("rna")
 
     def to_protein(self):
         """returns copy of self as an alignment of PROTEIN moltype seqs"""
-        data = [s.to_moltype("protein") for s in self.seqs]
-        new = self.__class__(
-            data=data, moltype=cogent3.PROTEIN, name=self.name, info=self.info
-        )
-        if isinstance(self, _Annotatable) and self.annotations:
-            new.annotations = self.annotations[:]
-        return new
+        return self.to_moltype("protein")
 
     def rc(self):
         """Returns the reverse complement alignment"""
@@ -1724,9 +1702,7 @@ class _SequenceCollectionBase:
         rc = self.__class__(
             data=seqs, names=self.names[:], name=self.name, info=self.info
         )
-        if isinstance(self, _Annotatable) and self.annotations:
-            self._annotations_nucleic_reversed_on(rc)
-        rc.annotation_db = self.annotation_db
+        rc.annotation_db = deepcopy(self.annotation_db)
         return rc
 
     def reverse_complement(self):
@@ -2264,13 +2240,6 @@ class SequenceCollection(_SequenceCollectionBase):
         yield a sequence segment that is consistently oriented irrespective
         of strand of the current instance.
         """
-        if self.annotation_db is None:
-            # todo gah this should not be required
-            # the annotation db attribute should be defined on creation,
-            # or when annotation_db is assigned
-            anno_db = merged_db_collection(self.seqs)
-            self.annotation_db = anno_db
-
         if self.annotation_db is None:
             return None
 
@@ -2285,8 +2254,10 @@ class SequenceCollection(_SequenceCollectionBase):
             seq = self.named_seqs[seqid]
             if isinstance(seq, Aligned):
                 start, end = seq.map.start, seq.map.end
+                offset = seq.data.annotation_offset
             else:
                 start, end = 0, len(seq)
+                offset = seq.annotation_offset
 
             for feature in self.annotation_db.get_features_matching(
                 seqid=seqid,
@@ -2297,6 +2268,8 @@ class SequenceCollection(_SequenceCollectionBase):
                 start=start,
                 end=end,
             ):
+                if offset:
+                    feature["spans"] = (array(feature["spans"]) - offset).tolist()
                 # passing self only used when self is an Alignment
                 yield seq.make_feature(feature, self)
 
@@ -2353,12 +2326,18 @@ class Aligned:
         a copy of self
         """
         new_seq = self.data.copy(exclude_annotations=exclude_annotations)
+        db = new_seq.annotation_db
         if sliced:
             span = self.map.get_covering_span()
             new_seq = type(new_seq)(
                 str(new_seq[span.start : span.end]), info=new_seq.info, name=self.name
             )
-            new_seq.annotation_db = None
+            new_seq.annotation_offset = self.map.start
+            if self.map.reverse or exclude_annotations:
+                new_seq.annotation_db = None
+            else:
+                new_seq.annotation_offset = self.map.start
+                new_seq.annotation_db = db
             new_map = self.map.zeroed()
         else:
             new_map = self.map
@@ -2450,10 +2429,7 @@ class Aligned:
             raise NotImplementedError
         start, end = coords[0]
         data = self.data[start:end]
-        # drop any lost spans
-        for i, a in enumerate(data.annotations):
-            data.annotations[i] = a.without_lost_spans()
-        data = data.to_rich_dict()
+        data = data.to_rich_dict(exclude_annotations=True)
         data["seq"] = str(self)
         data["version"] = __version__
         return data
@@ -3818,6 +3794,45 @@ class AlignmentI(object):
             width=width, height=height, wrap=wrap, vspace=vspace, colours=colours
         )
 
+    def get_drawable(self, width=600, vertical=False):
+        """returns Drawable instance"""
+        from cogent3.draw.drawable import Drawable
+
+        drawables = self.get_drawables()
+        if not drawables:
+            return None
+        # we order by tracks
+        top = 0
+        space = 0.25
+        annotes = []
+        for feature_type in drawables:
+            new_bottom = top + space
+            for i, annott in enumerate(drawables[feature_type]):
+                annott.shift(y=new_bottom - annott.bottom)
+                if i > 0:
+                    annott._showlegend = False
+                annotes.append(annott)
+
+            top = annott.top
+
+        top += space
+        height = max((top / len(self)) * width, 300)
+        xaxis = dict(range=[0, len(self)], zeroline=False, showline=True)
+        yaxis = dict(range=[0, top], visible=False, zeroline=True, showline=True)
+
+        if vertical:
+            all_traces = [t.T.as_trace() for t in annotes]
+            width, height = height, width
+            xaxis, yaxis = yaxis, xaxis
+        else:
+            all_traces = [t.as_trace() for t in annotes]
+
+        drawer = Drawable(
+            title=self.name, traces=all_traces, width=width, height=height
+        )
+        drawer.layout.update(xaxis=xaxis, yaxis=yaxis)
+        return drawer
+
 
 def _one_length(seqs):
     """raises ValueError if seqs not all same length"""
@@ -4668,7 +4683,7 @@ def make_gap_filter(template, gap_fraction, gap_run):
     return result
 
 
-class Alignment(_Annotatable, AlignmentI, SequenceCollection):
+class Alignment(AlignmentI, SequenceCollection):
     """An annotatable alignment class"""
 
     moltype = None  # note: this is reset to ASCII in moltype module
@@ -4702,7 +4717,7 @@ class Alignment(_Annotatable, AlignmentI, SequenceCollection):
     def _coerce_seqs(self, seqs, is_array):
         if any(isinstance(seq, ArraySequence) for seq in seqs):
             seqs = [self.moltype.make_seq(str(seq), name=seq.name) for seq in seqs]
-        elif not any(isinstance(seq, (_Annotatable, Aligned)) for seq in seqs):
+        elif not any(isinstance(seq, (Sequence, Aligned)) for seq in seqs):
             seqs = list(map(self.moltype.make_seq, seqs))
         return seqs
 
@@ -4756,10 +4771,18 @@ class Alignment(_Annotatable, AlignmentI, SequenceCollection):
 
     @annotation_db.setter
     def annotation_db(self, value):
-        from cogent3.core.annotation_db import SupportsFeatures
-
         if value and not isinstance(value, SupportsFeatures):
             raise TypeError
+
+        # Without knowing the contents of the db we cannot
+        # establish whether self.moltype is compatible, so
+        # we rely on the user to get that correct
+        # one approach to support validation might be to add
+        # to the SupportsFeatures protocol a is_nucleic flag,
+        # for both DNA and RNA. But if a user trys get_slice()
+        # on a '-' strand feature, they will get a TypError.
+        # I think that's enough.
+
         self._annotation_db = value
         for seq in self.seqs:
             seq.data.annotation_db = value
@@ -5108,19 +5131,14 @@ class Alignment(_Annotatable, AlignmentI, SequenceCollection):
     def to_moltype(self, moltype):
         """returns copy of self with moltype seqs"""
         new = super().to_moltype(moltype)
-        for annot in self.annotations:
-            new = annot.copy_annotations_to(new)
+        new.annotation_db = self.annotation_db
         return new
 
     def get_drawables(self):
         """returns a dict of drawables, keyed by type"""
         result = defaultdict(list)
-        for a in self.get_annotations_from_any_seq():
-            d = a.get_drawable()
-            if not d:
-                continue
-            result[a.type].append(d)
-
+        for f in self.get_features():
+            result[f.biotype].append(f.get_drawable())
         return result
 
     def add_feature(
diff --git a/src/cogent3/core/annotation.py b/src/cogent3/core/annotation.py
index 967a6e96c..680a1c8be 100644
--- a/src/cogent3/core/annotation.py
+++ b/src/cogent3/core/annotation.py
@@ -592,6 +592,12 @@ class Annotation(_AnnotationCore, _Serialisable):
             parent=self.parent, seqid=seqid, map=map, biotype=biotype, name=name
         )
 
+    def get_drawable(self):
+        """returns plotly trace"""
+        from cogent3.draw.drawable import make_shape
+
+        return make_shape(type_=self)
+
 
 # https://pythonspeed.com/products/filmemoryprofiler/
 
diff --git a/src/cogent3/core/annotation_db.py b/src/cogent3/core/annotation_db.py
index 322e19fd7..23ea5ac3d 100644
--- a/src/cogent3/core/annotation_db.py
+++ b/src/cogent3/core/annotation_db.py
@@ -921,11 +921,9 @@ class GenbankAnnotationDb(SqliteAnnotationDbMixin):
             store = {"seqid": seqid}
             # we create the location data directly
             if location := record.get("location", None):
-                store["spans"] = numpy.array(
-                    sorted([sorted((s.start, s.stop)) for s in location]), dtype=int
-                )
+                store["spans"] = numpy.array(location.get_coordinates(), dtype=int)
                 if strand := location.strand:
-                    store["strand"] = strand
+                    store["strand"] = "-" if strand == -1 else "+"
                 store["start"] = int(store["spans"].min())
                 store["end"] = int(store["spans"].max())
 
@@ -1059,12 +1057,12 @@ def convert_annotation_to_annotation_db(data: dict) -> dict:
 
 def _db_from_genbank(path, db):
     from cogent3 import open_
-    from cogent3.parse.genbank import MinimalGenbankParser
+    from cogent3.parse.genbank import RichGenbankParser
 
     with open_(path) as infile:
-        data = list(MinimalGenbankParser(infile))
+        data = list(RichGenbankParser(infile, db=db))[0][1]
 
-    return GenbankAnnotationDb(data=data[0]["features"], seqid=data[0]["locus"], db=db)
+    return data.annotation_db
 
 
 def _leave_attributes(*attrs):
diff --git a/src/cogent3/core/sequence.py b/src/cogent3/core/sequence.py
index b6e09188e..93519aef2 100644
--- a/src/cogent3/core/sequence.py
+++ b/src/cogent3/core/sequence.py
@@ -17,6 +17,7 @@ import os
 import re
 import warnings
 
+from collections import defaultdict
 from functools import singledispatch, total_ordering
 from operator import eq, ne
 from random import shuffle
@@ -37,7 +38,7 @@ from numpy import (
 from numpy.random import permutation
 
 from cogent3.core.alphabet import AlphabetError
-from cogent3.core.annotation import Annotation, _Annotatable
+from cogent3.core.annotation import Annotation
 from cogent3.core.annotation_db import (
     FeatureDataType,
     GenbankAnnotationDb,
@@ -144,7 +145,7 @@ class SequenceI(object):
             version=__version__,
         )
         if hasattr(self, "annotation_offset"):
-            offset = self.annotation_offset or self._seq.start
+            offset = self._seq.offset + self._seq.start
             data.update(dict(annotation_offset=offset))
 
         if (
@@ -793,7 +794,7 @@ class SequenceI(object):
 
 
 @total_ordering
-class Sequence(_Annotatable, SequenceI):
+class Sequence(SequenceI):
     """Holds the standard Sequence object. Immutable."""
 
     moltype = None  # connected to ACSII when moltype is imported
@@ -846,12 +847,7 @@ class Sequence(_Annotatable, SequenceI):
                 info.update(orig_seq.info)
         self.info = info
 
-        if isinstance(orig_seq, _Annotatable):
-            for ann in orig_seq.annotations:
-                ann.copy_annotations_to(self)
-
         self._repr_policy = dict(num_pos=60)
-
         self._annotation_db = None
         self.annotation_offset = annotation_offset
 
@@ -884,10 +880,16 @@ class Sequence(_Annotatable, SequenceI):
 
     @annotation_db.setter
     def annotation_db(self, value):
-        from cogent3.core.annotation_db import SupportsFeatures
-
         if value and not isinstance(value, SupportsFeatures):
             raise TypeError
+        # Without knowing the contents of the db we cannot
+        # establish whether self.moltype is compatible, so
+        # we rely on the user to get that correct
+        # one approach to support validation might be to add
+        # to the SupportsFeatures protocol a is_nucleic flag,
+        # for both DNA and RNA. But if a user trys get_slice()
+        # on a '-' strand feature, they will get a TypError.
+        # I think that's enough.
         self._annotation_db = value
 
     @deprecated_args(
@@ -1152,10 +1154,7 @@ class Sequence(_Annotatable, SequenceI):
         moltype.verify_sequence(s, gaps_allowed=True, wildcards_allowed=True)
         sv = SeqView(s)
         new = moltype.make_seq(sv, name=self.name, info=self.info)
-
-        new.clear_annotations()
-        for ann in self.annotations:
-            ann.copy_annotations_to(new)
+        new.annotation_db = self.annotation_db
         return new
 
     def _seq_filter(self, seq):
@@ -1192,7 +1191,10 @@ class Sequence(_Annotatable, SequenceI):
 
     def copy(self, exclude_annotations=False):
         """returns a copy of self"""
-        new = self.__class__(self._seq[:], name=self.name, info=self.info)
+        offset = self._seq.offset + self._seq.start
+        new = self.__class__(
+            self._seq[:], name=self.name, info=self.info, annotation_offset=offset
+        )
         db = None if exclude_annotations else copy.deepcopy(self.annotation_db)
         new._annotation_db = db
         return new
@@ -1286,12 +1288,11 @@ class Sequence(_Annotatable, SequenceI):
             yield from segment
 
     def gapped_by_map(self, map, recode_gaps=False):
+        # todo gah do we propagate annotations here?
         segments = self.gapped_by_map_segment_iter(map, True, recode_gaps)
         new = self.__class__(
             "".join(segments), name=self.name, check=False, info=self.info
         )
-        annots = self._sliced_annotations(new, map)
-        new.annotations = annots
         return new
 
     def _mapped(self, map):
@@ -1309,21 +1310,24 @@ class Sequence(_Annotatable, SequenceI):
         return f"{myclass}({seq})"
 
     def __getitem__(self, index):
+        preserve_offset = False
         if hasattr(index, "map"):
             index = index.map
 
-        # todo: kath preserve the offset?
-
         if isinstance(index, Map):
             new = self._mapped(index)
+            preserve_offset = not index.reverse
 
         elif isinstance(index, (int, slice)):
             new = self.__class__(
                 self._seq[index], name=self.name, check=False, info=self.info
             )
+            stride = getattr(index, "step", 1) or 1
+            preserve_offset = stride > 0
 
-        if self.annotation_db is not None:
+        if self.annotation_db is not None and preserve_offset:
             new.annotation_db = self.annotation_db
+            new.annotation_offset = self.annotation_offset
 
         if hasattr(self, "_repr_policy"):
             new._repr_policy.update(self._repr_policy)
@@ -1449,8 +1453,8 @@ class Sequence(_Annotatable, SequenceI):
         seq = self.__class__(
             "".join(gapless), name=self.get_name(), info=self.info, preserve_case=True
         )
-        if self.annotations:
-            seq.annotations = [a.remapped_to(seq, map) for a in self.annotations]
+        if self.annotation_db:
+            seq.annotation_db = self.annotation_db
         return (map, seq)
 
     def replace(self, oldchar, newchar):
@@ -1461,7 +1465,7 @@ class Sequence(_Annotatable, SequenceI):
     def is_annotated(self):
         """returns True if sequence has any annotations"""
         num = self.annotation_db.num_matches() if self.annotation_db else 0
-        return num != 0 or len(self.annotations) != 0
+        return num != 0
 
     def annotate_matches_to(self, pattern, annot_type, name, allow_multiple=False):
         """Adds an annotation at sequence positions matching pattern.
@@ -1506,24 +1510,55 @@ class Sequence(_Annotatable, SequenceI):
     def __add__(self, other):
         """Adds two sequences (other can be a string as well)"""
         new_seq = super(Sequence, self).__add__(other)
-        # Annotations which extend past the right end of the left sequence
-        # or past the left end of the right sequence are dropped because
-        # otherwise they will annotate the wrong part of the constructed
-        # sequence.
-        left = [
-            a for a in self._shifted_annotations(new_seq, 0) if a.map.end <= len(self)
-        ]
-        if hasattr(other, "_shifted_annotations"):
-            right = [
-                a
-                for a in other._shifted_annotations(new_seq, len(self))
-                if a.map.start >= len(self)
-            ]
-            new_seq.annotations = left + right
+        new_seq.annotation_db = None
+        # annotations are dropped in this case
+        return new_seq
+
+    def get_drawables(self):
+        """returns a dict of drawables, keyed by type"""
+        result = defaultdict(list)
+        for f in self.get_features():
+            result[f.biotype].append(f.get_drawable())
+        return result
+
+    def get_drawable(self, width=600, vertical=False):
+        """returns Drawable instance"""
+        from cogent3.draw.drawable import Drawable
+
+        drawables = self.get_drawables()
+        if not drawables:
+            return None
+        # we order by tracks
+        top = 0
+        space = 0.25
+        annotes = []
+        for feature_type in drawables:
+            new_bottom = top + space
+            for i, annott in enumerate(drawables[feature_type]):
+                annott.shift(y=new_bottom - annott.bottom)
+                if i > 0:
+                    annott._showlegend = False
+                annotes.append(annott)
+
+            top = annott.top
+
+        top += space
+        height = max((top / len(self)) * width, 300)
+        xaxis = dict(range=[0, len(self)], zeroline=False, showline=True)
+        yaxis = dict(range=[0, top], visible=False, zeroline=True, showline=True)
+
+        if vertical:
+            all_traces = [t.T.as_trace() for t in annotes]
+            width, height = height, width
+            xaxis, yaxis = yaxis, xaxis
         else:
-            new_seq.annotations = left
+            all_traces = [t.as_trace() for t in annotes]
 
-        return new_seq
+        drawer = Drawable(
+            title=self.name, traces=all_traces, width=width, height=height
+        )
+        drawer.layout.update(xaxis=xaxis, yaxis=yaxis)
+        return drawer
 
 
 class ProteinSequence(Sequence):
@@ -1751,7 +1786,7 @@ def _input_vals_neg_step(seqlen, start, stop, step):
 
 
 class SeqView:
-    __slots__ = ("seq", "start", "stop", "step", "offset")
+    __slots__ = ("seq", "start", "stop", "step", "_offset")
 
     def __init__(self, seq, *, start: int = None, stop: int = None, step: int = None):
         if step == 0:
@@ -1763,9 +1798,18 @@ class SeqView:
         self.seq = seq
         self.start = start
         self.stop = stop
-        self.offset = None  # todo gah this should default to 0
+        self._offset = 0
         self.step = step
 
+    @property
+    def offset(self) -> int:
+        return self._offset
+
+    @offset.setter
+    def offset(self, value: int):
+        self._offset = value or 0
+
+    # todo gah rename to reverse for compatability with Map
     @property
     def reversed(self):
         return self.step < 0
@@ -1791,7 +1835,7 @@ class SeqView:
         seq_index, _, _ = self._get_index(rel_index, include_boundary=include_boundary)
 
         # add offset and handle reversed views, now absolute relative to annotation coordinates
-        offset = 0 if self.offset is None else self.offset
+        offset = self.offset
         if self.reversed:
             abs_index = offset + len(self.seq) + seq_index + 1
         else:
@@ -1810,7 +1854,7 @@ class SeqView:
             raise IndexError("Index must be +ve and relative to the + strand")
 
         if self.reversed:
-            offset = 0 if self.offset is None else self.offset
+            offset = self.offset
 
             if (
                 tmp := ((len(self.seq) - (abs_index - offset)) + self.start + 1)
@@ -1820,7 +1864,7 @@ class SeqView:
                 rel_pos = (tmp // abs(self.step)) + 1
 
         else:
-            offset = self.start if self.offset is None else self.offset + self.start
+            offset = self.offset + self.start
 
             if (tmp := abs_index - offset) % self.step == 0 or stop:
                 rel_pos = tmp // self.step
diff --git a/src/cogent3/draw/drawable.py b/src/cogent3/draw/drawable.py
index 09c3e524b..f720f1e0d 100644
--- a/src/cogent3/draw/drawable.py
+++ b/src/cogent3/draw/drawable.py
@@ -776,21 +776,21 @@ class _MakeShape:
     def __call__(self, type_=None, name=None, coords=None, **kwargs):
         from cogent3.core.annotation import _Annotatable
 
-        if isinstance(type_, _Annotatable):
+        if hasattr(type_, "map"):
             if not type_.map.useful:
                 return None
 
             name = type_.name
             coords = type_.map.get_coordinates()
             reverse = type_.map.get_covering_span().reverse
-            type_ = type_.type
+            type_ = type_.biotype
         else:
             if coords[0][0] > coords[-1][1]:
                 reverse = True
             else:
                 reverse = False
             if coords is None:
-                raise Exception("No coordinates defined")
+                raise ValueError("No coordinates defined")
         kwargs.update(dict(reverse=reverse))
 
         klass = self._shapes.get(type_.lower(), Rectangle)
diff --git a/src/cogent3/parse/gbseq.py b/src/cogent3/parse/gbseq.py
index 71fe059d8..6b7cd32d8 100644
--- a/src/cogent3/parse/gbseq.py
+++ b/src/cogent3/parse/gbseq.py
@@ -74,13 +74,15 @@ def GbSeqXmlParser(doc):
         seq = alphabet.make_seq(raw_string, name=name)
 
         all = annotation.Map([(0, len(seq))], parent_length=len(seq))
-        seq.add_annotation(annotation.Source, all, name, all)
+        seq.add_feature(
+            biotype="source", name=name, spans=all.get_coordinates(), strand=all.reverse
+        )
 
         organism = str(
             record.getElementsByTagName("GBSeq_organism")[0].childNodes[0].nodeValue
         )
 
-        seq.add_annotation(annotation.Feature, "organism", organism, [(0, len(seq))])
+        seq.add_feature(biotype="organism", name=organism, spans=[(0, len(seq))])
 
         features = record.getElementsByTagName("GBFeature")
         for feature in features:
@@ -128,7 +130,7 @@ def GbSeqXmlParser(doc):
                         .childNodes[0]
                         .nodeValue
                     )
-            seq.add_annotation(annotation.Feature, key, feature_name, spans)
+            seq.add_feature(biotype=key, name=feature_name, spans=spans)
         yield (name, seq)
 
 
diff --git a/src/cogent3/parse/genbank.py b/src/cogent3/parse/genbank.py
index 162ce59ea..7464b61d9 100644
--- a/src/cogent3/parse/genbank.py
+++ b/src/cogent3/parse/genbank.py
@@ -1,5 +1,6 @@
-#!/usr/bin/env python
-from cogent3.core.annotation import Feature
+from typing import List, Optional, Tuple
+
+from cogent3.core.annotation_db import GenbankAnnotationDb
 from cogent3.core.genetic_code import GeneticCodes
 from cogent3.core.info import Info
 from cogent3.core.moltype import get_moltype
@@ -8,6 +9,7 @@ from cogent3.parse.record_finder import (
     DelimitedRecordFinder,
     LabeledRecordFinder,
 )
+from cogent3.util import warning as c3warn
 
 
 __author__ = "Rob Knight"
@@ -485,13 +487,8 @@ class LocationList(list):
     @property
     def strand(self):
         """Returns strand of components: 1=forward, -1=reverse, 0=both"""
-        curr = {}
-        for i in self:
-            curr[i.strand] = 1
-        if len(curr) >= 2:  # found stuff on both strands
-            return 0
-        else:
-            return list(curr.keys())[0]
+        curr = {i.strand: 1 for i in self}
+        return 0 if len(curr) >= 2 else list(curr.keys())[0]
 
     def __str__(self):
         """Returns (normalized) string representation of self."""
@@ -518,6 +515,10 @@ class LocationList(list):
             result.append(curr)
         return "".join(result)
 
+    def get_coordinates(self) -> List[Tuple[int, int]]:
+        """returns the segments in python coordinates"""
+        return sorted((i.start, i.stop + 1) for i in self)
+
 
 def parse_feature_table(lines):
     """Simple parser for feature table. Assumes starts with FEATURES line."""
@@ -678,8 +679,16 @@ def extract_nt_prot_seqs(rec, wanted=wanted_types):
         print("s :", seq)
 
 
+@c3warn.deprecated_args(
+    "2023.7", reason="no longer relevant", discontinued="add_annotation"
+)
 def RichGenbankParser(
-    handle, info_excludes=None, moltype=None, skip_contigs=False, add_annotation=None
+    handle,
+    info_excludes=None,
+    moltype=None,
+    skip_contigs=False,
+    add_annotation=None,
+    db: Optional[GenbankAnnotationDb] = None,
 ):
     """Returns annotated sequences from GenBank formatted file.
 
@@ -692,11 +701,9 @@ def RichGenbankParser(
     skip_contigs
         ignores records with no actual sequence data, typically
         a genomic contig.
-    add_annotation
-        a callback function to create an new annotation from a
-        GenBank feature. Function is called with the sequence, a feature dict
-        and the feature spans.
-
+    db
+        a GenbankAnnotationDb instance to which feature data will be
+        added
     """
     info_excludes = info_excludes or []
     moltype = get_moltype(moltype) if moltype else None
@@ -731,42 +738,8 @@ def RichGenbankParser(
                     yield rec["locus"], None
             continue
 
-        for feature in rec["features"]:
-            spans = []
-            reversed = None
-            if feature["location"] is None or feature["type"] in ["source", "organism"]:
-                continue
-            for location in feature["location"]:
-                (lo, hi) = (location.start, location.stop + 1)
-                if location.strand == -1:
-                    (lo, hi) = (hi, lo)
-                    assert reversed is not False
-                    reversed = True
-                else:
-                    assert reversed is not True
-                    reversed = False
-                # ensure we don't put in a span that starts beyond the sequence
-                if lo > len(seq):
-                    continue
-                # or that's longer than the sequence
-                hi = [hi, len(seq)][hi > len(seq)]
-                spans.append((lo, hi))
-
-            if add_annotation:
-                add_annotation(seq, feature, spans)
-            else:
-                for id_field in ["gene", "product", "clone", "note"]:
-                    if id_field in feature:
-                        name = feature[id_field]
-                        if not isinstance(name, str):
-                            name = " ".join(name)
-                        break
-                else:
-                    name = None
-                seq.add_annotation(Feature, feature["type"], name, spans)
-
+        db = getattr(db, "db", None)
+        seq.annotation_db = GenbankAnnotationDb(
+            data=rec["features"], seqid=rec["locus"], db=db
+        )
         yield (rec["locus"], seq)
-
-
-def parse(*args):
-    return RichGenbankParser(*args).next()[1]
diff --git a/src/cogent3/parse/tinyseq.py b/src/cogent3/parse/tinyseq.py
index 8af36c23c..2af805587 100644
--- a/src/cogent3/parse/tinyseq.py
+++ b/src/cogent3/parse/tinyseq.py
@@ -69,16 +69,12 @@ def TinyseqParser(doc):
 
         seq = alphabet.make_seq(raw_string, name=name)
 
-        seq.add_annotation(annotation.Feature, "genbank_id", name, [(0, len(seq))])
+        seq.add_feature(biotype="genbank_id", name=name, spans=[(0, len(seq))])
 
         organism = str(
             record.getElementsByTagName("TSeq_orgname")[0].childNodes[0].nodeValue
         )
 
-        seq.add_annotation(annotation.Feature, "organism", organism, [(0, len(seq))])
+        seq.add_feature(biotype="organism", name=organism, spans=[(0, len(seq))])
 
         yield (name, seq)
-
-
-def parse(*args):
-    return TinyseqParser(*args).next()[1]


Test patch:
diff --git a/tests/test_core/test_alignment.py b/tests/test_core/test_alignment.py
index f705742cb..d102c9817 100644
--- a/tests/test_core/test_alignment.py
+++ b/tests/test_core/test_alignment.py
@@ -47,6 +47,7 @@ from cogent3.core.alignment import (
     seqs_from_kv_pairs,
 )
 from cogent3.core.alphabet import AlphabetError
+from cogent3.core.annotation_db import GffAnnotationDb
 from cogent3.core.moltype import AB, ASCII, BYTES, DNA, PROTEIN, RNA
 from cogent3.core.sequence import ArraySequence, RnaSequence, Sequence
 from cogent3.maths.util import safe_p_log_p
@@ -1167,19 +1168,6 @@ class SequenceCollectionBaseTests(object):
         with self.assertRaises(AssertionError):
             seqs.dotplot(window=5, k=11)
 
-    def test_dotplot_annotated(self):
-        """exercising dotplot method with annotated sequences"""
-        seqs = self.Class(data={"Human": "CAGATTTGGCAGTT-", "Mouse": "CAGATTCAGCAGGTG"})
-
-        seqs = seqs.take_seqs(["Human", "Mouse"])
-
-        if type(self.Class) != ArrayAlignment:
-            # we annotated Human
-            seq = seqs.get_seq("Human")
-            _ = seq.add_feature(biotype="exon", name="fred", spans=[(10, 15)])
-
-        _ = seqs.dotplot(show_progress=False)
-
     def test_rename_seqs(self):
         """successfully rename sequences"""
         data = {"seq1": "ACGTACGTA", "seq2": "ACCGAA---", "seq3": "ACGTACGTT"}
@@ -2301,14 +2289,6 @@ class AlignmentBaseTests(SequenceCollectionBaseTests):
         self.assertTrue("-" not in found_motifs)
         self.assertEqual(lengths, {2})
 
-    def test_get_seq_entropy(self):
-        """ArrayAlignment get_seq_entropy should get entropy of each seq"""
-        seqs = [AB.make_seq(s, preserve_case=True) for s in ["abab", "bbbb", "abbb"]]
-        a = self.Class(seqs, alphabet=AB.alphabet)
-        entropy = a.entropy_per_seq()
-        e = 0.81127812445913283  # sum(p log_2 p) for p = 0.25, 0.75
-        assert_allclose(entropy, array([1, 0, e]))
-
     def test_seq_entropy_just_gaps(self):
         """ArrayAlignment get_seq_entropy should get entropy of each seq"""
         a = self.Class(dict(a="A---", b="----"), moltype=DNA)
@@ -2777,11 +2757,6 @@ class AlignmentTests(AlignmentBaseTests, TestCase):
         for name in rna.names:
             orig_seq = aln.get_seq(name)
             new_seq = rna.get_seq(name)
-            self.assertEqual(len(orig_seq.annotations), len(new_seq.annotations))
-            for src, dest in zip(orig_seq.annotations, new_seq.annotations):
-                self.assertEqual(src.get_coordinates(), dest.get_coordinates())
-                self.assertIsInstance(src, dest.__class__)
-                self.assertIs(dest.parent, new_seq)
         # check the sequence moltypes
         self.assertEqual({s.data.moltype.label for s in rna.seqs}, {"rna"})
         self.assertEqual(rna.moltype.label, "rna")
@@ -3220,9 +3195,11 @@ def test_get_annotations_from_any_seq(cls):
     """get_annotations_from_any_seq returns correct annotations"""
     data = {"seq1": "ACGTACGTA", "seq2": "ACCGAA---", "seq3": "ACGTACGTT"}
     seqs = cls(data, moltype=DNA)
-    seqs.get_seq("seq1").add_feature(biotype="exon", name="annotation1", spans=[(3, 8)])
-    seqs.get_seq("seq2").add_feature(biotype="exon", name="annotation2", spans=[(1, 2)])
-    seqs.get_seq("seq3").add_feature(biotype="exon", name="annotation3", spans=[(3, 6)])
+    db = GffAnnotationDb()
+    db.add_feature(seqid="seq1", biotype="exon", name="annotation1", spans=[(3, 8)])
+    db.add_feature(seqid="seq2", biotype="exon", name="annotation2", spans=[(1, 2)])
+    db.add_feature(seqid="seq3", biotype="exon", name="annotation3", spans=[(3, 6)])
+    seqs.annotation_db = db
     got = list(seqs.get_features())
     assert len(got) == 3
     assert 'exon "annotation1" at [3:8]/9' in str(got[0])
@@ -3246,6 +3223,7 @@ def test_annotate_matches_to():
     """Aligned.annotate_matches_to correctly delegates to sequence"""
 
     aln = Alignment(dict(x="TTCCACTTCCGCTT"), moltype="dna")
+    aln.annotation_db = GffAnnotationDb()
     seq = aln.named_seqs["x"]
     pattern = "CCRC"
     annot = seq.annotate_matches_to(
@@ -3262,6 +3240,7 @@ def test_annotate_matches_to():
 
     # handles regex from aa
     aln = Alignment(dict(x="TTCCACTTCCGCTT"), moltype="dna")
+    aln.annotation_db = GffAnnotationDb()
     gc = get_code(1)
     aa_regex = gc.to_regex("FHF")
     s = aln.named_seqs["x"].annotate_matches_to(
@@ -3320,59 +3299,37 @@ def test_copy_annotations_incompat_type_fails(seqcoll_db, gb_db):
         seqcoll_db.copy_annotations({"a": "ACGGT"})
 
 
-@pytest.mark.xfail(reason="todo gah to be implemented")
 def test_deepcopy_with_features():
     """correctly deepcopy Aligned objects in an alignment"""
     path = DATA_DIR / "brca1_5.paml"
     # generates an annotatable Alignment object
     aln = load_aligned_seqs(path, array_align=False, moltype="dna")
+    db = GffAnnotationDb()
     # when the annotation is outside(before) boundary of the slice
-    aln.add_feature(seqid="NineBande", biotype="exon", name="annot1", spans=[(0, 10)])
+    db.add_feature(seqid="NineBande", biotype="exon", name="annot1", spans=[(0, 10)])
     # when the annotation is across boundary of the slice
-    aln.add_feature(seqid="Mouse", biotype="exon", name="annot2", spans=[(10, 21)])
+    db.add_feature(seqid="Mouse", biotype="exon", name="annot2", spans=[(10, 21)])
     # when the annotation is within boundary of the slice
-    aln.add_feature(seqid="Human", biotype="exon", name="annot3", spans=[(20, 25)])
+    db.add_feature(seqid="Human", biotype="exon", name="annot3", spans=[(20, 25)])
     # when the annotation is across boundary of the slice
-    aln.add_feature(seqid="HowlerMon", biotype="exon", name="annot4", spans=[(25, 32)])
+    db.add_feature(seqid="HowlerMon", biotype="exon", name="annot4", spans=[(25, 32)])
     # when the annotation is outside(after) boundary of the slice
-    aln.add_feature(seqid="DogFaced", biotype="exon", name="annot5", spans=[(40, 45)])
+    db.add_feature(seqid="DogFaced", biotype="exon", name="annot5", spans=[(40, 45)])
+    aln.annotation_db = db
     aln = aln[20:30]
-
-    # for these species, each has an annotation spanning slice boundary or within it
-    for name in ["Mouse", "Human", "HowlerMon"]:
-        new_seq = aln.named_seqs[name].deepcopy(sliced=True)
-        seq = aln.named_seqs[name]
-        assert new_seq.map.parent_length, seq.map.parent_length
-
-        assert len(new_seq.data) == 10
-        assert new_seq.data.is_annotated()
-        assert len(new_seq.data.annotations) == 1
-        # tests the case when sliced argument if False
-        new_seq = aln.named_seqs[name].deepcopy(sliced=False)
-        assert new_seq.map.parent_length == seq.map.parent_length
-        assert len(new_seq.data) == len(aln.named_seqs[name].data)
-        assert new_seq.data.is_annotated()
-
-    # for these species, each has an annotation outside slice
-    for name in ["NineBande", "DogFaced"]:
-        new_seq = aln.named_seqs[name].deepcopy(sliced=True)
-        assert len(new_seq.data) == 10
-        assert new_seq.data.is_annotated()
-        # tests the case when sliced argument if False
-        new_seq = aln.named_seqs[name].deepcopy(sliced=False)
-        assert len(new_seq.data) == len(aln.named_seqs[name].data)
-        assert new_seq.data.is_annotated()
-        assert len(new_seq.data.annotations) == 1
-
-    # add another human annotation that is outside slice
-    aln.named_seqs["Human"].data.add_feature(
-        biotype="exon", name="annot6", spans=[(40, 45)]
-    )
-    # tests the case when sliced argument if False regarding the Human sequence
-    new_seq = aln.named_seqs["Human"].deepcopy(sliced=False)
-    assert len(new_seq.data) == len(aln.named_seqs["Human"].data)
-    assert new_seq.data.is_annotated()
-    assert len(new_seq.data.annotations) == 2
+    # no slice
+    copied = aln.deepcopy(sliced=False)
+    feats = list(copied.get_features(biotype="exon", allow_partial=True))
+    assert len(feats) == 3  # overlap is Mouse, Human, HowlerMon
+    # with slice
+    copied = aln.deepcopy(sliced=True)
+    feats = list(copied.get_features(biotype="exon", allow_partial=True))
+    assert len(feats) == 3
+    # rc drops annotations only when sliced is True
+    rced = aln.rc()
+    copied = rced.deepcopy(sliced=False)
+    feats = list(copied.get_features(biotype="exon", allow_partial=True))
+    assert len(feats) == 3  # overlap is Mouse, Human, HowlerMon
 
 
 @pytest.mark.parametrize("cls", (SequenceCollection, Alignment))
@@ -3449,3 +3406,25 @@ def test_to_rich_dict(cls, with_offset):
         "version": __version__,
     }
     assert got == expect
+
+
+@pytest.mark.parametrize("cls", (SequenceCollection, Alignment))
+def test_dotplot_annotated(cls):
+    """exercising dotplot method with annotated sequences"""
+    db = GffAnnotationDb()
+    db.add_feature(seqid="Human", biotype="exon", name="fred", spans=[(10, 15)])
+
+    seqs = cls(data={"Human": "CAGATTTGGCAGTT-", "Mouse": "CAGATTCAGCAGGTG"})
+    seqs.annotation_db = db
+    seqs = seqs.take_seqs(["Human", "Mouse"])
+    _ = seqs.dotplot(show_progress=False)
+
+
+@pytest.mark.parametrize("cls", (Alignment, ArrayAlignment))
+def test_get_seq_entropy(cls):
+    """ArrayAlignment get_seq_entropy should get entropy of each seq"""
+    seqs = [AB.make_seq(s, preserve_case=True) for s in ["abab", "bbbb", "abbb"]]
+    a = cls(seqs, alphabet=AB.alphabet)
+    entropy = a.entropy_per_seq()
+    e = 0.81127812445913283  # sum(p log_2 p) for p = 0.25, 0.75
+    assert_allclose(entropy, array([1, 0, e]))
diff --git a/tests/test_core/test_annotation.py b/tests/test_core/test_annotation.py
index 88abea06e..2c392b5aa 100644
--- a/tests/test_core/test_annotation.py
+++ b/tests/test_core/test_annotation.py
@@ -91,20 +91,6 @@ class TestAnnotations(unittest.TestCase):
             expected = seq_expecteds[annot_type]
             assert observed == expected
 
-    def test_feature_copy_annotations_to(self):
-        """test correct copy of annotations"""
-        orig = DnaSequence("TTTTTTTTTTAAAA", name="Orig")
-        annot = orig.add_annotation(Feature, "exon", "fred", [(0, 14)])
-        seq = RnaSequence("UUUUUUUUUUAAAA", name="Test")
-        got = annot.copy_annotations_to(seq)
-        self.assertEqual(len(orig.annotations), len(got.annotations))
-        for src, dest in zip(orig.annotations, got.annotations):
-            self.assertEqual(src.get_coordinates(), dest.get_coordinates())
-            self.assertIsInstance(src, dest.__class__)
-            self.assertIs(dest.parent, seq)
-        with self.assertRaises(AssertionError):
-            _ = annot.copy_annotations_to(seq[:-2])
-
 
 class TestMapSpans(unittest.TestCase):
     """Test attributes of Map & Spans classes critical to annotation
diff --git a/tests/test_core/test_annotation_db.py b/tests/test_core/test_annotation_db.py
index 3b6e5469c..bd366e4fc 100644
--- a/tests/test_core/test_annotation_db.py
+++ b/tests/test_core/test_annotation_db.py
@@ -222,10 +222,8 @@ def test_get_features_matching_matching_feature(seq, anno_db):
     assert len(got) == 1
 
 
-def test_gav():
-    # gah todo also need to check a slice with stride
+def test_feature_strand():
     from cogent3 import make_seq
-    from cogent3.core.annotation import Feature
 
     #            ++   ++++
     #              ---    --
@@ -261,24 +259,15 @@ def test_gav():
     minus = list(seq.get_features(name="minus"))[0]
     assert str(minus.get_slice()) == minus_seq
 
-    f = Feature(seq, "cds", "minus", reversed([reversed(s) for s in minus_spans]))
-    s = f.get_slice()
-    assert str(s) == minus_seq
-
     # now reverse complement the sequence
     rced = seq.rc()
-
-    rced_spans = [(5, 9), (12, 14)]
-    rf = Feature(rced, "cds", "minus", [reversed(s) for s in rced_spans])
-    got = rf.get_slice()
-
     plus = list(rced.get_features(name="plus"))[0]
     assert str(plus.get_slice()) == plus_seq
     minus = list(rced.get_features(name="minus"))[0]
     assert str(minus.get_slice()) == minus_seq
 
 
-def test_gav2():
+def test_feature_nucleic():
     from cogent3 import make_seq
     from cogent3.core import location as loc
 
diff --git a/tests/test_core/test_features.py b/tests/test_core/test_features.py
index 7bdb55c34..e0f8baea2 100755
--- a/tests/test_core/test_features.py
+++ b/tests/test_core/test_features.py
@@ -2,10 +2,8 @@ from unittest import TestCase
 
 import pytest
 
-from cogent3 import ASCII, DNA, make_aligned_seqs
-from cogent3.core.annotation import Feature, Variable
+from cogent3 import ASCII, DNA, get_moltype, make_aligned_seqs
 from cogent3.core.annotation_db import GffAnnotationDb
-
 # Complete version of manipulating sequence annotations
 from cogent3.util.deserialise import deserialise_object
 
@@ -83,200 +81,6 @@ class FeaturesTest(TestCase):
         shadow = exon1.union(exons).shadow()
         assert str(shadow.get_slice()) == expect
 
-    @pytest.mark.xfail(reason="todo gah update test to use latest API")
-    def test_feature_reverse(self):
-        """reverse complement of features"""
-
-        # When dealing with sequences that can be reverse complemented
-        # (e.g. DnaSequence) features are **not** reversed.
-        # Features are considered to have strand specific meaning
-        # (.e.g CDS, exons) and so stay on their original strands.
-        # We create a sequence with a CDS that spans multiple exons,
-        # and show that after getting the reverse complement we have
-        # exactly the same result from getting the CDS annotation.
-
-        plus = DNA.make_seq("AAGGGGAAAACCCCCAAAAAAAAAATTTTTTTTTTAAA", name="plus")
-        plus_cds = plus.add_annotation(
-            Feature, "CDS", "gene", [(2, 6), (10, 15), (25, 35)]
-        )
-        self.assertEqual(str(plus_cds.get_slice()), "GGGGCCCCCTTTTTTTTTT")
-        minus = plus.rc()
-        minus_cds = minus.get_features(biotype="CDS")[0]
-        self.assertEqual(str(minus_cds.get_slice()), "GGGGCCCCCTTTTTTTTTT")
-
-    @pytest.mark.xfail(reason="todo gah implement lost spans")
-    def test_lost_spans(self):
-        """features no longer included in an alignment represented by lost spans"""
-
-        # If the feature lies outside the sequence being copied to, you get a
-        # lost span
-
-        aln = make_aligned_seqs(
-            data=[["x", "-AAAA"], ["y", "TTTTT"]], array_align=False
-        )
-        seq = DNA.make_seq("CCCCCCCCCCCCCCCCCCCC", "x")
-        seq.add_feature("exon", "A", [(5, 8)])
-        aln.get_seq("x").copy_annotations(seq)
-        copied = list(aln.get_features(seqid="x", name="exon"))
-        self.assertEqual(str(copied), '[exon "A" at [5:5, -4-]/5]')
-        self.assertEqual(str(copied[0].get_slice()), ">x\n----\n>y\n----\n")
-
-    @pytest.mark.xfail(reason="todo gah update test to use latest API")
-    def test_annotated_separately_equivalence(self):
-        """allow defining features as a series or individually"""
-
-        # It shouldn't matter whether annotated coordinates are entered
-        # separately, or as a series.
-
-        data = [["human", "CGAAACGTTT"], ["mouse", "CTAAACGTCG"]]
-        as_series = make_aligned_seqs(data=data, array_align=False)
-        as_items = make_aligned_seqs(data=data, array_align=False)
-
-        # We add annotations to the sequences as a series.
-
-        self.assertEqual(
-            str(
-                as_series.get_seq("human").add_feature(
-                    biotype="cpgsite", name="cpg", spans=[(0, 2), (5, 7)]
-                )
-            ),
-            'cpgsite "cpg" at [0:2, 5:7]/10',
-        )
-        self.assertEqual(
-            str(
-                as_series.get_seq("mouse").add_feature(
-                    biotype="cpgsite", name="cpg", spans=[(5, 7), (8, 10)]
-                )
-            ),
-            'cpgsite "cpg" at [5:7, 8:10]/10',
-        )
-
-        # We add the annotations to the sequences one segment at a time.
-
-        self.assertEqual(
-            str(as_items.get_seq("human").add_feature("cpgsite", "cpg", [(0, 2)])),
-            'cpgsite "cpg" at [0:2]/10',
-        )
-        self.assertEqual(
-            str(as_items.get_seq("human").add_feature("cpgsite", "cpg", [(5, 7)])),
-            'cpgsite "cpg" at [5:7]/10',
-        )
-        self.assertEqual(
-            str(as_items.get_seq("mouse").add_feature("cpgsite", "cpg", [(5, 7)])),
-            'cpgsite "cpg" at [5:7]/10',
-        )
-        self.assertEqual(
-            str(as_items.get_seq("mouse").add_feature("cpgsite", "cpg", [(8, 10)])),
-            'cpgsite "cpg" at [8:10]/10',
-        )
-
-    @pytest.mark.xfail(
-        reason="todo gah update test to use latest API plus support serialisation"
-    )
-    def test_roundtripped_alignment(self):
-        """Alignment with annotations roundtrips correctly"""
-        # annotations just on member sequences
-        aln = make_aligned_seqs(
-            data=[["x", "-AAAAAAAAA"], ["y", "TTTT--TTTT"]], array_align=False
-        )
-        _ = aln.get_seq("x").add_annotation(Feature, "exon", "fred", [(3, 8)])
-        seq_exon = list(aln.get_features(seqid="x", name="exon"))[0]
-        expect = seq_exon.get_slice()
-
-        json = aln.to_json()
-        new = deserialise_object(json)
-        got_exons = list(new.get_features(seqid="x", name="exon"))[0]
-        self.assertEqual(got_exons.get_slice().to_dict(), expect.to_dict())
-
-        # annotations just on alignment
-        aln = make_aligned_seqs(
-            data=[["x", "-AAAAAGGGG"], ["y", "TTTT--CCCC"]], array_align=False
-        )
-        f = aln.add_annotation(Feature, "generic", "no name", [(1, 4), (6, 10)])
-        expect = f.get_slice().to_dict()
-        json = aln.to_json()
-        new = deserialise_object(json)
-        got = list(new.get_features(biotype="generic"))[0]
-        self.assertEqual(got.get_slice().to_dict(), expect)
-
-        # annotations on both alignment and sequence
-        aln = make_aligned_seqs(
-            data=[["x", "-AAAAAGGGG"], ["y", "TTTT--CCCC"]], array_align=False
-        )
-        f = aln.add_annotation(Feature, "generic", "no name", [(1, 4), (6, 10)])
-        _ = aln.get_seq("x").add_annotation(Feature, "exon", "1", [(3, 8)])
-        json = aln.to_json()
-        new = deserialise_object(json)
-        ## get back the exon
-        seq_exon = list(aln.get_features(seqid="x", name="exon"))[0]
-        expect = seq_exon.get_slice().to_dict()
-        got_exons = list(new.get_features(seqid="x", name="exon"))[0]
-        self.assertEqual(got_exons.get_slice().to_dict(), expect)
-        ## get back the generic
-        expect = f.get_slice().to_dict()
-        got = list(new.get_features(biotype="generic"))[0]
-        self.assertEqual(got.get_slice().to_dict(), expect)
-
-        # check masking of seq features still works
-        new = new.with_masked_annotations("exon", mask_char="?")
-        self.assertEqual(new[4:9].to_dict(), dict(x="?????", y="--CCC"))
-
-    @pytest.mark.xfail(reason="todo gah update test to use latest API")
-    def test_roundtripped_alignment2(self):
-        """Sliced Alignment with annotations roundtrips correctly"""
-        # annotations just on member sequences
-        aln = make_aligned_seqs(
-            data=[["x", "-AAAGGGGGAACCCT"], ["y", "TTTT--TTTTAGGGA"]], array_align=False
-        )
-        aln.get_seq("x").add_annotation(Feature, "exon", "E1", [(3, 8)])
-        aln.get_seq("x").add_annotation(Feature, "exon", "E2", [(10, 13)])
-        # at the alignment level
-        sub_aln = aln[:-3]
-        s = sub_aln.named_seqs["x"]
-        d = s.data[:11]
-        json = s.to_json()
-        new = deserialise_object(json)
-        gf1, gf2 = list(new.data.get_features(biotype="exon"))
-        self.assertEqual(str(gf1.get_slice()), "GGGGG")
-        self.assertEqual(str(gf2.get_slice()), "C")
-        # the sliced alignment
-        json = sub_aln.to_json()
-        got = deserialise_object(json)
-        x = got.named_seqs["x"]
-        self.assertEqual(str(x.data.annotations[0].get_slice()), "GGGGG")
-        self.assertEqual(str(x.data.annotations[1].get_slice()), "C")
-
-    @pytest.mark.xfail(reason="todo gah delete test not supporting Variable class")
-    def test_roundtrip_variable(self):
-        """should recover the Variable feature type"""
-        seq = DNA.make_seq("AAGGGGAAAACCCCCAAAAAAAAAATTTTTTTTTTAAA", name="plus")
-        xx_y = [[[2, 6], 2.4], [[10, 15], 5.1], [[25, 35], 1.3]]
-        y_valued = seq.add_annotation(Variable, "SNP", "freq", xx_y)
-        json = seq.to_json()
-        new = deserialise_object(json)
-        got = list(new.get_features(biotype="SNP"))[0]
-        # annoyingly, comes back as list of lists
-        self.assertEqual(got.xxy_list, [[list(xx), y] for xx, y in y_valued.xxy_list])
-
-    @pytest.mark.xfail(reason="todo gah change test to use latest API")
-    def test_nested_to_rich_dict(self):
-        """check the to_rich_dict method works with nested annotations"""
-        self.assertEqual(
-            self.exon1.to_rich_dict()["annotations"][0],
-            self.nested_feature.to_rich_dict(),
-        )
-
-    @pytest.mark.xfail(reason="todo gah update test to use latest API")
-    def test_nested_deserialise_annotation(self):
-        """nested annotations can be deserialised"""
-        got = self.s.to_json()
-        new = deserialise_object(got)
-        new_exon1 = new.annotations[0]
-        new_nested_feature = new_exon1.annotations[0]
-        self.assertEqual(
-            new_nested_feature.to_rich_dict(), self.nested_feature.to_rich_dict()
-        )
-
     def test_annotate_matches_to(self):
         """annotate_matches_to attaches annotations correctly to a Sequence"""
         seq = DNA.make_seq("TTCCACTTCCGCTT", name="x")
@@ -322,10 +126,10 @@ def test_feature_residue():
         moltype=DNA,
         array_align=False,
     )
+    db = GffAnnotationDb()
+    aln.annotation_db = db
     assert str(aln), ">x\nC-CCCAAAAA\n>y\n-T----TTTT\n"
-    exon = aln.get_seq("x").add_feature(biotype="exon", name="ex1", spans=[(0, 4)])
-    assert 'exon "ex1" at [0:4]/9' in str(exon)
-    assert str(exon.get_slice()), "CCCC"
+    db.add_feature(seqid="x", biotype="exon", name="ex1", spans=[(0, 4)])
     aln_exons = list(aln.get_features(seqid="x", biotype="exon"))
     assert 'exon "ex1" at [0:1, 2:5]/10' in str(aln_exons)
     assert aln_exons[0].get_slice().to_dict() == dict(x="CCCC", y="----")
@@ -577,9 +381,6 @@ def test_annotated_region_masks():
 
 def test_nested_annotated_region_masks():
     """masking a sequence with specific features when nested annotations"""
-    # gene = aln.get_seq("x").add_feature("gene", "norwegian", [(0, 4)])
-    # gene.add_feature("repeat", "blue", [(1, 3)])
-    # exon = aln.get_seq("y").add_feature("repeat", "frog", [(1, 4)])
     db = GffAnnotationDb()
     db.add_feature(seqid="x", biotype="gene", name="norwegian", spans=[(0, 4)])
     db.add_feature(seqid="x", biotype="repeat", name="blue", spans=[(1, 3)])
@@ -747,3 +548,121 @@ def test_roundtrip_json():
     assert new.annotation_offset == 3
     feat = list(new.get_features(biotype="exon", allow_partial=True))[0]
     assert str(feat.get_slice()) == "AA"
+
+
+def test_roundtripped_alignment():
+    """Alignment with annotations roundtrips correctly"""
+    # annotations just on member sequences
+    aln = make_aligned_seqs(
+        data=[["x", "-AAAAAAAAA"], ["y", "TTTT--TTTT"]], array_align=False
+    )
+    db = GffAnnotationDb()
+    db.add_feature(seqid="x", biotype="exon", name="fred", spans=[(3, 8)])
+    aln.annotation_db = db
+    seq_exon = list(aln.get_features(seqid="x", biotype="exon"))[0]
+    expect = seq_exon.get_slice()
+
+    json = aln.to_json()
+    new = deserialise_object(json)
+    got_exons = list(new.get_features(seqid="x", biotype="exon"))[0]
+    assert got_exons.get_slice().to_dict() == expect.to_dict()
+
+    # annotations just on alignment
+    aln = make_aligned_seqs(
+        data=[["x", "-AAAAAGGGG"], ["y", "TTTT--CCCC"]], array_align=False
+    )
+
+    f = aln.add_feature(biotype="generic", name="no name", spans=[(1, 4), (6, 10)])
+    expect = f.get_slice().to_dict()
+    json = aln.to_json()
+    new = deserialise_object(json)
+    got = list(new.get_features(biotype="generic"))[0]
+    assert got.get_slice().to_dict() == expect
+    # annotations on both alignment and sequence
+    aln = make_aligned_seqs(
+        data=[["x", "-AAAAAGGGG"], ["y", "TTTT--CCCC"]], array_align=False
+    )
+    db = GffAnnotationDb()
+    db.add_feature(
+        seqid=None,
+        biotype="generic",
+        name="no name",
+        spans=[(1, 4), (6, 10)],
+        on_alignment=True,
+    )
+    db.add_feature(seqid="x", biotype="exon", name="1", spans=[(3, 8)])
+    aln.annotation_db = db
+    json = aln.to_json()
+    new = deserialise_object(json)
+    ## get back the exon
+    seq_exon = list(aln.get_features(seqid="x", biotype="exon"))[0]
+    expect = seq_exon.get_slice().to_dict()
+    got_exons = list(new.get_features(seqid="x", biotype="exon"))[0]
+    assert got_exons.get_slice().to_dict() == expect
+    ## get back the generic
+    expect = f.get_slice().to_dict()
+    got = list(new.get_features(biotype="generic"))[0]
+    assert got.get_slice().to_dict() == expect
+
+    # check masking of seq features still works
+    new = new.with_masked_annotations("exon", mask_char="?")
+    assert new[4:9].to_dict() == dict(x="?????", y="--CCC")
+
+
+def test_feature_out_range():
+    """features no longer included in an alignment will not be returned"""
+    aln = make_aligned_seqs(data=[["x", "-AAAA"], ["y", "TTTTT"]], array_align=False)
+    db = GffAnnotationDb()
+    db.add_feature(seqid="x", biotype="exon", name="A", spans=[(5, 8)])
+    f = list(aln.get_features(seqid="x", biotype="exon"))
+    assert not f
+
+
+def test_roundtripped_alignment_with_slices():
+    """Sliced Alignment with annotations roundtrips correctly"""
+    # annotations just on member sequences
+    aln = make_aligned_seqs(
+        data=[["x", "-AAAGGGGGAACCCT"], ["y", "TTTT--TTTTAGGGA"]], array_align=False
+    )
+    db = GffAnnotationDb()
+    db.add_feature(seqid="x", biotype="exon", name="E1", spans=[(3, 8)])
+    db.add_feature(seqid="x", biotype="exon", name="E2", spans=[(10, 13)])
+    aln.annotation_db = db
+    # at the alignment level
+    sub_aln = aln[:-3]
+    new = deserialise_object(sub_aln.to_json())
+    gf1, gf2 = list(new.get_features(biotype="exon", allow_partial=True))
+    assert gf1.get_slice().to_dict() == {"x": "GGGGG", "y": "--TTT"}
+    assert gf2.get_slice().to_dict() == {"x": "C", "y": "G"}
+
+
+def test_feature_reverse():
+    """reverse complement of features"""
+
+    # When dealing with sequences that can be reverse complemented
+    # (e.g. DnaSequence) features are **not** reversed.
+    # Features are considered to have strand specific meaning
+    # (.e.g CDS, exons) and so stay on their original strands.
+    # We create a sequence with a CDS that spans multiple exons,
+    # and show that after getting the reverse complement we have
+    # exactly the same result from getting the CDS annotation.
+
+    plus = DNA.make_seq("AAGGGGAAAACCCCCAAAAAAAAAATTTTTTTTTTAAA", name="plus")
+    plus_cds = plus.add_feature(
+        biotype="CDS", name="gene", spans=[(2, 6), (10, 15), (25, 35)]
+    )
+    assert str(plus_cds.get_slice()) == "GGGGCCCCCTTTTTTTTTT"
+    minus = plus.rc()
+    minus_cds = list(minus.get_features(biotype="CDS"))[0]
+    assert str(minus_cds.get_slice()) == "GGGGCCCCCTTTTTTTTTT"
+
+
+@pytest.mark.parametrize("moltype", ("protein", "bytes", "text"))
+def test_rc_feature_on_wrong_moltype(moltype):
+    moltype = get_moltype(moltype)
+    seq = moltype.make_seq("AAGGGGAAAACCCCCAAAAAAAAAATTTTTTTTTTAAA", name="s1")
+    cds = seq.add_feature(
+        biotype="CDS", name="gene", spans=[(2, 6), (10, 15), (25, 35)], strand="-"
+    )
+    with pytest.raises(TypeError):
+        cds.get_slice()
diff --git a/tests/test_core/test_sequence.py b/tests/test_core/test_sequence.py
index 958380701..515f76ddd 100644
--- a/tests/test_core/test_sequence.py
+++ b/tests/test_core/test_sequence.py
@@ -162,38 +162,13 @@ class SequenceTests(TestCase):
     def test_sequence_to_moltype(self):
         """correctly convert to specified moltype"""
         s = Sequence("TTTTTTTTTTAAAA", name="test1")
-        annot1 = s.add_annotation(Feature, "exon", "fred", [(0, 10)])
-        annot2 = s.add_annotation(Feature, "exon", "trev", [(10, 14)])
+        s.add_feature(biotype="exon", name="fred", spans=[(0, 10)])
+        s.add_feature(biotype="exon", name="trev", spans=[(10, 14)])
         got = s.to_moltype("rna")
-        annot1_slice = str(annot1.get_slice())
-        annot2_slice = str(annot2.get_slice())
-        got1_slice = str(got.annotations[0].get_slice())
-        got2_slice = str(got.annotations[1].get_slice())
-        self.assertNotEqual(annot1_slice, got1_slice)
-        self.assertEqual(annot2_slice, got2_slice)
-        self.assertEqual(got.moltype.label, "rna")
-        self.assertEqual(got.name, "test1")
-
-        s = Sequence("AAGGGGAAAACCCCCAAAAAAAAAATTTTTTTTTTAAA", name="test2")
-        xx_y = [[[2, 6], 2.4], [[10, 15], 5.1], [[25, 35], 1.3]]
-        y_valued = s.add_annotation(Variable, "SNP", "freq", xx_y)
-        got = s.to_moltype("rna")
-        y_valued_slice = str(y_valued.get_slice())
-        got_slice = str(str(got.annotations[0].get_slice()))
-        self.assertNotEqual(y_valued_slice, got_slice)
-        self.assertEqual(got.moltype.label, "rna")
-        self.assertEqual(got.name, "test2")
-
-        s = Sequence("TTTTTTTTTTAAAAAAAAAA", name="test3")
-        data = [i for i in range(20)]
-        annot4 = s.add_annotation(SimpleVariable, "SNP", "freq", data)
-        got = s.to_moltype(RNA)
-        annot4_slice = str(annot4.get_slice())
-        got_slice = str(str(got.annotations[0].get_slice()))
-        self.assertNotEqual(annot4_slice[:10], got_slice[:10])
-        self.assertEqual(annot4_slice[10:20], got_slice[10:20])
-        self.assertEqual(got.moltype.label, "rna")
-        self.assertEqual(got.name, "test3")
+        fred = list(got.get_features(name="fred"))[0]
+        assert str(got[fred]) == "UUUUUUUUUU"
+        trev = list(got.get_features(name="trev"))[0]
+        assert str(got[trev]) == "AAAA"
 
         # calling with a null object should raise an exception
         with self.assertRaises(ValueError):
@@ -1922,21 +1897,14 @@ def test_seq_repr(one_seq, rc):
     assert expect.startswith(got), (expect, got)
 
 
-@pytest.mark.xfail(reason="todo: to be implemented")
-def test_db_bind_fails_with_incompatible_moltype(worm_seq_path, worm_gff_path):
-    from cogent3.core.annotation_db import load_annotations
-
-    seq = cogent3.load_seq(worm_seq_path, moltype="text")
-    db = load_annotations(worm_gff_path)
-    with pytest.raises(TypeError):
-        seq.annotation_db = db
-
-
-@pytest.mark.xfail(reason="todo: write test")
 def test_annotation_from_slice_with_stride():
-    raise NotImplementedError(
-        "check annotations get correct sequence when there's a stride"
-    )
+    seq = DNA.make_seq("AAACGCGCGAAAAAAA", name="s1")
+    seq.add_feature(biotype="exon", name="ex1", spans=[(3, 9)])
+    f = list(seq.get_features(name="ex1"))[0]
+    assert str(f.get_slice()) == "CGCGCG"
+    s1 = seq[1::2]
+    f = list(s1.get_features(name="ex1"))[0]
+    assert str(f.get_slice()) == "CCC"
 
 
 def test_absolute_position_base_cases(one_seq):
@@ -2040,7 +2008,7 @@ def test_relative_position_with_remainder(integer_seq):
 def test_absolute_relative_roundtrip(one_seq, value, offset, start, stop, step):
     # a round trip from relative to absolute then from absolute to relative, should return the same value we began with
     view = one_seq[start:stop:step]
-    view.annotation_offset = offset
+    view.annotation_offset = offset or 0
     abs_val = view._seq.absolute_position(value)
     rel_val = view._seq.relative_position(abs_val)
     assert rel_val == value
@@ -2056,10 +2024,10 @@ def test_absolute_relative_roundtrip_reverse(
 ):
     # a round trip from relative to absolute then from absolute to relative, should return the same value we began with
     view = integer_seq[start:stop:step]
-    view.offset = offset
+    view.offset = offset or 0
     abs_val = view.absolute_position(value)
     rel_val = view.relative_position(abs_val)
-    assert view.offset == offset
+    assert view.offset == (offset or 0)
     assert (view[rel_val]).value == view[value].value
 
 
@@ -2158,3 +2126,20 @@ def test_to_json(cls, with_offset):
         expect["annotation_offset"] = 0
 
     assert got == expect
+
+
+def test_offset_with_multiple_slices():
+    from cogent3.util.deserialise import deserialise_object
+
+    seq = DNA.make_seq("ACCCCGGAAAATTTTTTTTTAAGGGGGAAAAAAAAACCCCCCC", name="22")
+    gff3_path = DATADIR / "ensembl_sample.gff3"
+    seq.annotate_from_gff(gff3_path)
+    rd = seq[2:].to_rich_dict()
+    s1 = deserialise_object(rd)
+    assert s1.annotation_offset == 2
+    rd = s1[3:].to_rich_dict()
+    s2 = deserialise_object(rd)
+    assert s2.annotation_offset == 5
+    expect = {(f.seqid, f.biotype, f.name) for f in seq.get_features(start=5)}
+    got = {(f.seqid, f.biotype, f.name) for f in s2.get_features()}
+    assert got == expect
diff --git a/tests/test_draw/test_draw_integration.py b/tests/test_draw/test_draw_integration.py
index d99005616..787aaf587 100644
--- a/tests/test_draw/test_draw_integration.py
+++ b/tests/test_draw/test_draw_integration.py
@@ -4,6 +4,7 @@ import unittest
 from numpy.testing import assert_allclose
 
 from cogent3 import load_aligned_seqs, make_aligned_seqs, make_table
+from cogent3.core.annotation_db import GffAnnotationDb
 from cogent3.draw.drawable import AnnotatedDrawable, Drawable, get_domain
 from cogent3.util.union_dict import UnionDict
 
@@ -20,15 +21,27 @@ __status__ = "Alpha"
 
 def load_alignment(annotate1=False, annotate2=False):
     """creates an alignment with None, one or two sequences annotated"""
+    db = GffAnnotationDb()
+
     path = str(pathlib.Path(__file__).parent.parent / "data/brca1_5.paml")
     aln = load_aligned_seqs(path, array_align=False, moltype="dna")
     aln = aln.omit_gap_pos()
     if annotate1:
-        aln.get_seq(aln.names[0]).add_feature("gene", "abcde1", [(20, 50)])
-        aln.get_seq(aln.names[0]).add_feature("variation", "one", [(11, 12)])
+        db.add_feature(
+            seqid=aln.names[0], biotype="gene", name="abcde1", spans=[(20, 50)]
+        )
+        db.add_feature(
+            seqid=aln.names[0], biotype="variation", name="one", spans=[(11, 12)]
+        )
     if annotate2:
-        aln.get_seq(aln.names[1]).add_feature("gene", "abcde2", [(20, 50)])
-        aln.get_seq(aln.names[1]).add_feature("domain", "abcde2", [(10, 15)])
+        db.add_feature(
+            seqid=aln.names[1], biotype="gene", name="abcde2", spans=[(20, 50)]
+        )
+
+        db.add_feature(
+            seqid=aln.names[1], biotype="domain", name="abcde2", spans=[(10, 15)]
+        )
+    aln.annotation_db = db
     return aln
 
 
@@ -359,7 +372,7 @@ class AlignmentDrawablesTests(BaseDrawablesTests):
     def test_coevo_annotated(self):
         """coevolution on alignment with annotated seqs should add to heatmap plot"""
         aln = load_alignment(True, False)
-        aln = aln[30:]
+        aln = aln[:30]
         coevo = aln.coevolution(show_progress=False, drawable="heatmap")
         drawable = coevo.drawable
         self.assertIsInstance(drawable, AnnotatedDrawable)
@@ -382,10 +395,12 @@ class AlignmentDrawablesTests(BaseDrawablesTests):
         aln = make_aligned_seqs(
             data=dict(a="AAACGGTTT", b="CAA--GTAA"), array_align=False
         )
-        _ = aln.get_seq("b").add_feature("domain", "1", [(1, 5)])
-        _ = aln.get_seq("b").add_feature("variation", "1", [(1, 5)])
-        _ = aln.get_seq("b").add_feature("gene", "1", [(1, 5)])
-        _ = aln.get_seq("b").add_feature("gene", "1", [(5, 1)])
+        db = GffAnnotationDb()
+        db.add_feature(seqid="b", biotype="domain", name="1", spans=[(1, 5)])
+        db.add_feature(seqid="b", biotype="variation", name="1", spans=[(1, 5)])
+        db.add_feature(seqid="b", biotype="gene", name="1", spans=[(1, 5)])
+        db.add_feature(seqid="b", biotype="gene", name="1", spans=[(5, 1)])
+        aln.annotation_db = db
         aln.get_drawable()
 
 
diff --git a/tests/test_parse/test_cigar.py b/tests/test_parse/test_cigar.py
index 794fbcee4..3ac2bd115 100644
--- a/tests/test_parse/test_cigar.py
+++ b/tests/test_parse/test_cigar.py
@@ -1,7 +1,7 @@
-#!/usr/bin/env python
 import unittest
 
 from cogent3 import DNA, make_aligned_seqs
+from cogent3.core.annotation_db import GffAnnotationDb
 from cogent3.parse.cigar import (
     CigarParser,
     aligned_from_cigar,
@@ -73,12 +73,15 @@ class TestCigar(unittest.TestCase):
         aln = aln.to_type(array_align=False)
         assert aln == self.aln
         # test slice
+        db = GffAnnotationDb()
+        aln.annotation_db = db
         i = 1
         for start, end in self.slices:
-            self.aln.get_seq("FAKE01").add_feature(
-                biotype=f"annot{i}", name="annot", spans=[(start, end)]
+            db.add_feature(
+                seqid="FAKE01", biotype=f"annot{i}", name="annot", spans=[(start, end)]
             )
-            annot = list(self.aln.get_features(biotype=f"annot{i}"))
+
+            annot = list(aln.get_features(biotype=f"annot{i}"))
             annot = annot[0].union(annot[1:])
             slice_aln = annot.as_one_span().get_slice()
             i += 1
diff --git a/tests/test_parse/test_gbseq.py b/tests/test_parse/test_gbseq.py
index 060eb0958..b910a589e 100644
--- a/tests/test_parse/test_gbseq.py
+++ b/tests/test_parse/test_gbseq.py
@@ -1,4 +1,3 @@
-#!/usr/bin/env python
 import xml.dom.minidom
 
 from unittest import TestCase, main
@@ -185,10 +184,4 @@ class ParseGBseq(TestCase):
         ]:
             self.assertEqual(name, "AY286018.1")
             self.assertEqual(sample_seq, seq.to_fasta(block_size=len(sample_seq)))
-            self.assertEqual(str(seq.annotations), sample_annotations)
-
-    pass
-
-
-if __name__ == "__main__":
-    main()
+            self.assertEqual(seq.annotation_db.num_matches(), 4)
diff --git a/tests/test_parse/test_genbank.py b/tests/test_parse/test_genbank.py
index 74ce74196..6e9cf6ef9 100644
--- a/tests/test_parse/test_genbank.py
+++ b/tests/test_parse/test_genbank.py
@@ -1,7 +1,6 @@
-#!/usr/bin/env python
 """Unit tests for the GenBank database parsers.
 """
-from unittest import TestCase, main
+from unittest import TestCase
 
 from cogent3.parse.genbank import (
     Location,
@@ -433,20 +432,11 @@ ORIGIN
 
     def test_rich_parser(self):
         """correctly constructs +/- strand features"""
-        # a custom annotation function
-        from cogent3.core.annotation import Feature
-
-        def add_annotation(seq, feature, spans):
-            if feature["type"] != "CDS":
-                return
-            name = feature["locus_tag"][0]
-            seq.add_annotation(Feature, "CDS", name, spans)
-
-        infile = open("data/annotated_seq.gb")
-        parser = RichGenbankParser(infile, add_annotation=add_annotation)
+        with open("data/annotated_seq.gb") as infile:
+            parser = RichGenbankParser(infile)
+            seq = [s for l, s in parser][0]
 
-        seq = [s for l, s in parser][0]
-        cds = dict([(f.name, f) for f in seq.get_features_matching("CDS")])
+        cds = dict([(f.name, f) for f in seq.get_features(biotype="CDS")])
         expects = {
             "CNA00110": "MAGYDARYGNPLDPMSGGRPSPPETSQQDAYEYSKHGSSSGYLGQLPLGAD"
             "SAQAETASALRTLFGEGADVQALQEPPNQINTLAEGAAVAETGGVLGGDTTRSDNEALAIDPSL"
@@ -469,7 +459,6 @@ ORIGIN
         for locus in cds:
             got = cds[locus].get_slice().trim_stop_codon().get_translation()
             self.assertEqual(str(got), expects[locus])
-        infile.close()
 
     def test_rich_parser_moltype(self):
         """correctly handles moltypes"""
@@ -478,8 +467,9 @@ ORIGIN
             got_1 = [s for _, s in parser][0]
 
         # name formed from /product value
-        got = {f.name for f in got_1.get_features_matching("mRNA")}
-        self.assertEqual(got, {"conserved hypothetical protein", "chaperone, putative"})
+        feature_ids = {"CNA00110", "CNA00120"}
+        got = {f.name for f in got_1.get_features(biotype="mRNA")}
+        self.assertEqual(got, {"CNA00110", "CNA00120"})
 
         # the file defines itself as DNA
         self.assertEqual(got_1.moltype.label, "dna")
@@ -489,13 +479,12 @@ ORIGIN
             with open("data/annotated_seq.gb") as infile:
                 parser = RichGenbankParser(infile, moltype=moltype)
                 got_2 = [s for _, s in parser][0]
-
-            self.assertEqual(len(got_1.annotations), len(got_2.annotations))
-            self.assertEqual(got_2.moltype.label, moltype)
-            got = {f.name for f in got_1.get_features_matching("mRNA")}
             self.assertEqual(
-                got, {"conserved hypothetical protein", "chaperone, putative"}
+                got_1.annotation_db.num_matches(), got_2.annotation_db.num_matches()
             )
+            self.assertEqual(got_2.moltype.label, moltype)
+            got = {f.name for f in got_1.get_features(biotype="mRNA")}
+            self.assertEqual(got, feature_ids)
 
 
 class LocationTests(TestCase):
@@ -539,40 +528,32 @@ def test_Location_start():
 class LocationListTests(TestCase):
     """Tests of the LocationList class."""
 
-    def test_extract(self):
-        """LocationList extract should return correct sequence"""
-        l = Location(3)
-        l2_a = Location(5)
-        l2_b = Location(7)
-        l2 = Location([l2_a, l2_b], strand=-1)
-        l3_a = Location(10)
-        l3_b = Location(12)
-        l3 = Location([l3_a, l3_b])
-        ll = LocationList([l, l2, l3])
-        s = ll.extract("ACGTGCAGTCAGTAGCAT")
-        #               123456789012345678
-        self.assertEqual(s, "G" + "TGC" + "CAG")
-        # check a case where it wraps around
-        l5_a = Location(16)
-        l5_b = Location(4)
-        l5 = Location([l5_a, l5_b])
-        ll = LocationList([l5])
-        s = ll.extract("ACGTGCAGTCAGTAGCAT")
-        self.assertEqual(s, "CATACGT")
-
-
-if __name__ == "__main__":
-    from sys import argv
-
-    if len(argv) > 2 and argv[1] == "x":
-        filename = argv[2]
-        lines = open(filename)
-        for i in indent_splitter(lines):
-            print("******")
-            print(i[0])
-            for j in indent_splitter(i[1:]):
-                print("?????")
-                for line in j:
-                    print(line)
-    else:
-        main()
+
+def test_locationlist_extract():
+    """LocationList extract should return correct sequence"""
+    l = Location(3)
+    l2_a = Location(5)
+    l2_b = Location(7)
+    l2 = Location([l2_a, l2_b], strand=-1)
+    l3_a = Location(10)
+    l3_b = Location(12)
+    l3 = Location([l3_a, l3_b])
+    ll = LocationList([l, l2, l3])
+    s = ll.extract("ACGTGCAGTCAGTAGCAT")
+    #               123456789012345678
+    assert s == "G" + "TGC" + "CAG"
+    # check a case where it wraps around
+    l5_a = Location(16)
+    l5_b = Location(4)
+    l5 = Location([l5_a, l5_b])
+    ll = LocationList([l5])
+    s = ll.extract("ACGTGCAGTCAGTAGCAT")
+    assert s == "CATACGT"
+
+
+def test_location_list_get_coordinates():
+    l = "complement(join(5670..5918,5965..6126))"
+    l = location_line_tokenizer([l])
+    g = parse_location_line(l)
+    spans = g.get_coordinates()
+    assert spans == [(5669, 5918), (5964, 6126)]
diff --git a/tests/test_parse/test_tinyseq.py b/tests/test_parse/test_tinyseq.py
index 26f20018f..7d441162b 100644
--- a/tests/test_parse/test_tinyseq.py
+++ b/tests/test_parse/test_tinyseq.py
@@ -1,4 +1,3 @@
-#!/usr/bin/env python
 import xml.dom.minidom
 
 from unittest import TestCase, main
@@ -46,7 +45,7 @@ class ParseTinyseq(TestCase):
         ]:
             self.assertEqual(name, "AY286018.1")
             self.assertEqual(sample_seq, seq.to_fasta(block_size=len(sample_seq)))
-            self.assertEqual(str(seq.annotations), sample_annotations)
+            self.assertEqual(seq.annotation_db.num_matches(), 2)
 
     pass