Skip to content

Commit 5cbdf83

Browse files
committed
Update the cross_validate method to match the _make_rule method, for fitting kinetics.
I tried using a helper function that could be re-used, but it started to get too messy and complicated. This unfortunately requires continued maintenance keeping three different bits of code in sync doing the same thing, eg. rejecting BM fits when abs(n)>5 or E0<0.
1 parent 3d0b9b1 commit 5cbdf83

1 file changed

Lines changed: 10 additions & 8 deletions

File tree

rmgpy/data/kinetics/family.py

Lines changed: 10 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -3678,19 +3678,22 @@ def cross_validate(self, folds=5, template_rxn_map=None, test_rxn_inds=None, T=1
36783678

36793679
uncertainties[rxn] = self.rules.entries[entry.label][0].data.uncertainty
36803680

3681-
36823681
L = list(set(template_rxn_map[entry.label]) - set(rxns_test))
36833682

36843683
if L != []:
36853684
if isinstance(L[0].kinetics, Arrhenius):
36863685
kinetics = ArrheniusBM().fit_to_reactions(L, recipe=self.forward_recipe.actions)
3687-
if kinetics.E0.value_si < 0.0 or len(L) == 1:
3686+
# These conditions (eg. n>5) for rejecting the BM fit should correspond
3687+
# to those in _marke_rule, and those just below:
3688+
if len(L) == 1 or kinetics.E0.value_si < 0.0 or abs(kinetics.n.value_si) > 5.0:
36883689
kinetics = average_kinetics([r.kinetics for r in L])
36893690
else:
36903691
kinetics = kinetics.to_arrhenius(rxn.get_enthalpy_of_reaction(298.))
36913692
else:
36923693
kinetics = ArrheniusChargeTransferBM().fit_to_reactions(L, recipe=self.forward_recipe.actions)
3693-
if kinetics.E0.value_si < 0.0 or len(L) == 1:
3694+
# These conditions (eg. n>5) for rejecting the BM fit should correspond
3695+
# to those in _marke_rule, and those just above:
3696+
if len(L) == 1 or kinetics.E0.value_si < 0.0 or abs(kinetics.n.value_si) > 5.0:
36943697
kinetics = average_kinetics([r.kinetics for r in L])
36953698
else:
36963699
kinetics = kinetics.to_arrhenius_charge_transfer(rxn.get_enthalpy_of_reaction(298.))
@@ -4612,24 +4615,23 @@ def _make_rule(rr):
46124615
arr = ArrheniusChargeTransferBM
46134616
if n > 1:
46144617
kin = arr().fit_to_reactions(rs, recipe=recipe)
4618+
4619+
# If you update the following conditions, also update the cross_validate function.
46154620
if n == 1 or kin.E0.value_si < 0.0 or abs(kin.n.value_si) > 5.0:
4616-
46174621
if n > 1: #we have an ArrheniusBM
46184622
if kin.E0.value_si < 0.0:
46194623
reason = "E0<0"
46204624
elif abs(kin.n.value_si) > 5.0:
46214625
reason = "abs(n)>5"
46224626
else:
46234627
reason = "?"
4624-
4628+
46254629
# still run it through the averaging function when n=1 to standardize the units and run checks
46264630
kin = average_kinetics([r.kinetics for r in rs])
4627-
4628-
if n == 1:
46294631

4632+
if n == 1:
46304633
kin.uncertainty = RateUncertainty(mu=0.0, var=(np.log(fmax) / 2.0) ** 2, N=1, Tref=Tref, data_mean=data_mean, correlation=label)
46314634
else:
4632-
46334635
kin.comment = f"Blowers-Masel fit was bad ({reason}) so instead averaged from {n} reactions."
46344636
dlnks = np.array([
46354637
np.log(

0 commit comments

Comments
 (0)