Thank you for strongly championing our work! We are happy to further discuss the justification and other possible views for our method.

---

Justification for Random Placement of Small Molecules

Randomly placing solvent molecules around a solute is a well-established strategy in molecular simulations. For instance, protein–ligand docking protocols (e.g., Rosetta) often initialize ligands in random orientations relative to the protein before searching for binding modes. Similarly, Monte Carlo insertion methods like Widom’s test-particle approach randomly insert solvent molecules to explore solute–solvent configurations without bias. This stochastic initialization allows efficient sampling over a wide range of distances and orientations while remaining physically sound: in the limit of sufficient sampling, no unphysical configuration is favored, and the process mimics the early stages of solvation when solvent molecules approach from arbitrary directions.

Our approach constructs a “virtual interaction environment” by generating many such random solute–solvent geometries. This approximates the diversity of encounters seen in explicit solvation, while avoiding the high cost of molecular dynamics. It also exposes the model to both favorable and unfavorable orientations, encouraging generalization across interaction types. The method is computationally efficient and theoretically grounded, supported by literature precedent.

---

Local Relative Geometry as an Interaction Force

While we define local relative geometry for learning fine-grained interactions between molecules, we can view local relative geometry as an interaction force between molecules. More specifically, we define a local relative geometry using the normalized vector between a solute atom and its closest solvent atom.

This provides a physically motivated supervision signal rooted in classical intermolecular forces, many of which are central and act along the internuclear axis. For example, van der Waals interactions (described by the Lennard-Jones potential) exhibit repulsive or attractive forces directed along this axis. At short distances, repulsion dominates and aligns directly outward between nuclei.

This supervision scheme serves as a central-force approximation, consistent with classical force fields, and offers a lightweight surrogate for full force labels, which would require costly quantum chemistry or MD simulations. Notably, SchNet [1] demonstrated that even approximate force signals improve learning of molecular interactions. Our direction-based supervision enables the model to learn geometric features like hydrogen bond alignment or steric repulsion trajectories in an SE(3)-equivariant manner.

Since solvent atoms are placed near specific solute atoms, the dominant interaction direction aligns with the interatomic vector, making it a reasonable proxy for the net force axis. Thus, the unit direction vector serves as a pseudo-force label, conveying the primary interaction axis and encouraging the model to encode directionality of intermolecular interactions.

[1] Schütt, Kristof, et al. "Schnet: A continuous-filter convolutional neural network for modeling quantum interactions." Advances in neural information processing systems 30 (2017).

---

Limitations

Although the limitations of our work are addressed in Section 5.3, we will add a separate section in the camera-ready version to specifically detail concerns regarding the physicality of our virtual environment

Thank you for acknowledging our main contributions to the paper: replacing expensive quantum calculations with cost-effective 3D virtual environment simulation!

---

W1

We agree with the reviewer that relying on ETKDG-generated conformers, while computationally efficient, does introduce some inaccuracy. As you rightly noted, using experimental structures could indeed reduce this error; however, their acquisition comes with substantial expense and often limits the scalability of such methods.

Obtaining experimental 3D structures for molecules, especially for a large and diverse set, is a significant challenge. Techniques like X-ray crystallography or NMR spectroscopy are highly resource-intensive, requiring specialized equipment, significant time, and often large quantities of purified samples. Many molecules are difficult to crystallize or are unstable in conditions required for these experiments. This makes it impractical to generate experimental 3D structures for the vast number of molecules needed for large-scale machine learning model training and evaluation.

A core strength of our work, as you pointed out, is our aim to replace costly operations in traditional QM approaches with AI, thereby reducing overall experimental costs. Introducing the requirement for expensive experimental structures would, unfortunately, create another bottleneck, thereby limiting the practical applicability of our approach to a broader range of molecules and scenarios.

Nevertheless, we concur that incorporating experimental structures, despite the cost and difficulty, could further enhance our model's accuracy. This remains a valuable direction for future work, particularly in contexts where such high-fidelity data is available.

---

W2

We appreciate the reviewer's comment regarding the perceived lack of a formal mathematical proof for SE(3) equivariance's necessity. We want to clarify that the necessity of SE(3) equivariance in molecular modeling isn't about deriving a mathematical theorem, but rather about aligning our models with a fundamental and inherent property of the physical world.

Specifically, molecules exist in 3D space, and their intrinsic characteristics—such as energy, bond lengths, or interaction forces—are fundamentally independent of their arbitrary position or orientation. If we rotate or translate a molecule in space, its physical nature remains unchanged. For a computational model to accurately capture these true molecular properties, its internal representations and final predictions must transform predictably (or remain invariant, for scalar properties) when the input molecule undergoes such rigid transformations. This is precisely what SE(3) equivariance mathematically enforces and represents, as demonstrated in many previous works [1, 2].

[1] UNKE, Oliver, et al. SE (3)-equivariant prediction of molecular wavefunctions and electronic densities. NeurIPS 2021.

[2] HOOGEBOOM, Emiel, et al. Equivariant diffusion for molecule generation in 3d. ICML 2022.

---

Q1

Yes, we fully agree that incorporating molecular dynamics (MD) trajectories into our pre-training data would significantly enhance the realism of our virtual environments and further improve our work. Beyond just enriching the environment, this would also open avenues for new pre-training strategies, such as predicting molecular movement rather than just relative geometry.

While acknowledging that calculating extensive MD trajectories introduces a computational bottleneck (as we discussed in our response to W1), our immediate research interest is indeed focused on integrating this valuable MD data into our virtual environment framework. We are truly grateful for your insightful suggestion, which highlights a crucial direction for further enhancing our model's realism and capabilities.

We appreciate your strong endorsement of our research! It's rewarding to know that our virtual interaction environment construction, ablation studies, and sensitivity analyses resonated positively with you.

---

Q1

As detailed in Lines 189-191, our approach generates the target geometry anew in each epoch. Crucially, we utilize random Gaussian noise directly without any training on it. While learning this random noise could indeed be a promising research direction, it presents a significant challenge: we lack the ground truth geometry that would typically be obtained through methods like MD simulation. Without such a signal, the noise could easily degenerate to a zero vector, as our model is designed to predict this noise, making a zero vector an "easy target" for the local geometry learning described in Section 4.2.2.

However, if we could obtain MD trajectories of molecular interactions, we could first train a virtual interaction geometry generator to mimic these trajectories. This generator could then be used during our pre-training strategy instead of purely random generation, providing a more informed signal. We're truly grateful for your insightful question; it highlights a crucial direction for further enhancing our model's realism and capabilities.

---

Q2

As described in our answer to W1, we do not involve any training of noise during environment construction. This key feature means it scales effectively with all involved parameters, irrespective of their quantity.

While pre-training can encounter scaling issues as molecule and environment sizes increase, this concern is isolated to the training phase, not the environment construction itself.

---

Q3

Yes, we generally recommend setting the number of small molecules, $n$, to 5. As illustrated in Figure 4, we found that this value strikes an optimal balance between computational cost and model performance.

You're right that the optimal $n$ could depend on factors like the size of the molecules. We could dynamically control $n$ by calculating the size ratio between the larger and smaller molecules. However, doing so would introduce an additional hyperparameter that users would then need to tune, potentially adding complexity to the model's application.

However, we believe that exploring how to dynamically learn the environment based on experimental conditions or specific solvent types would be a very promising direction for future work. This would allow our model to adapt to a wider range of realistic scenarios, moving beyond static configurations to more accurately capture the complexities of molecular interactions in diverse experimental settings.

---

Q4

As detailed in Section 5.1, our model accounts for diverse interaction types, including solute-solvent, chromophore-solvent, and drug-drug interactions. Notably, 3DMRL also handles drug-drug interactions, where interacting molecules are of similar size. This suggests our method is well-suited for modeling interactions between molecules of comparable dimensions.

Currently, our experiments are limited to interactions between two partners. While we believe extending our virtual environment to accommodate more than two distinct interacting molecules (by arranging them without overlap) is straightforward, our 2D encoder parts primarily model pairwise interactions using cross-attention. To effectively model higher-order interactions involving more than two types of molecules, we would need to develop new methods.

Thank you for suggesting us highly promising future work, which is crucial for a deeper understanding of complex molecular systems!

---

Q5

We believe our method can conceptually be extended to large biomolecules like proteins. However, scaling to such large molecules presents practical challenges. Since proteins comprise a significantly greater number of atoms, the resulting virtual environments would become prohibitively large. To address this, we would need to develop specialized methods for constructing virtual environments using larger components than individual atoms, such as amino acid residues. This would be a crucial step for applying our approach to proteins and other large biomolecules.

---

Q6

Your question regarding the sufficiency of random virtual environments for modeling real-world interactions, and the performance delta compared to expensive QM methods, is insightful. As we discuss in Section 4.1, our goal was to mimic real-world molecular interactions through our virtual environment. The results in Table 3 demonstrate that our virtual environment strikes an effective balance between computational cost and the realism of the generated environment.

On the other hand, we expect that integrating our method with QM-derived data would lead to significant performance improvements. For instance, instead of purely random generation, we could use MD snapshots as initial configurations for our virtual environments, providing a more physically accurate starting point. Furthermore, leveraging entire MD trajectories as pre-training data would allow our model to learn from much richer and dynamic interaction landscapes. While these QM/MD-driven approaches inherently involve a higher computational cost, we believe they offer a path to even more realistically model interaction environments, thereby substantially enhancing our method's predictive power.

---

Limitations

Thank you for suggesting additional limitations of our work. We will make sure to incorporate the 1) physicality of the virtual environment and 2) computational complexity as additional limitations.

---

Paper Formatting Concerns

Thank you for your detailed suggestions! We will ensure that all suggestions are included in our camera-ready version.

Thank you for asking for this clarification. We apologize for any confusion our previous response may have caused.

Regarding Random Gaussian Vector

As the reviewer correctly pointed out, we construct our virtual environment with random Gaussian noise, i.e., $R^{2, I}$, which is also used to generate $\hat{f}$ in eq (4). However, the original geometry of larger molecule is ground truth geometry obtained from RDKit ETKDG algorithm without any added noise. Therefore, in Eq (5), our model learns to predict the optimal direction for the smaller molecules to move, simultaneously learning the geometry of larger molecules.

As reviewer mentioned, $f_{k}^{i}$ contains the $R^{2, i}$. However, the crucial point is that the model input is also constructed based on the same random noise. Therefore, our model is not predicting something random. Instead, it is learning the deterministic relationship between the noisy input and the corresponding noisy label, which is the geometry of larger molecules. This task is essentially a form of denoising, where the model learns to identify and predict the underlying ground truth structure by mapping a noisy input to its target.

Denoising in Machine Learning. Note that adding noise and reconstructing it is widely used in various machine learning field. For example, Diffusion models follow a 'forward process' where noise is gradually added to data to transform it into a simple Gaussian noise distribution. They then learn a 'reverse process' to denoise and restore the original data, thereby learning the complex data distribution. This process forces the model to understand the essential structure of the data from a noisy input and learn how to remove that noise.

Regarding ground truth geometry

Moreover, yes, we do not have access to ground truth geometry of whole interaction. That is, we only know the ground truth geometry of each individual molecule, but do not know how they are positioned in the space during the interaction. And in L252, we mentioned it as “ground truth relative geometry”, which is defined in L253-255. What we meant here was ground truth of relative geometry rather than the whole interaction geometry.

Thank you again for championing our work and for helping us clarify this important point. We will make sure to update the paper to explain this concept more clearly, and happy to discuss further if you have remaining questions!

Thank you for acknowledging our effort in extensive experiments of 40 downstream tasks and OOD generalization ability of 3DMRL! We're happy to address each of the weaknesses and questions that the reviewer raised.

---

W1 and Q1

We acknowledge the reviewer's concern that random Gaussian placement, even with collision checks, might lead to unrealistic interactions and that our virtual environment lacks physical validation. As we detailed in Section 5.3, we're aware that some molecule collisions can still occur.

Our approach prioritizes an efficient balance between computational complexity and performance, as demonstrated in Table 3. Regenerating the environment every time a collision occurs would increase computation time by over four times compared to 3DMRL. Furthermore, as shown in Figure 4, computational complexity significantly increases with the number of solvent molecules ($n$). Considering all possible realistic inter-atomic distances between molecules for every generated configuration would incur an even greater, potentially prohibitive, computational cost. Therefore, while we fully agree on the importance of a physically realistic environment, this was a necessary design choice given the computational constraints and the practicality of our proposed framework.

Regarding the question of how often generated configurations violate typical inter-molecular distance/angle distributions seen in MD trajectories, we couldn't quantitatively assess this as our datasets do not include MD trajectories. However, we truly appreciate you suggesting this promising future work. If we could obtain MD trajectories of molecular interactions, we would first train a virtual interaction geometry generator to mimic these trajectories. This generator could then be used in our pre-training strategy instead of purely random generation, providing a more physically informed signal. Your insightful question highlights a crucial direction for further enhancing our model's realism and capabilities.

---

W2

We agree with the reviewer that relying on ETKDG-generated conformers, while fast, introduces some inaccuracy. As discussed in Section 4.1, this inherent limitation of computationally generated conformers often leads to the underperformance of previous 3D encoder-based models when compared to methods like 3DMRL, as further evidenced in Appendix E.5.

However, it's crucial to recognize that this is a practical challenge faced by most 3D molecular modeling approaches: obtaining experimentally accurate 3D conformers for the larger molecules is prohibitively expensive and often unfeasible. Therefore, any evaluation of 3D encoder-based models resorts to the datasets that consist of small molecules such as MD17 and QM9.

In contrast, our approach explicitly addresses this challenge through its two-step design. While our pre-training leverages the 3D structural information, our fine-tuning step relies primarily on the more robust and readily available 2D molecular topology. This design makes our model inherently less sensitive to the imperfections of computationally generated 3D conformers during the critical downstream task, allowing for more generalizable and practical applications.

---

W3

We agree with the reviewer that our global contrastive loss shares similarities with Infomax. However, Infomax is designed exclusively for single-molecule property prediction and cannot learn interaction geometry. To overcome this, we developed virtual interaction geometry, an approach that, to our knowledge, is entirely novel. Beyond global contrastive loss, we also introduce novel local relative geometry learning. This aims to capture localized geometry, which is crucial for understanding how molecules interact in diverse environments. We want to emphasize that the true novelty of our work lies in these innovative methods for modeling molecular interactions.

---

W4 and Q2

We understand your concerns about the efficiency of our method due to its inherent complexity. However, we'd like to clarify that certain aspects of our approach are highly efficient, while others represent known bottlenecks, primarily from components borrowed from prior work.

As described in Lines 189-191, our virtual environment construction is fast because it relies solely on efficient matrix operations for transitions and rotations. This real-time generation capability is a key strength. Furthermore, as shown in Table 3, 3DMRL maintains a similar computational time across various virtual interaction environment construction variants while consistently demonstrating superior performance on downstream tasks.

On the other hand, the pre-training phase does see an increase in complexity as the size of the virtual environment grows. Specifically, Figure 4(a) illustrates that computational time increases with the number of smaller molecules ($n$) attached to the target atoms. This is primarily because the $f_{3D}$​ encoder model, which we utilize, requires more computation as the input size expands.

To further validate this, we performed time complexity analyses on ablated models. Our findings, presented in the table below, show that the "Only Glob." model, which requires 3D encoding of virtual geometry, demands significantly more training time than the "Only Local" model. This indicates that the 3D encoder is indeed the primary computational bottleneck of our approach.

In conclusion, our proposed components, such as virtual geometry construction and local geometry learning, do not introduce additional computational cost. The main efficiency challenge arises from the 3D encoder, a component integrated from previous works.