v. 0.0.3.16: fixed bug with changing kcat values for catalytic reacti…

…ons with underscores in reaction id

src/PAModelpy/CatalyticEvent.py

@@ -11,6 +11,7 @@
 from typing import Optional, Dict, Union
 from warnings import warn
 from copy import copy, deepcopy
+import re
 
 
 class CatalyticEvent(Object):
@@ -393,6 +394,27 @@ def change_kcat_values(self, enzyme_kcat_dict : dict):
                 enzyme_var.kcats[ce_reaction.id][direction] = kcat
                 self._model._change_kcat_in_enzyme_constraint(ce_reaction, enzyme, direction, kcat)
 
+    @staticmethod
+    def _extract_reaction_id_from_catalytic_reaction_id(input_str: str) -> str:
+        # Define the regex pattern for protein IDs, obtained from UniProtKB, 2024-08-07
+        # https://www.uniprot.org/help/accession_numbers
+        protein_id_pattern = r'(?:[OPQ][0-9][A-Z0-9]{3}[0-9]|[A-NR-Z][0-9]([A-Z][A-Z0-9]{2}[0-9]){1,2})'
+
+        # Remove the 'CE_' prefix if it exists
+        if input_str.startswith('CE_'):
+            input_str = input_str[3:]
+
+        # Define the regex pattern to match protein IDs
+        protein_id_regex = re.compile(protein_id_pattern)
+
+        # split off all protein ids from the reaction
+        reaction_id = protein_id_regex.split(input_str)[0]
+
+        # Remove any trailing or leading underscores that might remain
+        reaction_id = reaction_id.strip('_')
+
+        return reaction_id
+
     def __copy__(self) -> "CatalyticEvent":
         """
         Copy the CatalyticEvent.

src/PAModelpy/Enzyme.py

@@ -245,7 +245,12 @@ def change_kcat_values(self, rxn2kcat: Dict):
 
         # update the enzyme variables
         for rxn_id, kcats in rxn2kcat.items():
-            catalytic_event_id = self.catalytic_event_id.format(rxn_id)
+
+            if 'CE_' not in rxn_id:
+                catalytic_event_id = self.catalytic_event_id.format(rxn_id)
+            else:
+                catalytic_event_id = f"CE_{CatalyticEvent._extract_reaction_id_from_catalytic_reaction_id(rxn_id)}"
+
             # change rxn2kcat dictionary
             self.rxn2kcat[rxn_id] = kcats
             # is there already a link between enzyme and reaction?

src/PAModelpy/PAModel.py

@@ -15,6 +15,7 @@
 import pandas as pd
 from copy import copy, deepcopy
 import inspect
+import re
 
 # sys.path.append('../')
 from .EnzymeSectors import (
@@ -1374,13 +1375,14 @@ def change_kcat_value(self, enzyme_id:str, kcats:dict):
             # also change the active enzyme sector
             active_enzyme = self.sectors.get_by_id("ActiveEnzymeSector")
             for rxn, kcat_f_b in kcats.items():
-                # if a catalytic reaction is given, then extract the actual reaction id from it
+                # if a catalytic reaction is given, then extract the actual reaction id from it using the protein id convention from uniprot
                 if 'CE' in rxn:
-                    rxn = rxn.split('_')[1]
+                    rxn = CatalyticEvent._extract_reaction_id_from_catalytic_reaction_id(rxn)
                 active_enzyme.rxn2protein[rxn][enzyme_id] = kcat_f_b
         else:
             warnings.warn(f'The enzyme {enzyme_id} does not exist in the model. The kcat can thus not be changed.')
 
+
     def _change_kcat_in_enzyme_constraint(self, rxn:Union[str, cobra.Reaction], enzyme_id: str,
                                                 direction: str, kcat: float):
         constraint_id = f'EC_{enzyme_id}_{direction}'

src/PAModelpy/__init__.py

@@ -1,4 +1,4 @@
-print('Loading PAModelpy modules version 0.0.3.15')
+print('Loading PAModelpy modules version 0.0.3.16')
 
 from .Enzyme import *
 from .EnzymeSectors import *

src/pyproject.toml

@@ -9,7 +9,7 @@ find = {}  # Scan the project directory with the default parameters
 
 [project]
 name = 'PAModelpy'
-version = '0.0.3.15'
+version = '0.0.3.16'
 authors = [{name='Samira van den Bogaard', email = '[email protected]'}]
 description = 'Python framework for building and analysing protein allocation models'
 readme = 'README.md'

tests/unit_tests/test_pamodel/test_pamodel.py

@@ -1,12 +1,13 @@
 import pytest
 from cobra.io import load_json_model
 
-from src.PAModelpy import PAModel,Config,ActiveEnzymeSector, UnusedEnzymeSector, TransEnzymeSector
+from src.PAModelpy import PAModel,Config,ActiveEnzymeSector, UnusedEnzymeSector, TransEnzymeSector, CatalyticEvent
+from Scripts.pam_generation_uniprot_id import set_up_ecoli_pam
 from tests.unit_tests.test_pamodel.test_pam_generation import set_up_toy_pam_with_isozymes_and_enzymecomplex
 
 def test_if_pamodel_change_kcat_function_works():
     #arrange
-    sut = build_toy_pam()
+    sut = build_toy_pam(sensitivity=False)
     input_kcat = 10
     enzyme_id = 'E1'
     rxn= sut.reactions.get_by_id('CE_R1_'+enzyme_id)
@@ -26,6 +27,42 @@ def test_if_pamodel_change_kcat_function_works():
     assert input_kcat == pytest.approx(model_kcat_b, 1e-4)
     assert input_kcat == pytest.approx(model_kcat_f, 1e-4)
 
+def test_if_pamodel_gets_reaction_id_from_complex_catalytic_reaction_id():
+    # Arrange
+    sut = build_toy_pam(sensitivity=False)
+    ce_reaction_id = "CE_CYTBO3_4pp_P0ABJ1_P0ABJ5_P0ABJ7_P0ABJ8"
+    reaction_id = "CYTBO3_4pp"
+
+    # Act
+    parsed_reaction_id = CatalyticEvent._extract_reaction_id_from_catalytic_reaction_id(ce_reaction_id)
+
+    # Assert
+    assert parsed_reaction_id == reaction_id
+
+def test_if_pamodel_change_kcat_function_works_with_catalytic_reactions():
+    #arrange
+    sut = set_up_ecoli_pam(sensitivity=False)
+    input_kcat = 10
+    enzyme_id = 'P0ABJ1'
+    rxn_id = "CYTBO3_4pp"
+    ce_rxn= sut.reactions.query(f'CE_{rxn_id}_{enzyme_id}')[0]
+    enzyme_complex_id = "_".join(ce_rxn.id.split("_")[3:])
+    constraint_name = f'EC_{enzyme_complex_id}_'
+
+    #act
+    sut.change_kcat_value(enzyme_complex_id, kcats ={ce_rxn.id:{'f': input_kcat, 'b': input_kcat}})
+    coeff_b = sut.constraints[constraint_name+'b'].get_linear_coefficients([ce_rxn.reverse_variable])[ce_rxn.reverse_variable]
+    #/(3600*1e-6) to correct for dimension modifications in the model
+    model_kcat_b = 1/coeff_b/(3600*1e-6)
+
+    coeff_f = sut.constraints[constraint_name+'f'].get_linear_coefficients([ce_rxn.forward_variable])[ce_rxn.forward_variable]
+    # /(3600*1e-6) to correct for dimension modifications in the model
+    model_kcat_f = 1/coeff_f/(3600*1e-6)
+
+    #assert
+    assert input_kcat == pytest.approx(model_kcat_b, 1e-4)
+    assert input_kcat == pytest.approx(model_kcat_f, 1e-4)
+
 def test_if_pamodel_change_reaction_bounds_function_works_with_sensitivity_constraints():
     #arrange
     toy_pam = build_toy_pam(sensitivity=True)