Improve group additivity ranking for aromatic species#1731
Merged
Conversation
Group values are often fitted to the most aromatic resonance structure of a molecule. This can cause issues if the thermo value for a less representative structure ends up being lower than that of the aromatic structure, even though it has a more accurate group estimate. This approach ensures more predictable behavior after adding new group values for aromatic species.
Codecov Report
@@ Coverage Diff @@
## master #1731 +/- ##
==========================================
+ Coverage 32.57% 32.64% +0.06%
==========================================
Files 87 87
Lines 26090 26170 +80
Branches 6866 6877 +11
==========================================
+ Hits 8500 8542 +42
- Misses 16627 16656 +29
- Partials 963 972 +9
Continue to review full report at Codecov.
|
alongd
reviewed
Sep 27, 2019
alongd
left a comment
Member
There was a problem hiding this comment.
Looks straight forward and OK.
I added one commentt re commit messages.
Would you like to add an actual H298 GAV test? Maybe for the species in your example?
|
|
||
| cpdef identify_ring_membership(self) | ||
|
|
||
| cpdef int count_aromatic_rings(self) |
Member
There was a problem hiding this comment.
Could you change the commit message to lower_case_with_underscores to be consistent with the code?
(same comment for the Minor code improvements to ThermoDatabase.prioritizeThermo commit)
mliu49
force-pushed
the
aromatic_gav
branch
3 times, most recently
from
September 27, 2019 21:05
701fb2f to
276a91c
Compare
Merged
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Motivation or Problem
Estimating thermo for aromatic species using group additivity can sometimes give weird results because we will get an estimate for every resonance structure and then take the one with the lowest H298. In some cases, this can lead to severe underprediction of enthalpy if there is an inappropriate group which results in a lower enthalpy than a more appropriate group.
For example:
This species has resonance structures where the radical can delocalize into the aromatic ring. One of those structures gives the estimate with the lowest H298 of 61.4 kcal/mol:
Thermo group additivity estimation: group(Cs-(Cds-Cds)CbHH) + group(Cb-Cs) + group(Cb-(Cds-Cds)) + group(Cds-Cds(Cds-Cds)Cs) + group(Cb-H) + group(Cds-CdsCbH) + group(Cb-H) + group(Cds-Cds(Cds-Cds)H) + group(Cb-H) + group(Cb-H) + group(Cds-CdsHH) + group(Cdd-CdsCds) + polycyclic(s2_6_6_ben_ene_1) + radical(Benzyl_S_dihydronaphthalene)However the calculated H298 of the molecule is 101.7 kcal/mol. If we were to use GAV to estimate the thermo of resonance structure depicted above, then we would get 107.3 kcal/mol with the following groups:
Thermo group additivity estimation: group(Cbf-CbCbCbf) + group(Cbf-CbCbCbf) + group(Cb-(Cds-Cds)) + group(Cb-H) + group(Cb-H) + group(Cb-H) + group(Cb-H) + group(Cb-H) + group(Cb-H) + group(Cb-H) + group(Cds-CdsCbH) + group(Cds-CdsHH) + polycyclic(s2_6_6_naphthalene) + radical(Cds_S)Description of Changes
This PR adds an additional parameter, number of aromatic rings, by which to rank GAV estimates, in addition to H298.
The number of aromatic rings is determined by the actual bond orders in the resonance structure, rather than true aromaticity detection.
The result is that estimates for resonance structures with more aromatic rings is prioritized over those with fewer, even if it results in a higher H298.
Testing
Besides testing with the species mentioned above, I ran comparisons against the SABIC_aromatics_1dHR library (not currently on master), which contains ~350 aromatic/near-aromatic species calculated using CBS-QB3 with 1D hindered rotors.
The H298 MAE on master is 5.263 kcal/mol, while the MAE on this branch is 5.228 kcal/mol. This is not a substantial difference, mainly because this change does not affect many species. However, it does suggest that this change does not have noticeable negative effects.