25816358406

Committed 13 May 2026 05:46PM UTC coverage: 85.426% (+0.5%) from 84.971%

Build # 25816358406

Build Type

Pull #355

github

Committed by

iskandr

Commit Message

Address PR #355 review: drop pickle schema bump, surface fasta_version

* Reverted FASTA_PICKLE_SCHEMA_VERSION + the .v2.pickle filename suffix.
  v1 (bare-keyed) pickles load cleanly under the new code path - the
  rebuilt _stripped_index stays empty, but lookups still resolve via
  the version-stripped retry path so Ensembl-style callers see no
  regression. New test test_old_bare_keyed_pickle_still_loads_under_new_code
  pins this.

* Added Transcript.fasta_version - returns the int version that the
  cDNA FASTA header carried for this transcript, or None when the
  genome has no transcript FASTA / the header was bare.

* Added Protein.fasta_version - same shape but for the protein FASTA.
  Protein now carries an optional genome= reference so the view object
  can consult genome.protein_sequences.fasta_version(). Default is
  None for back-compat with callers constructing Protein outside of
  Transcript.

* Opted out of gtfparse 2.7.0's cast_version_columns=True in
  database.py's read_gtf call. When the *_version columns are missing
  on some rows (e.g. start_codon rows don't carry transcript_version),
  pandas' nullable Int64 routes through float on the sqlite write path
  and the value lands as '7.0' text, which broke pyensembl's
  int(result[0]) parse. Keeping strings and doing the int conversion
  on the property side as before.

* New tests for both fasta_version accessors, including a "GTF claims
  version 5, FASTA claims version 3" disagreement case so downstream
  tools can detect mismatched GTF/FASTA pairs.

Pull Request Pull Request #355: Preserve FASTA-header versions in SequenceData (closes #351)

Coverage Stats

60 of 61 new or added lines in 5 files covered. (98.36%)

1653 of 1935 relevant lines covered (85.43%)

3.42 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

86.78

/pyensembl/sequence_data.py

# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

from os import remove
from os.path import dirname, exists, abspath, split, join

import datacache
import logging
from collections import Counter
import pickle
from .common import load_pickle, dump_pickle
from .fasta import _split_ens_version, parse_fasta_dictionary


logger = logging.getLogger(__name__)


def lookup_sequence_with_version_fallback(sequence_data, identifier):
    """
    Look up ``identifier`` in ``sequence_data``, tolerating ENS ``.N``
    version-suffix mismatches in either direction.

    Both Ensembl and GENCODE work with versioned IDs (e.g.
    ``ENSP00000123456.3``); the two formats just split that information
    differently. Ensembl GTFs keep the bare ID in ``protein_id`` /
    ``transcript_id`` and put the version in a separate ``*_version``
    attribute, while GENCODE GTFs embed the version directly in the ID
    attribute. Whichever convention the FASTA on disk used, the
    sequence dict now keys on the form that header carried — so either
    a versioned or a bare caller-supplied ID may need a fallback.

    Resolution order:
      1. literal lookup of ``identifier``
      2. if ``identifier`` is a versioned ENS ID, strip ``.N`` and retry
         (handles caller-supplied versioned ID against a bare FASTA)
      3. consult ``sequence_data._stripped_index`` for a versioned alias
         of a bare caller-supplied ID (handles bare ID against a
         versioned FASTA, the GENCODE case)
      4. None

    Non-Ensembl IDs (e.g. TAIR ``AT1G01010.1`` where ``.N`` is an
    isoform suffix, not a version) skip step (2): stripping there is
    semantically wrong.
    """
    if not identifier:
        return None
    sequence = sequence_data.get(identifier)
    if sequence is not None:
        return sequence
    if identifier.startswith("ENS") and "." in identifier:
        bare, version = _split_ens_version(identifier)
        if version is not None:
            sequence = sequence_data.get(bare)
            if sequence is not None:
                return sequence
    stripped_index = getattr(sequence_data, "_stripped_index", None)
    if stripped_index:
        versioned = stripped_index.get(identifier)
        if versioned is not None:
            return sequence_data.get(versioned)
    return None


class SequenceData(object):
    """
    Container for reference nucleotide and amino acid sequenes.
    """

    def __init__(self, fasta_paths, cache_directory_path=None):
        if type(fasta_paths) is str:
            fasta_paths = [fasta_paths]

        self.fasta_paths = [abspath(path) for path in fasta_paths]
        self.fasta_directory_paths = [split(path)[0] for path in self.fasta_paths]
        self.fasta_filenames = [split(path)[1] for path in self.fasta_paths]
        if cache_directory_path:
            self.cache_directory_paths = [cache_directory_path] * len(self.fasta_paths)
        else:
            self.cache_directory_paths = self.fasta_directory_paths
        for path in self.fasta_paths:
            if not exists(path):
                raise ValueError("Couldn't find FASTA file %s" % (path,))
        self.fasta_dictionary_filenames = [
            filename + ".pickle" for filename in self.fasta_filenames
        ]
        self.fasta_dictionary_pickle_paths = [
            join(cache_path, filename)
            for cache_path, filename in zip(
                self.cache_directory_paths, self.fasta_dictionary_filenames
            )
        ]
        self._init_lazy_fields()

    def _init_lazy_fields(self):
        self._fasta_dictionary = None
        self._fasta_keys = None
        # Maps a bare ENS ID -> the versioned form actually keyed in
        # _fasta_dictionary, so callers passing the unversioned form can
        # still resolve to a versioned FASTA entry (GENCODE case).
        self._stripped_index = None
        # Maps a sequence identifier -> the integer version suffix the
        # FASTA header carried, or absent for headers without a version.
        # Lets callers sanity-check FASTA / GTF alignment.
        self._versions = None

    def clear_cache(self):
        self._init_lazy_fields()
        for path in self.fasta_dictionary_pickle_paths:
            if exists(path):
                remove(path)

    def __str__(self):
        return "SequenceData(fasta_paths=%s)" % (self.fasta_paths,)

    def __repr__(self):
        return str(self)

    def __contains__(self, sequence_id):
        if self._fasta_keys is None:
            self._fasta_keys = set(self.fasta_dictionary.keys())
        return sequence_id in self._fasta_keys

    def __eq__(self, other):
        # test to see if self.fasta_paths and other.fasta_paths contain
        # the same list of paths, regardless of order
        return (other.__class__ is SequenceData) and Counter(
            self.fasta_paths
        ) == Counter(other.fasta_paths)

    def __hash__(self):
        return hash(self.fasta_paths)

    def _add_to_fasta_dictionary(self, fasta_dictionary_tmp):
        for identifier, sequence in fasta_dictionary_tmp.items():
            if identifier in self._fasta_dictionary:
                logger.warn(
                    "Sequence identifier %s is duplicated in your FASTA files!"
                    % identifier
                )
                continue
            self._fasta_dictionary[identifier] = sequence
            bare, version = _split_ens_version(identifier)
            if version is not None:
                if bare in self._stripped_index:
                    # Two FASTA entries collide on the same bare ENS ID with
                    # different versions. The newer of the two is kept as
                    # the canonical alias since version numbers grow
                    # monotonically.
                    existing = self._stripped_index[bare]
                    _, existing_version = _split_ens_version(existing)
                    if existing_version is None or version > existing_version:
                        self._stripped_index[bare] = identifier
                else:
                    self._stripped_index[bare] = identifier
                self._versions[identifier] = version

    def _load_or_create_fasta_dictionary_pickle(self):
        self._fasta_dictionary = dict()
        self._stripped_index = dict()
        self._versions = dict()
        for fasta_path, pickle_path in zip(
            self.fasta_paths, self.fasta_dictionary_pickle_paths
        ):
            if exists(pickle_path):
                # try loading the cached file
                # but we'll fall back on recreating it if loading fails
                try:
                    fasta_dictionary_tmp = load_pickle(pickle_path)
                    self._add_to_fasta_dictionary(fasta_dictionary_tmp)
                    logger.info("Loaded sequence dictionary from %s", pickle_path)
                    continue
                except (pickle.UnpicklingError, AttributeError):
                    # catch either an UnpicklingError or an AttributeError
                    # resulting from pickled objects refering to classes
                    # that no longer exists
                    logger.warn(
                        "Failed to load %s, attempting to read FASTA directly",
                        pickle_path,
                    )
            logger.info("Parsing sequences from FASTA file at %s", fasta_path)

            fasta_dictionary_tmp = parse_fasta_dictionary(fasta_path)
            self._add_to_fasta_dictionary(fasta_dictionary_tmp)
            logger.info("Saving sequence dictionary to %s", pickle_path)
            datacache.ensure_dir(dirname(pickle_path))
            dump_pickle(fasta_dictionary_tmp, pickle_path)

    def index(self, overwrite=False):
        if overwrite:
            self.clear_cache()
        self._load_or_create_fasta_dictionary_pickle()

    @property
    def fasta_dictionary(self):
        if not self._fasta_dictionary:
            self._load_or_create_fasta_dictionary_pickle()
        return self._fasta_dictionary

    def get(self, sequence_id):
        """Get sequence associated with given ID or return None if missing"""
        return self.fasta_dictionary.get(sequence_id)

    def fasta_version(self, sequence_id):
        """
        Return the integer FASTA-header version for ``sequence_id``, or
        ``None`` if the header didn't carry a version (e.g. older Ensembl
        releases) or the ID isn't in the FASTA.

        Accepts either the versioned or bare form of an ENS ID — the
        bare form is resolved through ``_stripped_index``.
        """
        if not sequence_id:
            return None
        # ensure lazy load
        _ = self.fasta_dictionary
        version = self._versions.get(sequence_id)
        if version is not None:
            return version
        versioned = self._stripped_index.get(sequence_id)
        if versioned is not None:
            return self._versions.get(versioned)
        return None

1	# Licensed under the Apache License, Version 2.0 (the "License");
2	# you may not use this file except in compliance with the License.
3	# You may obtain a copy of the License at
4	#
5	# http://www.apache.org/licenses/LICENSE-2.0
6	#
7	# Unless required by applicable law or agreed to in writing, software
8	# distributed under the License is distributed on an "AS IS" BASIS,
9	# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
10	# See the License for the specific language governing permissions and
11	# limitations under the License.
12
13	from os import remove	4✔
14	from os.path import dirname, exists, abspath, split, join	4✔
15
16	import datacache	4✔
17	import logging	4✔
18	from collections import Counter	4✔
19	import pickle	4✔
20	from .common import load_pickle, dump_pickle	4✔
21	from .fasta import _split_ens_version, parse_fasta_dictionary	4✔
22
23
24	logger = logging.getLogger(__name__)	4✔
25
26
27	def lookup_sequence_with_version_fallback(sequence_data, identifier):	4✔
28	"""
29	Look up ``identifier`` in ``sequence_data``, tolerating ENS ``.N``
30	version-suffix mismatches in either direction.
31
32	Both Ensembl and GENCODE work with versioned IDs (e.g.
33	``ENSP00000123456.3``); the two formats just split that information
34	differently. Ensembl GTFs keep the bare ID in ``protein_id`` /
35	``transcript_id`` and put the version in a separate ``*_version``
36	attribute, while GENCODE GTFs embed the version directly in the ID
37	attribute. Whichever convention the FASTA on disk used, the
38	sequence dict now keys on the form that header carried — so either
39	a versioned or a bare caller-supplied ID may need a fallback.
40
41	Resolution order:
42	1. literal lookup of ``identifier``
43	2. if ``identifier`` is a versioned ENS ID, strip ``.N`` and retry
44	(handles caller-supplied versioned ID against a bare FASTA)
45	3. consult ``sequence_data._stripped_index`` for a versioned alias
46	of a bare caller-supplied ID (handles bare ID against a
47	versioned FASTA, the GENCODE case)
48	4. None
49
50	Non-Ensembl IDs (e.g. TAIR ``AT1G01010.1`` where ``.N`` is an
51	isoform suffix, not a version) skip step (2): stripping there is
52	semantically wrong.
53	"""
54	if not identifier:	4✔
55	return None	4✔
56	sequence = sequence_data.get(identifier)	4✔
57	if sequence is not None:	4✔
58	return sequence	4✔
59	if identifier.startswith("ENS") and "." in identifier:	4✔
60	bare, version = _split_ens_version(identifier)	4✔
61	if version is not None:	4✔
62	sequence = sequence_data.get(bare)	4✔
63	if sequence is not None:	4✔
64	return sequence	4✔
65	stripped_index = getattr(sequence_data, "_stripped_index", None)	4✔
66	if stripped_index:	4✔
67	versioned = stripped_index.get(identifier)	4✔
68	if versioned is not None:	4✔
69	return sequence_data.get(versioned)	4✔
70	return None	4✔
71
72
73	class SequenceData(object):	4✔
74	"""
75	Container for reference nucleotide and amino acid sequenes.
76	"""
77
78	def __init__(self, fasta_paths, cache_directory_path=None):	4✔
79	if type(fasta_paths) is str:	4✔
80	fasta_paths = [fasta_paths]	×
81
82	self.fasta_paths = [abspath(path) for path in fasta_paths]	4✔
83	self.fasta_directory_paths = [split(path)[0] for path in self.fasta_paths]	4✔
84	self.fasta_filenames = [split(path)[1] for path in self.fasta_paths]	4✔
85	if cache_directory_path:	4✔
86	self.cache_directory_paths = [cache_directory_path] * len(self.fasta_paths)	4✔
87	else:
88	self.cache_directory_paths = self.fasta_directory_paths	×
89	for path in self.fasta_paths:	4✔
90	if not exists(path):	4✔
91	raise ValueError("Couldn't find FASTA file %s" % (path,))	×
92	self.fasta_dictionary_filenames = [	4✔
93	filename + ".pickle" for filename in self.fasta_filenames
94	]
95	self.fasta_dictionary_pickle_paths = [	4✔
96	join(cache_path, filename)
97	for cache_path, filename in zip(
98	self.cache_directory_paths, self.fasta_dictionary_filenames
99	)
100	]
101	self._init_lazy_fields()	4✔
102
103	def _init_lazy_fields(self):	4✔
104	self._fasta_dictionary = None	4✔
105	self._fasta_keys = None	4✔
106	# Maps a bare ENS ID -> the versioned form actually keyed in
107	# _fasta_dictionary, so callers passing the unversioned form can
108	# still resolve to a versioned FASTA entry (GENCODE case).
109	self._stripped_index = None	4✔
110	# Maps a sequence identifier -> the integer version suffix the
111	# FASTA header carried, or absent for headers without a version.
112	# Lets callers sanity-check FASTA / GTF alignment.
113	self._versions = None	4✔
114
115	def clear_cache(self):	4✔
116	self._init_lazy_fields()	4✔
117	for path in self.fasta_dictionary_pickle_paths:	4✔
118	if exists(path):	4✔
119	remove(path)	4✔
120
121	def __str__(self):	4✔
122	return "SequenceData(fasta_paths=%s)" % (self.fasta_paths,)	×
123
124	def __repr__(self):	4✔
125	return str(self)	×
126
127	def __contains__(self, sequence_id):	4✔
128	if self._fasta_keys is None:	×
129	self._fasta_keys = set(self.fasta_dictionary.keys())	×
130	return sequence_id in self._fasta_keys	×
131
132	def __eq__(self, other):	4✔
133	# test to see if self.fasta_paths and other.fasta_paths contain
134	# the same list of paths, regardless of order
135	return (other.__class__ is SequenceData) and Counter(	×
136	self.fasta_paths
137	) == Counter(other.fasta_paths)
138
139	def __hash__(self):	4✔
140	return hash(self.fasta_paths)	×
141
142	def _add_to_fasta_dictionary(self, fasta_dictionary_tmp):	4✔
143	for identifier, sequence in fasta_dictionary_tmp.items():	4✔
144	if identifier in self._fasta_dictionary:	4✔
145	logger.warn(	×
146	"Sequence identifier %s is duplicated in your FASTA files!"
147	% identifier
148	)
149	continue	×
150	self._fasta_dictionary[identifier] = sequence	4✔
151	bare, version = _split_ens_version(identifier)	4✔
152	if version is not None:	4✔
153	if bare in self._stripped_index:	4✔
154	# Two FASTA entries collide on the same bare ENS ID with
155	# different versions. The newer of the two is kept as
156	# the canonical alias since version numbers grow
157	# monotonically.
158	existing = self._stripped_index[bare]	4✔
159	_, existing_version = _split_ens_version(existing)	4✔
160	if existing_version is None or version > existing_version:	4✔
161	self._stripped_index[bare] = identifier	4✔
162	else:
163	self._stripped_index[bare] = identifier	4✔
164	self._versions[identifier] = version	4✔
165
166	def _load_or_create_fasta_dictionary_pickle(self):	4✔
167	self._fasta_dictionary = dict()	4✔
168	self._stripped_index = dict()	4✔
169	self._versions = dict()	4✔
170	for fasta_path, pickle_path in zip(	4✔
171	self.fasta_paths, self.fasta_dictionary_pickle_paths
172	):
173	if exists(pickle_path):	4✔
174	# try loading the cached file
175	# but we'll fall back on recreating it if loading fails
176	try:	4✔
177	fasta_dictionary_tmp = load_pickle(pickle_path)	4✔
178	self._add_to_fasta_dictionary(fasta_dictionary_tmp)	4✔
179	logger.info("Loaded sequence dictionary from %s", pickle_path)	4✔
180	continue	4✔
181	except (pickle.UnpicklingError, AttributeError):	×
182	# catch either an UnpicklingError or an AttributeError
183	# resulting from pickled objects refering to classes
184	# that no longer exists
185	logger.warn(	×
186	"Failed to load %s, attempting to read FASTA directly",
187	pickle_path,
188	)
189	logger.info("Parsing sequences from FASTA file at %s", fasta_path)	4✔
190
191	fasta_dictionary_tmp = parse_fasta_dictionary(fasta_path)	4✔
192	self._add_to_fasta_dictionary(fasta_dictionary_tmp)	4✔
193	logger.info("Saving sequence dictionary to %s", pickle_path)	4✔
194	datacache.ensure_dir(dirname(pickle_path))	4✔
195	dump_pickle(fasta_dictionary_tmp, pickle_path)	4✔
196
197	def index(self, overwrite=False):	4✔
198	if overwrite:	4✔
199	self.clear_cache()	×
200	self._load_or_create_fasta_dictionary_pickle()	4✔
201
202	@property	4✔
203	def fasta_dictionary(self):	4✔
204	if not self._fasta_dictionary:	4✔
205	self._load_or_create_fasta_dictionary_pickle()	4✔
206	return self._fasta_dictionary	4✔
207
208	def get(self, sequence_id):	4✔
209	"""Get sequence associated with given ID or return None if missing"""
210	return self.fasta_dictionary.get(sequence_id)	4✔
211
212	def fasta_version(self, sequence_id):	4✔
213	"""
214	Return the integer FASTA-header version for ``sequence_id``, or
215	``None`` if the header didn't carry a version (e.g. older Ensembl
216	releases) or the ID isn't in the FASTA.
217
218	Accepts either the versioned or bare form of an ENS ID — the
219	bare form is resolved through ``_stripped_index``.
220	"""
221	if not sequence_id:	4✔
NEW 222	return None	×
223	# ensure lazy load
224	_ = self.fasta_dictionary	4✔
225	version = self._versions.get(sequence_id)	4✔
226	if version is not None:	4✔
227	return version	4✔
228	versioned = self._stripped_index.get(sequence_id)	4✔
229	if versioned is not None:	4✔
230	return self._versions.get(versioned)	4✔
231	return None	4✔

openvax / pyensembl / 25816358406

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous