EGraFFBench: Evaluation of Equivariant Graph Neural Network Force Fields for Atomistic Simulations

Q1. I have serious doubts about the experimental results in Table 1. The MD17 results of NequIP and Equiformer are much worse than reported in their original papers. Although BOTNet and MACE report rMD17 results instead of MD17 in their papers, it doesn’t make any sense that their performance on MD17 will be such bad. I don’t think authors run their code correctly.

Response: We thank the reviewer for the careful reading and critical comments. We appreciate the opportunity to address your concerns, which has indeed enhanced the quality of our work. Please find below our detailed response to the issues raised:

1. Rectification of Methodological Issues: We are grateful for your observation regarding the high MAE values. We carefully looked at the codes, log files, and the datasets, and identified an issue related to the MACE model for the MD17 dataset. Specifically, while training the models on MD17, the periodic boundary condition was kept as 'True' while keeping the box large enough to avoid any self interaction. However, the internal architectural setup of MACE (and BOTNet) requires the lattice parameter to be specifically mentioned when the PBC is set to "True", even for small molecules. To address this concern, we have now re-trained all the models on the MD17 dataset without PBC and updated all the Tables and results associated with it. We believe these revised results align more closely with the expected outcomes based on the original papers of the compared models.

2. Comparison with literature data: Further to compare our reported results with existing literature studies, we have included a new Table 19 in our revised manuscript. This table provides a direct comparison of our results with those reported in the existing literature. It can be observed that the values obtained in the present work are comparable with those available in the literature. We believe this addition will offer a clearer context and a more comprehensive understanding of our experimental outcomes.

However, it should be noted that the results cannot be directly compared with the literature in many cases as the literature does not clearly provide the exact training-dev-test dataset. To address this concern, we have included the exact training-dev-test datasets, codes, hyperparameters, config files, and log files in the repository associated with the present work. For example, DimeNET and PaiNN model in Fu et al. was trained on 9500 configurations. In our study, we have used a training set of 950 configurations across all the models for a fair comparison. The dataset size was limited due to the high computational cost in transformer models.

Overall, with the release of all hyper-parameters, train-dev-test splits, etc., we believe we present the most transparent and reproducible evaluation of GNN force fields. Nonetheless, if there exist any outstanding concerns regarding codes or results, we request the reviewer to raise any of those.

Q2. Since the molecular dynamic simulations are based on frozen models that are trained to predict energy and force, the quality of Table 1 makes the MD simulation not convincing.

Response: We agree with the reviewer that higher energy and force error in training will make the MD simulation unrealistic. As mentioned in the previous comment, we have now retrained the models for the MD17 datasets, and performed MD simulation with the new trained models and corrected the results accordingly in Tables 1, 2 and 3.

Moreover, as mentioned earlier, all the codes and datasets to reproduce to the present work are made available in the GitHub (https://anonymous.4open.science/r/MDBENCHGNN-BF68). The logfiles are also included for easy analysis. With these additional data, we hope the reviewer can evaluate the results and authenticity of the present work in a transparent fashion.

Nevertheless, if there are any outstanding concerns, we request the reviewer to raise those.

Dear Reviewer Kn6d,

Since we are into the last few hours of author-reviewer discussion, we are keenly awaiting your feedback on the responses to your comments. Specifically, we have

• retrained all the models on MD17,
• compared the results with the literature,
• included two additional models, namely DimeNet++ and PaiNN, and performed all the tasks on these models,
• released all the log files, codes, and datasets and
• addressed all other concerns raised.

With these changes, we hope the reviewer finds the manuscript suitable for ICLR as a benchmarking paper. If there are any outstanding, we request you to raise those. Otherwise, we request you support the manuscript by increasing the score. Looking forward to your response.

Q1. I believe this would be difficult to do given the amount of work it would take, but it would be really valuable to have more of an ablation style study of what particular architectural choices help / hurt in different metrics. Even with the present benchmarks, there are many architectural difference between these models that the takeaways are a bit binary -- e.g. this model is or is not enough for this task -- rather than -- e.g. this architectural choice seems to help with X. If you think that your results can support such guidance I think that would be extremely useful to the community.

Response: Indeed, we observe that models perform different across different atomic systems. This could be attributed to the differences in the nature of bonding, size of the molecules or atomic systems, and the nature of the training set in terms of how much scenario it covers. A detailed analysis of which features associated with each architecture benefits a given dataset is proposed as a future work.

To address this comment, we have now included additional text in the Limitations and future work section of the main manuscript as follows.

Another interesting aspect is the empirical evaluation of which particular architectural feature of a model helps in giving a superior performance for a given dataset or system (defined by the type of bonding, number of atoms, crystalline vs disordered, etc.). Such a detailed analysis can be a guide to designing improved architecture while also providing thumb rules toward the use of an appropriate architecture for a given system.

Q2. For each of the models, you can predict forces either as a direct prediction or via backprop to atomic coordinates. From 2.1 it seems that backprop was always used for forces for this benchmark (which makes sense given the desire for conservative forces). However, it would still be interesting to contrast the evaluation efficiency gain vs. stability loss in this setting and whether this changed substantially between methods. Do you have any runs that would give insight into this?

Response: We agree with the reviewer that a study evaluating the direct force architecture vs backprop force architecture would be insightful. However, we did not employ this approach due to two reasons as follows.

1. Being a benchmarking work, we adhered to the architectures as reported in the original works.
2. These models specifically considered backpropagation for force calculations to make them conservative as the reviewer rightly pointed out. This property is crucial to run stable molecular dynamics simulations.

Indeed, to address this comment, we attempted to directly predict forces from the model. However, several major architectural modifications would be required to achieve this in some of the models. For instance, allegro predicts the energy directly from the edge embedding, which is then summed up to obtain the total energy. Thus, obtaining the per atom force from this architecture requires a non-trivial modification of the architecture. Nevertheless, we thank the reviewer for this suggestion. We will consider this approach as part of a future study.

Q1. The work neglects some useful metrics from the previous works, like stability, RDF, diffusivity, etc in [1]. It would be better to include the performance of latest models on these metrics.

Response: Thank you for raising this concern. Indeed, previous work referred to by the reviewer uses stability, RDF, and diffusivity. Among these metrics, while stability is an evaluation of dynamics, RDF is an evaluation of the structure. In the present work, we have proposed metrics that capture both structure and dynamics in a more systematic fashion. Specifically, to evaluate the structure, we have already added two RDF-based metrics, namely JS Divergence and Wright Factor. While the first metric captures how close the distributions of predicted and actual RDF is, the second metric represents the relative error between the predicted and actual RDF. Similarly, stability of the dynamics is captured through energy and force violation. Specifically, in these two metrics, we compare the energy and force on a molecular dynamics (MD) simulation performed with the machine learned potential with respect to the ground truth trajectory. Thus, the force and energy violation error clearly evaluate the divergence of force and energy of the simulated trajectory with respect to the ground truth. These metrics, in turn, give a clear indication of the stability of the simulation.

The third metric proposed, namely, diffusion is not necessarily metric; it is rather a property of the system which is derived from the dynamics. Note that there are several such properties which can be obtained from the MD simulations. Nevertheless, in order to address the concern, we have now included RMSD (note that diffusion constant is the slope of MSD with respect to time) of the trajectories for all the datasets and models in Appendix A.9 (see Figs. 8 and 9).

Q2. Though new datasets and metrics are proposed. The major conclusion that low energy/force doesn't guarantee performing well in molecular simulations is not fresh, as pointed out in [1].

Response: We agree with the reviewer that the low energy or force doesn't inherently ensure superior performance in molecular simulations, which has also been highlighted previously [1]. However, we do not consider this to be the major contribution of the work.

The main contributions of our work are in the novel tasks proposed to evaluate equivariant graph neural network forcefields as detailed below.

1. Complex tasks: We introduce complex tasks that evaluate ML potentials for realistic applications. This involves training on diverse liquid systems with varying chemical compositions (LiPS20, see App. A.2, Table 11) and assessing performance on the dynamic simulations of unseen crystalline and disordered structures. To the best of our knowledge, this is the first time such a benchmarking task has been performed on equivariant graph neural network force fields (see Sec 4.4). Such tasks can be a critical one for evaluating any new forcefield for materials discovery.
2. Equivariant graph neural network forcefield: This is the first work that evaluates state-of-the-art equivariant graph neural network forcefields in a consistent fashion on the same dataset and on challenging tasks with rigorous structure and dynamics metrics. To further emphasize this, we show table 21 and table 10, the datasets on which previous works have evaluated the equivariant GNN forcefields. It can be observed that there is a clear lack of systematic evaluation of these forcefields.

We also thank the reviewer for noting the contribution of new datasets and metrics in our work. In addition, the key insights obtained from the work are included in the conclusion section in a point-by-point fashion.

Q3. Some benchmarking settings may not be convincing and not reflect the scenario in real applications. More discussions are included in the following Questions section.

Response: In our assessment, we've thoroughly considered the questions raised regarding the benchmarking settings. Our approach aims to encompass a diverse range of tasks and metrics to gauge performance comprehensively. We have addressed the concerns raised by the reviewer in the Questions section. We're open to further discussions to refine our benchmarking methodology and ensure a more accurate representation aligned with real applications.

I appreciate the authors' efforts in answering my questions. However, there are still some questions:

1. Indeed diffusivity is a property from MD simulations, but it could still be an important metric to evaluate ML force field. Since in most cases for MD simulations, we care about such statistical properties rather than the accuracy of each timestep.
2. I thank the authors for including the hyperparameter in A.10. However, as reviewer Kn6d pointed out, the reproduced results (as shown in Table 19) are usually worse than the number reported in original and other works, especially for the new added PaiNN and DimeNet++. So I still have the concern of underfitting.