Hybrid machine-learning/molecular-mechanics (ML/MM) methods extend the classical QM/MM paradigm by replacing the quantum description with neural network interatomic potentials trained to reproduce accurately quantum-mechanical (QM) results. By describing only the chemically active region with ML and the surrounding environment with molecular mechanics (MM), ML/MM models achieve near-QM/MM fidelity at a fraction of the computational cost, enabling routine simulation of reaction mechanisms, vibrational spectra, and binding free energies in complex biological or condensed-phase environments. The key challenge lies in coupling the ML and MM regions, a task addressed through three main strategies: (1) mechanical embedding (ME), where ML regions interact with fixed MM charges via classical electrostatics; (2) polarization-corrected mechanical embedding (PCME), where a vacuum-trained ML potential is supplemented post hoc with electrostatic corrections; and (3) environment-integrated embedding (EIE), where ML potentials are trained with explicit inclusion of MM-derived fields, enhancing accuracy but requiring specialized data. Since ML/MM builds on the scaffolding of QM/MM, most proposed coupling strategies rely heavily on electrostatics, polarization, and other physicochemical concepts, and the development and analysis of ML/MM schemes sits naturally at the intersection of physical chemistry and modern data science. This review surveys the conceptual foundations of ML/MM schemes, classifies existing implementations, and highlights key applications and open challenges, providing a critical snapshot of the current state-of-the-art and positioning ML/MM not merely as a computational alternative but as the natural evolution of QM/MM toward data-driven, scalable multiscale modeling.

We introduce an innovative machine learning (ML)-based framework for multiscale molecular modeling in which the ML subsystem is treated as an electrostatic entity interacting with its molecular mechanics (MM) environment through classical electrostatics. The integration of ML accuracy with multiscale modeling is accomplished by leveraging the capabilities of the ANI neural networks to predict geometry-dependent atomic partial charges at the minimal basis iterative stockholder (MBIS) level, going beyond static mechanical embedding. This ML/MM approach can closely approximate state-of-the-art multiscale quantum-classical (QM/MM) methods while significantly lowering computational requirements, thereby facilitating more efficient and precise simulations in computational chemistry. The method requires no additional training beyond the initial model setup and is integrated into Amber, one of the most widely used software suites for molecular modeling, ensuring accessibility to the broader community. We validate its performance across a variety of challenging applications, including the solvation structure, vibrational spectra, torsion free energy profiles, and protein-ligand interactions, achieving excellent agreement with QM/MM benchmarks. This framework not only advances the frontiers of multiscale modeling but also showcases the potential of machine learning to achieve quantum-level accuracy with exceptional efficiency for complex chemical systems.

In this work we introduce TorchANI-Amber, an interface for routine molecular dy- namics simulations of biomolecular systems using ANI-style machine learning poten- tials. TochANI-Amber incorporates the ANI neural network potentials into the Amber software suite, supporting Amber’s two engines: sander and pmemd. In addition to implementing all published ANI models, the interface is extensible to other energy predicting potentials through a simple mechanism requiring no knowledge of Amber’s codebase. To illustrate this versatility, we implement extensions to the AIMNet2 and Nutmeg potentials. The interface is integrated with Amber’s neighborlists, and it also supports an optimized CUDA implementation for computing the ANI models’ feature vectors, enabling simulations of systems with hundreds of thousands of atoms at the neural network’s level of theory (approaching DFT accuracy). The interface is designed so that all amber capabilities can be used with ANI potentials instead of force fields. To evaluate the energy conservation, stability, and performance of the ANI potentials as used through the interface, we run MD simulations on different biomolecular systems, including ubiquitin and Trp-cage proteins in explicit solvent. Additionally, we demon- strate the use of these potentials in the context of enhanced sampling techniques, such as temperature replica-exchange molecular dynamics.

Coenzyme A (CoA) is a key cellular metabolite which participates in diverse metabolic pathways, regulation of gene expression and the antioxidant defense mechanism. Human NME1 (hNME1), which is a moonlighting protein, was identified as a major CoA-binding protein. Biochemical studies showed that hNME1 is regulated by CoA through both covalent and non-covalent binding, which leads to a decrease in the hNME1 nucleoside diphosphate kinase (NDPK) activity. In this study, we expanded the knowledge on previous findings by focusing on the non-covalent mode of CoA binding to the hNME1. With X-ray crystallography, we solved the CoA bound structure of hNME1 (hNME1-CoA) and determined the stabilization interactions CoA forms within the nucleotide-binding site of hNME1. A hydrophobic patch stabilizing the CoA adenine ring, while salt bridges and hydrogen bonds stabilizing the phosphate groups of CoA were observed. With molecular dynamics studies, we extended our structural analysis by characterizing the hNME1-CoA structure and elucidating possible orientations of the pantetheine tail, which is absent in the X-ray structure due to its flexibility. Crystallographic studies suggested the involvement of arginine 58 and threonine 94 in mediating specific interactions with CoA. Site-directed mutagenesis and CoA-based affinity purifications showed that arginine 58 mutation to glutamate (R58E) and threonine 94 mutation to aspartate (T94D) prevent hNME1 from binding to CoA. Overall, our results reveal a unique mode by which hNME1 binds CoA, which differs significantly from that of ADP binding: the α- and β-phosphates of CoA are oriented away from the nucleotide-binding site, while 3′-phosphate faces catalytic histidine 118 (H118). The interactions formed by the CoA adenine ring and phosphate groups contribute to the specific mode of CoA binding to hNME1.

We present a novel integration of the ANI neural networks into the Amber software suite, offering a sophisticated machine learning/molecular mechanics (ML/MM) framework. The implementation is designed as a general-purpose tool for the simulation of neutral organic molecules, requiring no additional training for its use beyond the initial setup. The framework leverages a new ANI potential that accurately predicts geometry-dependent atomic partial charges at the Minimal Basis Iterative Stockholder (MBIS) level, enhancing the modeling of electrostatic interactions within ML/MM systems. Additionally, we incorporate a polarization correction to address the distortion effects on the ML subsystem from MM point charges. Our approach is validated through simulations of solvation profiles, vibrational spectra, and torsion free energy profiles of small molecules in aqueous environments, as well as protein-ligand interactions. Our findings demonstrate that this ML/MM framework can approximate QM/MM electrostatic embedding with significantly reduced computational demands, paving the way for more efficient and accurate simulations in computational chemistry.

Machine Learning (ML) methods have reached high accuracy levels for the prediction of in vacuo molecular properties. However, the simulation of large systems through solely ML methods (like those based on neural network potentials) is still a challenge. In this context, one of the most promising frameworks for integrating ML schemes in the simulation of complex molecular systems are the so-called ML/MM methods. These multiscale approaches combine ML methods with classical forcefields (MM), in the same spirit as the succesful hybrid quantum mechanics-molecular mechanics methods (QM/MM). The key issue for such ML/MM methods is the adequate description of the coupling between the region of the system described by ML and the region described at the MM level. In the context of QM/MM schemes, the main ingredient of the interaction is electrostatic, and the state of the art is the so called electrostatic-embedding. In this study, we analyze the quality of simpler mechanical embedding-based approaches, specifically focusing on their application within a ML/MM framework utilizing atomic partial charges derived in vacuo. Taking as reference electrostatic embedding calculations performed at a QM(DFT)/MM level, we explore different atomic charges schemes, as well as a polarization correction computed using atomic polarizabilites. Our benchmark data set comprises a set of about 80k small organic structures from the ANI-1x database, solvated in water. The results suggest that the MBIS atomic charges yield the best agreement with the reference coupling energy. Remarkable enhancements are achieved by including a simple polarization correction.

Magnesium (Mg2+), the second most abundant intracellular cation, plays a crucial role in cellular functions. In this study, we investigate the interaction between Mg2+ and coenzyme A (CoA), a thiol-containing cofactor central to cellular metabolism also involved in protein modifications. Isothermal titration calorimetry revealed a 1:1 binding stoichiometry between Mg2+ and free CoA under biologically relevant conditions. Association constants of (537 ± 20) M-1 and (312 ± 7) M-1 were determined at 25 °C and pH 7.2 and 7.8, respectively, suggesting that a significant fraction of CoA is likely bound to Mg2+ both in the cytosol and in the mitochondrial matrix. Additionally, the process is entropically-driven, and our results support that the origin of the entropy gain is solvent-related. On the other hand, the combination of 1- and 2-dimensional nuclear magnetic resonance spectroscopy with molecular dynamics simulations and unsupervised learning demonstrate a direct coordination between Mg2+ and the phosphate groups of the 4-phosphopantothenate unit and bound to position 5' of the adenosine ring. Interestingly, the phosphate in position 3' only indirectly contributes to Mg2+ coordination. Finally, we discuss how the binding of Mg2+ to CoA perturbates the chemical environment of different CoA atoms, regardless of their apparent proximity to the coordination site, through the modulation of the CoA conformational landscape. This insight holds implications for understanding the impact on both CoA and Mg2+ functions in physiological and pathological processes.

We present a comprehensive theoretical examination of the structural properties of dianionic polysulfides [S-n](2-) (n = 2-6), their conjugated monoacids [HSn](-) (n = 2-6), and a selection of 1e(-)-oxidized radical anions [S-n](center dot-) (n = 2-4), in aqueous and dimethyl sulfoxide (DMSO) solutions. We investigated the structures and stabilities of various conformational isomers within these families of compounds by employing Quantum Mechanics-Molecular Mechanics (QM-MM) Molecular Dynamics (MD) simulations. The explicit inclusion of solvent molecules in the calculations revealed stable conformational structures that were previously unreported and might have appreciable concentrations in real systems. The interconversions between the isomeric structures proceed on the order of hundreds of picoseconds and are energetically similar to the isomerization processes in substituted cyclohexanes. We also conducted a detailed analysis of the stability of different isomers of the radical anion [S-4](center dot-) in solution. Our findings highlight the significant influence of the solvent on the isomerizations, a result that could be particularly relevant for enhancing the performance of metal-sulfur batteries.

The oxidation of Met to methionine sulfoxide (MetSO) by oxidants such as hydrogen peroxide, hypochlorite, or peroxynitrite has profound effects on protein function. This modification can be reversed by methionine sulfoxide reductases (msr). In the context of pathogen infection, the reduction of oxidized proteins gains significance due to microbial oxidative damage generated by the immune system. For example, Mycobacterium tuberculosis (Mt) utilizes msrs (MtmsrA and MtmsrB) as part of the repair response to the host-induced oxidative stress. The absence of these enzymes makes Mycobacteria prone to increased susceptibility to cell death, pointing them out as potential therapeutic targets. This study provides a detailed characterization of the catalytic mechanism of MtmsrA using a comprehensive approach, including experimental techniques and theoretical methodologies. Confirming a ping-pong type enzymatic mechanism, we elucidate the catalytic parameters for sulfoxide and thioredoxin substrates (k(cat)/K-M = 2656 +/- 525 M-1 s(-1) and 1.7 +/- 0.8 x 10(6) M-1 s(-1), respectively). Notably, the entropic nature of the activation process thermodynamics, representing similar to 85% of the activation free energy at room temperature, is underscored. Furthermore, the current study questions the plausibility of a sulfurane intermediate, which may be a transition-state-like structure, suggesting the involvement of a conserved histidine residue as an acid-base catalyst in the MetSO reduction mechanism. This mechanistic insight not only advances our understanding of Mt antioxidant enzymes but also holds implications for future drug discovery and biotechnological applications.

Peroxynitrite is a very reactive species implicated in a variety of pathophysiological cellular processes. Particularly, peroxynitrite-mediated oxidation of cellular thiol-containing compounds such as cysteine residues is a key process which has been extensively studied. Cysteine plays roles in many redox biochemistry pathways. In contrast, selenocysteine, the 21st amino acid, is only present in 25 human proteins. Investigating the molecular basis of selenocysteine's reactivity may provide insights into its unique role in these selenocysteine-containing proteins. The two-electron oxidation of thiols or selenols by peroxynitrite is a process that is carried out by the thiolate/selenate forms and peroxynitrous acid.In this work, we shed light on the molecular basis of the differential reactivity of both species towards peroxynitrite by means of state-of-the-art computer simulations. We performed electronic structure calculations of the reaction in the methanethiolate and methaneselenolate model systems with peroxynitrous acid at different levels of theory using an implicit solvent scheme. In addition, we employed a multi-scale quantum mechanics/molecular mechanics approach for obtaining free energy profiles of these chemical reactions in aqueous solution, which enabled the comparison between the simulations and the available experimental data. Our results suggest that the larger reactivity observed in the selenocysteine case at physiological pH is mainly due to the lower pKa, which affords a larger fraction of the reactive anionic species in these conditions, and in a second place to a slightly enhanced intrinsic reactivity of the selenate form due to its larger nucleophilicity.

The determination of minimum free energy pathways (MFEP) is one of the most widely used strategies to study reactive processes. For chemical reactions in complex environments, the combination of quantum mechanics (QM) with a molecular mechanics (MM) representation is usually necessary in a hybrid QM/MM framework. However, even within the QM/MM approximation, the affordable sampling of the phase space is, in general, quite restricted. To reduce drastically the computational cost of the simulations, several methods such as umbrella sampling require performing a priori a selection of a reaction coordinate. The quality of the computed results, in an affordable computational time, is intimately related to the reaction coordinate election which is, in general, a nontrivial task. In this work, we provide an approach to model reactive processes in complex environments that does not require the a priori selection of a reaction coordinate. The proposed methodology combines QM/MM simulations with an extrapolation of the nudged elastic bands (NEB) method to the free energy surface (FENEB). We present and apply our own FENEB scheme to optimize MFEP in different reactive processes, using QM/MM frameworks at semiempirical and density functional theory levels. Our implementation is based on performing the FENEB optimization by uncoupling the optimization of the band in a perpendicular and tangential direction. In each step, a full optimization with the spring force is performed, which guarantees that the images remain evenly distributed. The robustness of the method and the influence of sampling on the quality of the optimized MFEP and its associated free energy barrier are studied. We show that the FENEB method provides a good estimation of the reaction barrier even with relatively short simulation times, supporting that its combination with QM/MM frameworks provides an adequate tool to study chemical processes in complex environments.

Challenging the basis of our chemical intuition, recent experimental evidence reveals the presence of a new type of intrinsic fluorescence in biomolecules that exists even in the absence of aromatic or electronically conjugated chemical compounds. The origin of this phenomenon has remained elusive so far. In the present study, we identify a mechanism underlying this new type of fluorescence in different biological aggregates. By employing non-adiabatic ab initio molecular dynamics simulations combined with a data-driven approach, we characterize the typical ultrafast non-radiative relaxation pathways active in non-fluorescent peptides. We show that the key vibrational mode for the non-radiative decay towards the ground state is the carbonyl elongation. Non-aromatic fluorescence appears to emerge from blocking this mode with strong local interactions such as hydrogen bonds. While we cannot rule out the existence of alternative non-aromatic fluorescence mechanisms in other systems, we demonstrate that this carbonyl-lock mechanism for trapping the excited state leads to the fluorescence yield increase observed experimentally, and set the stage for design principles to realize novel non-invasive biocompatible probes with applications in bioimaging, sensing, and biophotonics.

It is well established that proteins and peptides can release sulfur under alkaline treatment, mainly through the beta-elimination of disulfides with the concomitant formation of persulfides and dehydroalanine derivatives. In this study, we evaluated the formation of glutathione persulfide (GSSH/GSS-) by exposure of glutathione disulfide (GSSG) to alkaline conditions. The kinetics of the reaction between GSSG and HO- was investigated by UV-Vis absorbance, reaction with 5,5'-dithio-bis-(2-nitrobenzoic acid) (DTNB), and cold cyanolysis, obtaining an apparent second-order rate constant of similar to 10(-3) M-1 s(-1) at 25 degrees C. The formation of GSSH and the dehydroalanine derivative was confirmed by HPLC and/or mass spectrometry. However, the mixtures did not equilibrate in a timescale of hours, and additional species, including thiol and diverse sulfane sulfur compounds were also formed, probably through further reactions of the persulfide. Cold cyanolysis is frequently used to quantify persulfides, since it measures sulfane sulfur. This method involves a step in which the sample to be analyzed is incubated with cyanide at alkaline pH. When cold cyanolysis was applied to samples containing GSSG, sulfane sulfur products that were not present in the original sample were measured. Thus, our results reveal the risk of overestimating the amount of sulfane sulfur compounds in samples that contain disulfides due to their decay to persulfides and other sulfane sulfur compounds at alkaline pH. Overall, our study highlights that the beta-elimination of disulfides is a potential source of persulfides, although we do not recommend the preparation of GSSH from incubation of GSSG in alkali. Our study also highlights the importance of being cautious when doing and interpreting cold cyanolysis experiments.

The mechanism of the metal centered reduction of metmyoglobin (MbFeIII) by sulfide species (H2S/HS-) under an argon atmosphere has been studied by a combination of spectroscopic, kinetic, and computational methods. Asymmetric S-shaped time-traces for the formation of MbFeII at varying ratios of excess sulfide were observed at pH 5.3 < pH < 8.0 and 25 °C, suggesting an autocatalytic reaction mechanism. An increased rate at more alkaline pHs points to HS- as relevant reactive species for the reduction. The formation of the sulfanyl radical (HS•) in the slow initial phase was assessed using the spin-trap phenyl N-tert-butyl nitrone. This radical initiates the formation of S-S reactive species as disulfanuidyl/ disulfanudi-idyl radical anions and disulfide (HSSH•-/HSS•2- and HSS-, respectively). The autocatalysis has been ascribed to HSS-, formed after HSSH•-/HSS•2- disproportionation, which behaves as a fast reductant toward the intermediate complex MbFeIII(HS-). We propose a reaction mechanism for the sulfide-mediated reduction of metmyoglobin where only ferric heme iron initiates the oxidation of sulfide species. Beside the chemical interest, this insight into the MbFeIII/sulfide reaction under an argon atmosphere is relevant for the interpretation of biochemical aspects of ectopic myoglobins found on hypoxic tissues toward reactive sulfur species.

Persulfides (RSSH/RSS-) are species closely related to thiols (RSH/RS-) and hydrogen sulfide (H2S/HS-), and can be formed in biological systems in both low and high molecular weight cysteine-containing compounds. They are key intermediates in catabolic and biosynthetic processes, and have been proposed to participate in the transduction of hydrogen sulfide effects. Persulfides are acidic, more acidic than thiols, and the persulfide anions are expected to be the predominant species at neutral pH. The persulfide anion has high nucleophilicity, due in part to the alpha effect, i.e., the increased reactivity of a nucleophile when the neighboring atom has high electron density. In addition, persulfides have electrophilic character, a property that is absent in both thiols and hydrogen sulfide. In this article, the biochemistry of persulfides is described, and the possible ways in which the formation of a persulfide could impact on the properties of the biomolecule involved are discussed.

Coenzyme A (CoA) is a key cellular metabolite known for its diverse functions in metabolism and regulation of gene expression. CoA was recently shown to play an important antioxidant role under various cellular stress conditions by forming a disulfide bond with proteins, termed CoAlation. Using anti-CoA antibodies and liquid chromatography tandem mass spectrometry (LC-MS/MS) methodologies, CoAlated proteins were identified from various organisms/tissues/cell-lines under stress conditions. In this study, we integrated currently known CoAlated proteins into mammalian and bacterial datasets (CoAlomes), resulting in a total of 2093 CoAlated proteins (2862 CoAlation sites). Functional classification of these proteins showed that CoAlation is widespread among proteins involved in cellular metabolism, stress response and protein synthesis. Using 35 published CoAlated protein structures, we studied the stabilization interactions of each CoA segment (adenosine diphosphate (ADP) moiety and pantetheine tail) within the microenvironment of the modified cysteines. Alternating polar-non-polar residues, positively charged residues and hydrophobic interactions mainly stabilize the pantetheine tail, phosphate groups and the ADP moiety, respectively. A flexible nature of CoA is observed in examined structures, allowing it to adapt its conformation through interactions with residues surrounding the CoAlation site. Based on these findings, we propose three modes of CoA binding to proteins. Overall, this study summarizes currently available knowledge on CoAlated proteins, their functional distribution and CoA–protein stabilization interactions.

The role of inorganic sulfur species in biological systems has gained considerable interest since the recognition of sulfanes, particularly dihydrogen sulfide or sulfane, H 2 S, disulfane, HSSH, trisulfane, HSSSH, and their conjugate bases, as endogenous species and mediators of signaling functions in different tissues. The one-electron oxidation of H 2 S/HS − has been assigned as the onset of signaling processes or oxidative detoxification mechanisms. These varied sulfur containing inorganic species are, together with organic counterparts, reunited as reactive sulfur species (RSS). In order to shed light on this rich and still not completely explored chemistry, we have performed electronic structure calculations at different levels of theory, to provide estimations and the molecular basis of the pK a values of the polysulfides HSSH and HSSSH and of the radical HS • . In addition, we also reported the characterization of selected inorganic RSS including both radical and non-radical species with different protonation states with the intention of assisting the interpretation of chemical/biochemical experiments involving these species.

The optimization of minimum free energy pathways (MFEP) is one of the most widely used strategies to study activated processes. For chemical reactions, this requires the use of quantum mechanics. Using quantum mechanics molecular mechanics (QM- MM) Hamiltionians allows the simulation of reactive processes in complex environments by treating with quantum mechanics only the chemically relevant part of the system. However, even within this approximation, the affordable simulation lenghts of QM-MM simulations is in general, quite limited. Free energy methods based on the sampling of the potential energy surface require long simulations times to provide converged and accurate results. As consequence, the combination of QM-MM methods and free energy calculations is computationally expensive. Moreover, the user usually needs to perform an a priori selection of the reaction coordinate. This may be not trivial for the general case. One of the most established methods for finding potential energy profiles without selecting a reaction coordinate is the nudged elastic band method (NEB). In this work, we used the extension of this method to the exploration of the free energy surface for finding MFEP (FENEB). We present and apply to reactive systems an improved version of the basic optimization scheme of FENEB that increases its robustness, and is based on decoupling the optimization of the band in the perpendicular direction to the band, from the optimization of the tangential direction. In each optimization step, a full optimization with the spring force is performed, in order to keep the images evenly distributed. Additionally, we evaluate the influence of sampling in the quality of the optimized MFEP and free energy barrier computed from it. We show and discuss that the FENEB method provides a good estimation of the reaction barrier even with relatively short simulations lenghts and that it scales better than umbrella sampling both with simulation lenght and with dimensionality. Overall, our results support that the combination of QM-MM methods and the FENEB provides an adequate tool study chemical processes in complex environments.