EnzymeML / PyEnzyme

🧬 - Data management and modeling framework based on EnzymeML.
BSD 2-Clause "Simplified" License
21 stars 9 forks source link

TL_Pysces.write does not write out global parameters #24

Closed fbergmann closed 2 years ago

fbergmann commented 2 years ago

I think the current implementation of the pysces thinlayer does not update global parameters.

JR-1991 commented 2 years ago

This is not necessary since global parameters used within a reaction's rate law are references to global parameters found in the global_parameters attribute of the EnzymeMLDocument instance. Thus, when the write-method assigns estimated parameters from lmfit's Parameters instance, these will be updated simultaneously.

Ran a quick test to check, whether that's the case using the example from our manuscript with no assigned values yet:

enzmldoc = EnzymeMLDocument.fromFile("Model_4.omex")
for reaction in enzmldoc.reaction_dict.values():
    if reaction.model:
        for param in reaction.model.parameters:
            param.value = 10000

which resulted in the expected output

<listOfParameters>
      <parameter id="K_si" value="10000" constant="false"/>
      <parameter id="K_n" value="10000" constant="false"/>
      <parameter id="v_r" value="10000" constant="true"/>
</listOfParameters>

This was done to ensure global parameters behave global in the PyEnzyme environment as it is expected in the SBML document.

I find this is a quite handy way and eliminates possible boilerplate, but if a global parameter is not present in any equation this might be problematic. Is this the case or can we be sure that any global parameter is part of at least a pair of equations?

fbergmann commented 2 years ago

for the copasi class, i added:

           if reaction is not None:
                enz_reaction = nu_enzmldoc.getReaction(reaction.getSBMLId())
                if enz_reaction:
                    parameter = enz_reaction.model.getParameter(name)
                    parameter.value = value
            else:
                p = nu_enzmldoc.global_parameters.get(name)
                p.value = value

to ensure that global parameters would be written out. That should always work, right?

JR-1991 commented 2 years ago

That should work in most cases, but if a reaction is not found in the document a SpeciesNotFoundError will be raised. Hence, I would either recommend a try/except block or using the get-method of built-in dict which will result in a None if the reaction is not present:

           if reaction is not None:
                enz_reaction = nu_enzmldoc.reaction_dict.get(reaction_id)
                if enz_reaction:
                    parameter = enz_reaction.model.getParameter(name)
                    parameter.value = value
            else:
                p = nu_enzmldoc.global_parameters.get(name)
                p.value = value
fbergmann commented 2 years ago

that is checked before, since copasi reactions will only exist for reactions that are already in the enzyme ml document.