Symmetry-Driven Discovery of Dynamical Variables in Molecular Simulations

We thank you for your review and we hope we can address some of your concerns below :

Quantitative Metrics and Baselines

We understand your concerns regarding the rigor of the evaluation criteria. Please refer to the response for Reviewer 9HNJ to view a table of quantitative results that might address some of the concerns you raised. As seen in the table, the grid of initial configurations given by the discovered DoF provides coverage of the entire conformational landscape. As a result, we are able to find at least some points within each conformational basin within only 10 simulation steps.

To address your second point about the lack of comparison with other enhanced sampling methods, we find other method to be orthogonal in nature to some of the other sampling methods available. Most enhanced sampling methods like umbrella sampling and metadynamics based methods require specified CVs in order to sample the space. As a result these methods, work on a reduced low-dimensional space unlike the our proposed method which directly works on the all atom model. While replica exchange methods do directly sample from the high-dimensional all-atom positions, they do not involve learning any local symmetries. Furthermore the replica exchange methods require communication between different replicas and cannot explore different parts of the energy landscape simultaneously. We provide a highly parallelizable method which uses local symmetry information to efficiently index and explore conformational basins. Using parallel exploration we only need to simulate <10 steps/starting configuration (total 1000 starting configuration) to get to all the known conformational basins for the explored systems. As a result, the method cannot be directly compared to the sequential algorithms like metadynamics and replica exchange models.

Scaling to Larger Systems

As mentioned in the paper, the direct optimization based methods can be applied to any system immersed in a force field. Even in the presence of explicit solvent, the procedure only needs to access the force exerted on the atoms in the molecule for given configurations. As correctly stated by the reviewer, the paper is introducing the idea of discovering local DoF, we tried to focus mostly on small systems with well-known and interpretable DoF.

Concerns Regarding Applicability to ML

We believe the symmetry discovery/optimization part of the procedure as well as the problem of transfer and generalizability of learnt DoF is an ML algorithm, an ML venue is much better qualified to judge its validity than a bio venue. Additionally, comp bio, drug discovery and physical simulation are among the fastest growing and most important use cases of ML.

Thank you for your thorough and thoughtful review. We address some points below:

Quantitative Metrics & Baselines

Please refer to our response to Reviewer 9HNJ for additional results. Although there are existing methods for CV discovery, our method considers the slightly orthogonal task of discovering local DoF. Unlike CVs, the discovered DoF are local in nature and do not necessarily correspond to any global Collective Variable. Moreover, we show that we can sample all the conformers for these molecules within 20 fs (10 steps) of starting the simulation from the sampled gridpoints. This setup is very different from the standard CV setup in metadynamics or replica exchange where the initial configuration is fixed.

Presentation

Thank you for the pointing out the notational inconsistencies. We will add a glossary of the relevant notation in the supplementary material. Furthermore, we will also add a section on illustrative examples applying the proposed method on small systems. Please refer to our response to Reviewer p5r2 for a discussion of a small synthetic system similar to the one you proposed. As the Hessian based method in

Ablation Studies

We will add the suggested ablation studies in the supplementary material. It is not always possible to plot the trajectories for the highest valued eigenvectors of the Hessian as the transformations can lead to unnatural states that lead to errors in the openmm simulation. However, we will add a more detailed discussion in the paper.

Local CV and Relation to other CV discovery work

We thank the reviewer for the reference [1], but the CVs discussed do not align very well with the present work. Although one could define local variables $z_i = L_i x_0$, these resulting variables behave differently from collective variables. We treat all configurations that can be reached from $x_0$ as an equivalence class; thus, all the z_i's would form a basis for this subspace. Thus, all the local $z_i$'s can be used to index the local conformation basin, but it is unclear if these can be patched together or if patching these would give a good global CV.

Clarifications for Typos

• Page 1: "both computationally efficient and physically insightful." -> what do you mean by physically insightful? We mean it provides insight into the effective dynamics of the system.
• Page 3: I find the notation for general linear group strange. … $GL(\mathcal{X})$ doesn't make sense to me, … not consistent with the subsequent usage of (special) Euclidean groups. Actually, $GL(V)$ is the standard way the general linear group for a vector space $V$ is denoted, though it may not look as familiar as $GL(n,R)$. We will explain and add the definition of $SE(n)$.
• Page 3: The Orbit definition seems trivial as written … Do you mean to consider a subgroup? You are correct for $\mathcal{X} \sim R^n$. We are indeed usually concerned with a subgroup (the approximate symmetries of $E$). Will reword it.
• Page 3: symbols in "Group parameters as DOF" is quite unclear to me.
  • Lie algebra basis $L_a$, is this a single basis element, or does $a$ run through … the Lie algebra? $L_a$ is a single basis element.
  • For the notation $\theta \cdot L$: Yes, it is $\sum_a \theta_a L_a$
  • If θ is a vector, the Taylor expansion error term should be O(|θ|2) Any combination $\theta_a \theta_b$ is second order, but, yes, one can pull out the norm of $\theta = |\theta| \hat{\theta}$ and write the error as $O(\theta^2)$.
  • Since … always GL, what's the point of this abstraction to Lie groups/algebras? Why not just identify the Lie algebra with $M_n(R)$? They are matrices, yes, but we want a subgroup
  • Boldfaced $I$ vs normal. Any difference? No, will fix it.
  • Overall, I suggest … properly define the notation. Thank you, will do.
• Page 4: Equation 1 is not clear, is $S$ considered a subset of $GL(n,R)$ or $M_n(R)$? It should be the latter… We should clarify in the paper that we are working with the canonical representation of $GL(n,R)$ which maps elements to a subset of $M_n(R)$. So, yes, all in $M_n(R)$.
  Also, what's the dependence of δE on DOF? Why not just use $E_{barrier}$ which is already used in the previous page and not introduce an extra η? Eq 3, $\delta E = \theta \nabla E \cdot Lx$ (here $\theta \in R$) shows dependence of $\delta E$ on DOF $L$. We can replace $\eta$ with $E_{barrier} in eq 1.
• Page 4: the variables $\theta_l,\theta_s$ do not appear in eqs 1 and 2. Sorry, we changed it to $\epsilon_l, \epsilon_s$, will fix.
• g≈I+θL now takes θ as a number instead of a vector. A bit of abuse of notation here, will fix. Here $L$ denotes an arbitrary vector in the Lie algebra instead of specific basis elements $L_a$. Then, $\theta \in R$ is a small $\epsilon$ making clear that $g$ is near identity.
• P6: Eq 15, $K= n_L$ as in eq (14)? Yes, will fix.

[1] Understanding recent deep-learning techniques for identifying collective variables of molecular dynamics

We thank the reviewer for their comments, and we try to address some of the concerns below:

More Quantitative Results

As pointed out by the reviewer, some of the results presented in the paper might look a bit qualitative. Here, we provide some more quantitative statements about the results shown in the paper. We intend to add this table in the main paper in order to better illustrate the effectiveness of the discovered DOF.

Long Simulation Results: The table below shows the frequency of observing the molecule in a specific conformation as a function of the simulation time ( simulation is run with a step size of 2 femtoseconds). We consider a molecule within a confirmation basin if the dihedral angles are within the known range for the specific conformer.

DOF-assisted Simulation Results: Under our method, we only run a simulation of 2 picosecond duration for every starting configuration provided by our DOF discovery algorithm. Thus, the conformers do not have enough time to escape from one conformational basin to another. Thus, the initial distribution of the conformers stays the same throughout the short simulation. Here, we provide the frequency of finding a molecule within a given conformational basin (within 20 ps of simulation) given the initial 1000 starting configurations.

Simulating the algorithm on additional smaller systems

As suggested by the reviewer, we intend to add some example applications of the algorithm for smaller analytical systems. Please check bullet point 2 in our response to Reviewer p5r2 for further details and a small example system.

Literature Review for Hessian-based methods

We thank the reviewer for pointing us toward this literature. We will include a detailed overview of the work in the literature and a comparison with the current paper. For geometric structure optimization, most of the methods in the literature build upon the idea of the Newton-Raphson method to accelerate gradient descent using second-order Hessian information to provide quadratic convergence near local minima by providing faster updates without calculating the full Hessian. On the other hand, transition state-finding methods utilize information about the eigenvector corresponding to a negative eigenvalue of the Hessian to find saddle points. Some synchronous transit methods like LST and QST also use Hessian information to better estimate the energy landscape near intermediate states between two final reaction states. While conformer discovery shares some of the same objectives as geometric structure optimization and transition state finding, unlike geometric structure optimization, we need to escape local energy basins to find a diverse set of conformers. Unlike the transition state finding task, we do not have a specified final reaction state that we intend to reach. Our primary objective is to discover certain local degrees of freedom to explore the energy landscape and find more conformers efficiently. In our method, we discover local DoF for a given conformer and use the discovered DoF to explore the basin near the conformer without any additional knowledge of CVs.

Better Presentation of the Methodology

As suggested by the reviewer, some of the methodology is quite hard to understand from the current presentation. We plan to reduce some of the jargon used in the paper and add a table of references in supplementary materials to help with the notation used in the paper. Additionally, we will add a section with some illustrative examples of applying our proposed methods to small systems. We present the following 1-d, 2-particle system suggested by Reviewer zCrM as an example. Due to space constraints, we only show the discovered local DoF using the Full Hessian method. We leave a more detailed discussion for larger systems in the paper. Consider the energy function: $E(x) = F(|x_1| - x_2) + \epsilon \sin\left(\frac{x_2}{\epsilon}\right)$ for $x \in \mathbb{R}^2$, where $F$ is a 1-D potential. In this case, we can use the first order optimality condition to see that at any minima, $F(|x_1| - x_2) = 0$ and $\cos\left(\frac{x_2}{\epsilon}\right) = 0$. Then, we see that the Hessian at any such minima, $x'$, can be given as