Tautomer Search

Many medicinally relevant molecules exist as multiple tautomers, and understanding which tautomer predominates can be key to subsequent computational tasks: for instance, Hu et al. found in 2016 that correct assignment of tautomeric state dramatically improved relative binding affinity predictions. Rowan's tautomer prediction workflow uses machine-learned interatomic potentials to enable fast and minimally empirical prediction of the relative stability of different tautomers.

Rowan's tautomer workflow proceeds in three stages: enumeration, screening, and refinement.

Settings

The conformer generation and refinement steps can be set independently, exactly as in the standalone conformer search and multistage optimization workflows:

• Conformer generation settings control how trial conformers are generated for each significant tautomer. Rowan supports several conformer generators, including ETKDG (RDKit), iMTD-GC (CREST metadynamics), and openconf (a hybrid ETKDG/torsional Monte Carlo scheme), each with its own tunable parameters. By default, Rowan uses openconf with max 20 conformers per tautomer.
• Multistage optimization (MSO) settings control how those conformers are optimized and scored during refinement: one or more optimization steps, optionally at increasing levels of theory or in solvent, followed by a final single-point energy. By default, Rowan optimizes with AIMNet2 and implicit solvent water (using ALPB for the optimization and CPCM-X for the single-point energy). A g-xTB // GFN2-xTB preset is also available, and any other combination of levels of theory can be selected manually.
• Final correction optionally replaces the MSO solvation treatment with a COSMO-RS solvation free-energy correction. COSMO-RS requires a DFT surface calculation and is therefore more expensive than the default CPCM-X-based treatment.

Enumeration

Rowan enumerates tautomers using the tautomer enumeration functions in RDKit and a separate function that looks for same-net-charge protonation microstates (including zwitterions) that RDKit's enumerator can miss.

Screening

Each proposed tautomer is first optimized with GFN2-xTB in implicit solvent (chosen to match the MSO settings) and then scored with a single-point energy using the final level of theory in the MSO settings. Tautomers more than 10 kcal/mol above the lowest-scoring candidate are discarded, and at most 50 "significant" tautomers are carried forward.

Refinement

For each significant tautomer, Rowan runs a conformer search using the requested conformer generation and MSO settings. A vibrational frequency calculation is performed on each optimized conformer to obtain a thermal free-energy correction.

When COSMO-RS is selected as the final correction, Rowan applies it to at most the three lowest-MSO-free-energy conformers for each tautomer. The correction uses the solvent selected for the MSO single-point calculation, or water if that calculation is gas phase. Each corrected conformer free energy is

The DFT COSMO surface supplies the solvation correction; the electronic-energy anchor remains the final single-point method selected in the MSO settings. The conformers are then Boltzmann-weighted into a single free energy for each tautomer, and all tautomers are Boltzmann-weighted against one another at 298.15 K to obtain relative energies and populations.

Accuracy

(The following benchmark was run using our old default settings: an ETKDG conformer search and AIMNet2/CPCM-X(Water) // GFN2-xTB refinement.)

On the aqueous subset of the TautoBase benchmark set, Rowan's tautomer workflow displays a mean absolute error of 2.10 kcal/mol and a root mean squared error of 2.99 kcal/mol. This is comparable to the performance of high-level quantum chemical methods reported by Chodera and co-workers: B3LYP/cc-pVTZ/SMD(water) was reported to give an RMSE of 3.1 kcal/mol vs. TautoBase (on a slightly smaller subset).

A more relevant benchmark for real-world usage is classification accuracy—how much of the time can Rowan predict the correct lowest-energy tautomer? On the full dataset, Rowan predicts the correct tautomer 89% of the time. Some of these comparisons are not particularly challenging: for compounds with an experimental ∆∆G of less than 3 kcal/mol (shown in red), Rowan is still correct 77% of the time.