Learning Over Molecular Conformer Ensembles: Datasets and Benchmarks

Response to Reviewer pme5 (Part 1/3)

Thanks for your detailed and constructive reviews! We would like to make the following clarifications regarding your main concerns. Due to the length constraint on the comment box, we split our response into several parts. *W1. In the four tasks developed in this study, the objective is defined as predicting the Boltzmann average of various properties over multiple conformers. Under this goal setting, it seems intuitive that using information from multiple conformers would naturally improve prediction accuracy. Therefore, it has not been proven that ‘considering multiple conformers contributes to machine learning predictions of real data (e.g., actual experimental measurements of molecules rather than computed values).’ In this sense, the utility of this benchmark remains a bit artificial, and practical values would be unclear.* A: Regarding the practical relevance of considering multiple conformers, we note that to the best of our knowledge, the utility of encoding conformer ensembles in deep learning models for real-world molecular property prediction tasks has not *yet* been convincingly demonstrated. For instance, Axelrod and Gómez-Bombarelli [1] trained multiple models that encode conformer ensembles from the GEOM-DRUGS dataset in order to predict experimental protein-ligand biological activity. However, their ensemble models did not improve significantly upon models that encode a single 3D conformer. **This emphasizes the need for a machine learning benchmark that is explicitly designed to help conformer ensemble models reach their full potential, so that they can eventually be deployed in real-world tasks relevant to drug discovery, computational chemistry, etc.** **Although our benchmark includes simulated data, we view this as necessary to ensure the quality of the benchmark so that model analysis is not obfuscated by the noise in real-world tasks involving heterogeneous experimental data. Our benchmark is thus crafted to curate a set of well-defined, carefully controlled tasks to enable rigorous model benchmarking.** *Crucially*, the failure of prior deep learning models to effectively exploit the extra structural information contained in conformer ensembles *does not* mean that conformer ensembles are not useful for real-world tasks. We would like to emphasize that even *simulated* Boltzmann-averaged properties, while not typically included in GNN benchmarks (e.g., MoleculeNet), are still of substantial utility for the drug discovery and computational chemistry communities. For instance, simulated Boltzmann-averaged properties are often employed in computational chemistry to predict chemical reactivity [2-3], to aid in the discovery of new catalysts [4], and to computationally approximate protein-ligand binding affinities [5]. [1] Molecular Machine Learning with Conformer Ensembles. https://iopscience.iop.org/article/10.1088/2632-2153/acefa7/meta [2] AARON: An Automated Reaction Optimizer for New Catalysts. https://pubs.acs.org/doi/10.1021/acs.jctc.8b00578 [3] Multi-Instance Learning Approach to the Modeling of Enantioselectivity of Conformationally Flexible Organic Catalysts. https://pubs.acs.org/doi/10.1021/acs.jcim.3c00393 [4] Conformational Effects on Physical-Organic Descriptors: The Case of Sterimol Steric Parameters. https://pubs.acs.org/doi/pdf/10.1021/acscatal.8b04043 [5] End-Point Binding Free Energy Calculation with MM/PBSA and MM/GBSA: Strategies and Applications in Drug Design. https://pubs.acs.org/doi/10.1021/acs.chemrev.9b00055

Response to Reviewer pme5 (Part 2/3)

*W2. It is possible that machine learning predictions based solely on the ground-state structure, as traditionally done, are already practically useful. Regarding how considering multiple conformers contributes to molecular machine learning, the contribution of this study might be limited.* A: Thank you for your valuable comments. We agree that traditional machine learning predictions based on the ground-state structure have been useful, especially for benchmarking 3D graph neural networks (e.g., with the QM9 dataset). However, we emphasize that few *experimentally observable* properties are exactly dependent on a single static conformer: most experimental measurements are implicitly a thermodynamic average over the accessible conformer distribution. Moreover, it is often not practical to precisely determine the global minimum energy conformer of relatively complex molecules in high-throughput property prediction scenarios, since determining the lowest-energy conformer of all but the tiniest molecules requires extensive geometry optimizations across a large conformational space. In contrast, averaging properties over a representative conformer ensemble can be more robust (in addition to being more accurate) than relying on the lowest-energy conformer, as missing the global-minimum conformer will have less of an influence on the Boltzmann-average. Therefore, our approach can better align with real-world use-cases, where predicting properties based on conformer ensembles more holistically represents physical phenomena and better represents the feasibility of real-world studies. *Q1. Additionally, since all four datasets prepared are secondary data from referenced primary data, it’s unclear how challenging it would be for researchers to prepare them on their own. While the paper mentions quality control and the removal of redundancies, it’s unclear whether there is any original information added to this study.* *While Drugs-75K, a subset of the GEOM-Drugs [27], Kraken [33], EE [34], BDE [36] all have associated citations, indicating they are secondary data, it was unclear whether the presented datasets are a simple curated version of these primary data, or if any original information was generated in this study. If there is any original information, it would be very helpful to be clarified.* A: Thank you for your constructive reviews. First of all, we would like to clarify that historically it is not yet shown that conformer ensembles could improve performance on real-world molecular tasks. Therefore, our benchmark should consist of well-defined, simulated tasks without additional confounding variables to ensure a fair comparison with existing baseline models, and to assist in the development of new models that can more effectively make use of extra structural information contained in conformer ensembles. Secondly, regarding the added information to the datasets, the Drugs-75K dataset requires significant computational effort. The original GEOM-Drugs dataset was constructed using semi-empirical Density Functional Theory (DFT) methods, which is less accurate than full DFT. To curate the Drugs-75K subset, we generate the conformer ensembles with Auto3D and compute their corresponding energies using AIMNET-NSE. This process took approximately 600 CPU hours on an AMD EPYC 7763 server. Last, the datasets compiled in MARCEL have been reorganized and standardized to form a unified benchmark, requiring substantial effort as well. We have developed a Python interface for easy access of these datasets into PyTorch, facilitating their use and extension by the research community.

Response to Reviewer pme5 (Part 3/3)

*Q2. Regarding ‘Dataset preparation’ in the Supplementary Material, why are different methods used to generate conformers depending on the data? (Auto3D for Drugs-75k, Q2MM for EE, Open Babel with DTF?) Is this point not problematic for benchmarking several methods?* A: Thank you for your careful reviews! In our benchmark, the choice of method for generating conformers was carefully considered to best suit the specific characteristics of each dataset. For example, Auto3D was chosen for Drugs-75K because of its efficiency in generating high-quality conformers. Similarly, Q2MM was employed for the EE dataset due to its ability to accurately generate Transition State Force Fields (TSFFs). It is important to note that these tools are not only well-established but have also been extensively validated within the chemistry community. It is also important to note that many real-world tasks in the chemistry community employ different conformer generation workflows, and hence our benchmark properly reflects this diversity of conformer generation strategies. Also, all baseline models are evaluated using the same set of datasets. This ensures a fair and direct comparison across different methods. The diversity in conformer generation does not detract from the comparability of results but rather enriches the benchmark by covering a broader spectrum of realistic scenarios in molecular machine learning. *Q3. While I understand that the ‘two conformer ensemble learning strategies’ are useful for predicting the ‘Boltzmann-averaged value of each property across the conformer ensemble,’ can they be said to be generally useful for predicting molecular properties? Could you provide any supported evidence for this claim?* A: Please refer to our response to W1. *Q4. Is the BDE task also about predicting the Boltzmann-averaged value?* A: Thank you for your valuable feedback regarding our paper. You raised an important question about whether the Binding Energy Difference (BDE) task in our study also involves predicting the Boltzmann-averaged value. Indeed, a Boltzmann-averaged BDE would more holistically capture physical reality, as it would account for the full distribution of thermodynamically-accessible conformers under experimental conditions. However, in our case, the BDE dataset was computed with a simpler approach to manage the computational cost of DFT calculations. In the BDE dataset, binding energy difference is computed from single-point energy calculations (at a high level of DFT) between the two lowest-energy conformations of the unbound/bound complexes. These lowest-energy structures (one for each of the unbound and bound structures) were identified following a conformer search with OpenBabel. The lowest energy conformers then underwent further geometry optimization with DFT at the B3LYP/3-21G and B3LYP/def2-SVP levels, followed by a single-point energy calculation at the more expensive def2-TZVP level. The regression label is the (single-point) energy difference between the unbound and bound structures. The input conformers encoded by our models, however, are the original OpenBabel conformers prior to conformer selection and geometry optimization. Please note that the BDE is still a function (albeit indirect) of this conformer ensemble, and thus nicely fits within the MARCEL benchmark. The use of OpenBabel conformers as input to the models more realistically represents the real-world scenario during inference in which expensive DFT calculations are not practical. We will revise the relevant sections to ensure this distinction is explicitly stated and understood, thus preventing any potential confusion.