From a6c3f5ffefb3e7e002208b5ce3e5158134340e3a Mon Sep 17 00:00:00 2001 From: Yin Hoon Chew Date: Wed, 4 Mar 2020 16:56:09 -0500 Subject: [PATCH] Account for cases where forward rate law may not exist and fix impute kinetic constant for cases where the min and max fluxes are equivalent --- tests/eukaryote/test_metabolism.py | 26 ++++++++++++++- wc_model_gen/eukaryote/metabolism.py | 49 +++++++++++++++------------- 2 files changed, 51 insertions(+), 24 deletions(-) diff --git a/tests/eukaryote/test_metabolism.py b/tests/eukaryote/test_metabolism.py index 92e02ea..4ba81ef 100644 --- a/tests/eukaryote/test_metabolism.py +++ b/tests/eukaryote/test_metabolism.py @@ -752,11 +752,16 @@ def test_impute_kinetic_constant(self): model = self.model gen = self.gen - gen.clean_and_validate_options() + gen.clean_and_validate_options() bound_values = {'ex_m1': (10, 10), 'ex_m2': (10, 10), 'ex_m3': (0, 0), 'r1': (0, 1), 'r2': (-2, 3), 'r3': (0, 3), 'r4': (0, None), 'biomass_reaction': (0, None)} gen.impute_kinetic_constant(bound_values) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r1-forward')], 1) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r1-backward')], 0) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r2-forward')], 3) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r2-backward')], 2) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r3-forward')], 3) self.assertEqual(model.parameters.get_one(id='k_cat_r1_forward_enzyme3').value, 150.) self.assertEqual(model.parameters.get_one(id='k_cat_r1_forward_enzyme3').comments, 'Value imputed as the median of measured k_cat values') self.assertEqual(model.parameters.get_one(id='k_cat_r1_backward_enzyme3').value, 150.) @@ -771,3 +776,22 @@ def test_impute_kinetic_constant(self): self.assertEqual(model.parameters.get_one(id='k_cat_r2_backward_enzyme2').comments, 'Value imputed based on FVA bound value') self.assertEqual(model.parameters.get_one(id='k_cat_r3_forward_enzyme2').value, 300.) self.assertEqual(model.parameters.get_one(id='k_cat_r3_forward_enzyme2').comments, 'Measured value adjusted to relax bound') + + bound_values = {'ex_m1': (10, 10), 'ex_m2': (10, 10), 'ex_m3': (0, 0), 'r1': (1, 1), 'r2': (-2, -2), 'r3': (1, 3), 'r4': (0, None), 'biomass_reaction': (0, None)} + gen.impute_kinetic_constant(bound_values) + + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r1-forward')], 1) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r1-backward')], 0) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r2-forward')], 0) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r2-backward')], 2) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r3-forward')], 3) + + bound_values = {'ex_m1': (10, 10), 'ex_m2': (10, 10), 'ex_m3': (0, 0), 'r1': (0, 0), 'r2': (-4, -2), 'r3': (3, 3), 'r4': (0, None), 'biomass_reaction': (0, None)} + gen.impute_kinetic_constant(bound_values) + + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r1-forward')], 0) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r1-backward')], 0) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r2-forward')], 0) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r2-backward')], 4) + self.assertEqual(gen._law_bound_pairs[model.rate_laws.get_one(id='r3-forward')], 3) + \ No newline at end of file diff --git a/wc_model_gen/eukaryote/metabolism.py b/wc_model_gen/eukaryote/metabolism.py index d67b186..ad89964 100644 --- a/wc_model_gen/eukaryote/metabolism.py +++ b/wc_model_gen/eukaryote/metabolism.py @@ -450,16 +450,19 @@ def calibrate_submodel(self): elif reaction.rate_laws: for_ratelaw = reaction.rate_laws.get_one(direction=wc_lang.RateLawDirection.forward) - if all(i.distribution_init_concentration.mean==0. for i in for_ratelaw.expression.species): - max_constr = 0. - elif not any(numpy.isnan(p.value) for p in for_ratelaw.expression.parameters): - max_constr = for_ratelaw.expression._parsed_expression.eval({ - wc_lang.Species: {i.id: i.distribution_init_concentration.mean \ - for i in for_ratelaw.expression.species} - }) * scale_factor - upper_bound_adjustable.append(reaction.id) + if for_ratelaw: + if all(i.distribution_init_concentration.mean==0. for i in for_ratelaw.expression.species): + max_constr = 0. + elif not any(numpy.isnan(p.value) for p in for_ratelaw.expression.parameters): + max_constr = for_ratelaw.expression._parsed_expression.eval({ + wc_lang.Species: {i.id: i.distribution_init_concentration.mean \ + for i in for_ratelaw.expression.species} + }) * scale_factor + upper_bound_adjustable.append(reaction.id) + else: + max_constr = None else: - max_constr = None + max_constr = None rev_ratelaw = reaction.rate_laws.get_one(direction=wc_lang.RateLawDirection.backward) if rev_ratelaw: @@ -785,24 +788,24 @@ def impute_kinetic_constant(self, bound_values): median_kcat = numpy.median([k.value for i in submodel.reactions for j in i.rate_laws \ for k in j.expression.parameters if k.value and not numpy.isnan(k.value)]) - law_bound_pairs = {} - for reaction in submodel.reactions: + self._law_bound_pairs = {} + for reaction in submodel.reactions: if reaction.rate_laws: bounds = bound_values[reaction.id] - if len(reaction.rate_laws)==1: - law_bound_pairs[reaction.rate_laws[0]] = bounds[1] - else: - if bounds[0]<=0 and bounds[1]>=0: - law_bound_pairs[reaction.rate_laws[0]] = bounds[1] - law_bound_pairs[reaction.rate_laws[1]] = -bounds[0] - elif bounds[0]>=0 and bounds[1]>=0: - law_bound_pairs[reaction.rate_laws[0]] = bounds[1] - law_bound_pairs[reaction.rate_laws[1]] = 0. + for_ratelaw = reaction.rate_laws.get_one(direction=wc_lang.RateLawDirection.forward) + rev_ratelaw = reaction.rate_laws.get_one(direction=wc_lang.RateLawDirection.backward) + if for_ratelaw: + if bounds[1] >= 0: + self._law_bound_pairs[for_ratelaw] = bounds[1] + else: + self._law_bound_pairs[for_ratelaw] = 0. + if rev_ratelaw: + if bounds[0] <= 0: + self._law_bound_pairs[rev_ratelaw] = -bounds[0] else: - law_bound_pairs[reaction.rate_laws[0]] = 0. - law_bound_pairs[reaction.rate_laws[1]] = -bounds[0] + self._law_bound_pairs[rev_ratelaw] = 0. - for law, value in law_bound_pairs.items(): + for law, value in self._law_bound_pairs.items(): complexes = law.expression.species if value!=0.0: if len(complexes) == 1: