RESID in ProForma

AUTHORS

@@ -18,6 +18,7 @@ Ralf Gabriels [RalfG]
 Dikshant jha [dikshant182004]
 Copilot
 Yuhui Hong [josiehong]
+Anthony Cesnik [acesnik]
 
 
 On Bitbucket (before March 1, 2020)

CHANGELOG

@@ -1,3 +1,15 @@
+dev
+---
+
+  - Support RESID in ProForma, adding :py:class:`~.ResidResolver` and
+    :py:class:`~.ResidModification` (`#214 <https://github.com/levitsky/pyteomics/pull/214>`_
+    by Anthony Cesnik). Both the ``RESID`` and ``R`` prefixes resolve, completing the
+    ProForma 2.0 controlled-vocabulary list, and resolving unprefixed RESID names
+    and accessions such as ``[AA0038]`` through the generic resolver chain.
+
+    .. note ::
+      ProForma parsing now requires psims 1.4.0 or newer.
+
 5.0.1
 -----
 

doc/source/api/proforma.rst

@@ -83,7 +83,7 @@ coupled with minimal information about mass and position data.
 Dependencies
 ------------
 
-To resolve PSI-MOD, XL-MOD, and GNO identifiers, :mod:`psims` is required. By default,
+To resolve PSI-MOD, XL-MOD, GNO and RESID identifiers, :mod:`psims` is required. By default,
 :mod:`psims` retrieves the most recent version of each controlled vocabulary from the internet, but
 includes a fall-back version to use when the network is unavailable. It can also create
 an application cache on disk.
@@ -134,7 +134,7 @@ These features are independent from each other:
 
 3. Top Down Extensions
 
-    - [ ] Additional CV/ontologies for protein modifications: RESID (the prefix R MUST be used for RESID CV/ontology term names)
+    - [x] Additional CV/ontologies for protein modifications: RESID (the prefix R MUST be used for RESID CV/ontology term names)
     - [x] Chemical formulas (this feature occurs in two places in this list).
 
 4. Cross-Linking Extensions
@@ -198,6 +198,8 @@ Modification Tags
 
 .. autoclass:: PSIModModification
 
+.. autoclass:: ResidModification
+
 .. autoclass:: XLMODModification
 
 .. autoclass:: GNOmeModification
@@ -241,4 +243,6 @@ Modification Resolvers
 
 .. autoclass:: XLMODResolver
 
+.. autoclass:: ResidResolver
+
 .. autoclass:: GNOResolver

pyproject.toml

@@ -34,10 +34,12 @@ TDA = ["numpy"]
 graphics = ["matplotlib"]
 DF = ["pandas>=0.17"]
 Unimod = ["lxml", "sqlalchemy >= 1.4"]
-numpress = ["pynumpress"]
+numpress = ["pynumpress < 0.1.2; python_version < '3.11'", "pynumpress; python_version >= '3.11'"]
 mzMLb = ["h5py", "hdf5plugin"]
-proforma = ["psims > 0.1.42"]
-all = ["lxml", "numpy", "psims>0.1.42", "matplotlib", "pandas>=0.17", "sqlalchemy>=1.4", "pynumpress", "h5py", "hdf5plugin", "spectrum_utils>=0.4"]
+proforma = ["psims >= 1.4.0"]
+all = ["lxml", "numpy", "psims >= 1.4.0", "matplotlib", "pandas>=0.17", "sqlalchemy>=1.4",
+       "pynumpress < 0.1.2; python_version < '3.11'", "pynumpress; python_version >= '3.11'",
+       "h5py", "hdf5plugin", "spectrum_utils>=0.4"]
 
 [project.urls]
 Documentation = "http://pyteomics.readthedocs.io"

pyteomics/auxiliary/psims_util.py

@@ -1,13 +1,20 @@
 from functools import partial
 
+
+def _needs_psims(name):
+    raise ImportError("Loading %s requires the `psims` library. To access it, please install `psims`" % name)
+
+
+def _needs_newer_psims(name):
+    raise ImportError("Loading %s requires a newer version of the `psims` library. "
+                      "To access it, please upgrade `psims`" % name)
+
+
 try:
     from psims.controlled_vocabulary.controlled_vocabulary import (load_psimod, load_xlmod, load_gno, obo_cache, load_unimod, load_psims)
     from psims.controlled_vocabulary.relationship import HasValueTypeRelationship
     _has_psims = True
 except ImportError:
-    def _needs_psims(name):
-        raise ImportError("Loading %s requires the `psims` library. To access it, please install `psims`" % name)
-
     load_psimod = partial(_needs_psims, 'PSIMOD')
     load_xlmod = partial(_needs_psims, 'XLMOD')
     load_gno = partial(_needs_psims, 'GNO')
@@ -16,3 +23,12 @@ def _needs_psims(name):
     obo_cache = None
     HasValueTypeRelationship = None
     _has_psims = False
+
+try:
+    # Imported on its own rather than with the others, so that a psims
+    # predating RESID support loses only RESID. Folding it into the tuple above
+    # would turn one missing name into "psims is not installed" and disable
+    # every controlled vocabulary.
+    from psims.controlled_vocabulary.controlled_vocabulary import load_resid
+except ImportError:
+    load_resid = partial(_needs_newer_psims if _has_psims else _needs_psims, 'RESID')

pyteomics/proforma.py

@@ -25,7 +25,8 @@
 from .mass import Composition, std_aa_mass, Unimod, nist_mass, calculate_mass, std_ion_comp, mass_charge_ratio, std_aa_comp
 from .auxiliary import PyteomicsError, BasicComposition
 from .auxiliary.utils import add_metaclass, memoize
-from .auxiliary.psims_util import load_psimod, load_xlmod, load_gno, obo_cache, _has_psims
+from .auxiliary.psims_util import (load_psimod, load_xlmod, load_gno, load_resid, obo_cache,
+                                   _has_psims)
 
 try:
     import numpy as np
@@ -152,6 +153,7 @@ class TagTypeEnum(Enum):
     comkp = 13
     limit = 14
     custom = 15
+    resid = 16
 
     group_placeholder = 999
 
@@ -274,6 +276,7 @@ def is_modification(self) -> bool:
             TagTypeEnum.massmod,
             TagTypeEnum.psimod,
             TagTypeEnum.custom,
+            TagTypeEnum.resid,
         )
 
     def find_modification(self) -> Optional["TagBase"]:
@@ -872,6 +875,58 @@ def _resolve_impl(self, name=None, id=None, **kwargs):
         return rec
 
 
+class ResidResolver(ModificationResolver):
+    """Resolve RESID accessions and names.
+
+    RESID identifiers are not numeric (``AA0038``), so unlike the other
+    prefixed vocabularies they always arrive as a name rather than an id.
+    """
+
+    def __init__(self, **kwargs):
+        super(ResidResolver, self).__init__('resid', **kwargs)
+        self._database = kwargs.get("database")
+
+    def load_database(self):
+        return load_resid()
+
+    def parse_identifier(self, identifier: str):
+        """Strip the CV prefix, tolerating the space the spec writes after it.
+
+        The ProForma specification gives ``EM[R: L-methionine sulfone]...`` as
+        an example, so the remainder needs stripping before lookup.
+        """
+        name, id = super(ResidResolver, self).parse_identifier(identifier)
+        if name is not None:
+            name = name.strip()
+        return name, id
+
+    def _resolve_impl(self, name=None, id=None, **kwargs):
+        if name is not None:
+            defn = self.database[name]
+        elif id is not None:
+            defn = self.database['AA{:04d}'.format(id)]
+        else:
+            raise ValueError("Must provide one of `name` or `id`")
+        mass = defn.monoisotopic_mass
+        if mass is None:
+            # RESID lists no single mass change for this entry, either because
+            # it has none or because the change depends on which residue was
+            # modified. Rather than pick one, point at where the alternatives
+            # are: each carries the parent residue it applies to.
+            raise ModificationMassNotFoundError(
+                "RESID gives no single mass for %r; it lists %d corrections, each for a "
+                "different parent residue. See `.resolver.database[%r].corrections`." % (
+                    (name or id), len(defn.corrections), defn.id))
+        return {
+            'mass': mass,
+            'composition': defn.composition,
+            'name': defn.name,
+            'id': defn.id,
+            'provider': self.name,
+            "source": self
+        }
+
+
 class GenericResolver(ModificationResolver):
 
     def __init__(self, resolvers, **kwargs):
@@ -1422,6 +1477,15 @@ class XLMODModification(ModificationBase):
     _tag_type = TagTypeEnum.xlmod
 
 
+class ResidModification(ModificationBase):
+    __slots__ = ()
+
+    resolver = ResidResolver()
+    prefix_name = "RESID"
+    short_prefix = 'R'
+    _tag_type = TagTypeEnum.resid
+
+
 class CustomModification(ModificationBase):
     __slots__ = ()
 
@@ -1445,6 +1509,10 @@ class GenericModification(ModificationBase):
         PSIModModification.resolver,
         XLMODModification.resolver,
         GNOmeModification.resolver,
+        # After the vocabularies that share names with RESID, so anything they
+        # already answer keeps its existing answer, and before the non-strict
+        # pass so a RESID name is preferred over a fuzzy Unimod match.
+        ResidModification.resolver,
         # Some really common names aren't actually found in the XML exactly, so default
         # to non-strict matching now to avoid masking other sources here.
         partial(UnimodModification.resolver, strict=False)

test-requirements.txt

@@ -7,6 +7,7 @@ psycopg[binary]
 cython
 h5py
 hdf5plugin
-pynumpress
-psims > 0.1.42
+pynumpress < 0.1.2; python_version < '3.11'
+pynumpress; python_version >= '3.11'
+psims >= 1.4.0
 spectrum_utils

tests/test_proforma.py

@@ -10,7 +10,7 @@
 from pyteomics.proforma import (
     Chimeric, PSIModModification, ProForma, TaggedInterval, parse, MassModification, ProFormaError, TagTypeEnum,
     ModificationRule, StableIsotope, GenericModification, Composition, to_proforma, ModificationMassNotFoundError,
-    UnimodModification, ModificationTarget,
+    UnimodModification, ModificationTarget, ResidModification,
     AdductParser, ChargeState, proteoforms, _coerce_string_to_modification,
     std_aa_comp, obo_cache, process_tag_tokens, peptidoforms)
 from pyteomics import mass
@@ -450,7 +450,7 @@ def test_from_spec(self):
             "<[TMT6plex]@K,N-term:A,N-term:B>ATPEILTCNSIGCLK",
             "EM[Oxidation]EVEES[Phospho]PEK",
             "EM[L-methionine sulfoxide]EVEES[O-phospho-L-serine]PEK",
-            # "EM[R: L-methionine sulfone]EVEES[O-phospho-L-serine]PEK", # don't support RESID
+            "EM[R: L-methionine sulfone]EVEES[O-phospho-L-serine]PEK",
             "EMEVTK[X:DSS#XL1]SESPEK",
             "NEEYN[GNO:G59626AS]K",
             "NEEYN[G:G59626AS]K",
@@ -461,7 +461,7 @@ def test_from_spec(self):
             "EM[Oxidation]EVE[Cation:Mg[II]]ES[Phospho]PEK",
             "EM[MOD:00719]EVEES[MOD:00046]PEK",
             "EM[UNIMOD:35]EVEES[UNIMOD:56]PEK",
-            # "EM[RESID:AA0581]EVEES[RESID:AA0037]PEK", # don't support RESID
+            "EM[RESID:AA0581]EVEES[RESID:AA0037]PEK",
             "EMEVTK[XLMOD:02001#XL1]SESPEK[#XL1]",
             "EMK[XLMOD:02000#XL1]EVTKSE[XLMOD:02010#XL2]SK[#XL1]PEK[#XL2]AR",
             "EMEVTK[XLMOD:02001#XL1]SESPEK",
@@ -800,6 +800,68 @@ def test_unknown_mass(self):
         self.assertRaises(ModificationMassNotFoundError, lambda: state.resolve())
 
 
+class ResidModificationResolverTest(unittest.TestCase):
+    def test_resolve_by_accession(self):
+        mod = ResidModification("RESID:AA0038").resolve()
+        self.assertEqual(mod['id'], 'AA0038')
+        self.assertEqual(mod['name'], 'O-phospho-L-threonine')
+        self.assertEqual(mod['provider'], 'resid')
+        # The modification mass, not the mass of the modified residue, which
+        # would be 181.014009.
+        self.assertAlmostEqual(mod['mass'], 79.966331, 6)
+
+    def test_resolve_by_name(self):
+        mod = ResidModification("L-methionine sulfone").resolve()
+        self.assertEqual(mod['id'], 'AA0251')
+        self.assertAlmostEqual(mod['mass'], 31.989829, 6)
+
+    def test_short_prefix_tolerates_the_space_the_spec_writes(self):
+        mod = ResidModification("R: L-methionine sulfone").resolve()
+        self.assertEqual(mod['id'], 'AA0251')
+
+    def test_agrees_with_psimod(self):
+        for resid_id, psimod_id in [("RESID:AA0038", "MOD:00047"),
+                                    ("RESID:AA0055", "MOD:00064")]:
+            self.assertAlmostEqual(ResidModification(resid_id).mass,
+                                   PSIModModification(psimod_id).mass, 4)
+
+    def test_bare_accession_resolves_through_the_generic_chain(self):
+        # `[AA0038]` with no prefix. RESID sits after the vocabularies that
+        # share names with it, so this is additive: nothing that resolved
+        # before changes, and 619 bare accessions that resolved to nothing now
+        # resolve.
+        mod = GenericModification("AA0038").resolve()
+        self.assertEqual(mod["provider"], "resid")
+        self.assertEqual(mod["name"], "O-phospho-L-threonine")
+        self.assertAlmostEqual(GenericModification("AA0038").mass, 79.966331, 6)
+
+    def test_shared_names_keep_their_existing_resolver(self):
+        # RESID and PSI-MOD both call MOD:00046/AA0037 O-phospho-L-serine.
+        # Placing RESID after PSI-MOD means the existing answer stands.
+        self.assertEqual(
+            GenericModification("O-phospho-L-serine").resolve()["provider"], "psimod")
+
+    def test_unknown_mass(self):
+        # 2-pyrrolidone-5-carboxylic acid loses ammonia when formed from
+        # glutamine and water when formed from glutamate, so RESID lists no
+        # single delta and none should be invented.
+        state = ResidModification("RESID:AA0031")
+        self.assertRaises(ModificationMassNotFoundError, lambda: state.resolve())
+
+    def test_ambiguous_entry_still_offers_every_delta(self):
+        # Refusing to pick is only defensible if the alternatives are reachable
+        # and right, so check the two that AA0031 actually has.
+        database = ResidModification("RESID:AA0031").resolver.database
+        by_parent = {database[c.parents[0]].name: c.monoisotopic_mass
+                     for c in database['AA0031'].corrections}
+        self.assertAlmostEqual(by_parent['L-glutamine'], -17.026549, 6)
+        self.assertAlmostEqual(by_parent['L-glutamic acid'], -18.010565, 6)
+
+    def test_parses_and_totals(self):
+        seq = ProForma.parse("EM[RESID:AA0581]EVEES[RESID:AA0037]PEK")
+        self.assertAlmostEqual(seq.mass, 1301.4734302, 4)
+
+
 class ModificationTest(unittest.TestCase):
     def test_mass_modification_hashable(self):
         mod = MassModification(57.08)