Modeling Chemical Reactions by QM/MM Calculations: The Case of the Tautomerization in Fireflies Bioluminescent Systems

PyDFT-QMMM: A modular, extensible software framework for DFT-based QM/MM molecular dynamics

PyDFT-QMMM is a Python-based package for performing hybrid quantum mechanics/molecular mechanics (QM/MM) simulations at the density functional level of theory. The program is designed to treat short-range and long-range interactions through user-specified combinations of electrostatic and mechanical embedding procedures within periodic simulation domains, providing necessary interfaces to external quantum chemistry and molecular dynamics software. To enable direct embedding of long-range electrostatics in periodic systems, we have derived and implemented force terms for our previously described QM/MM/PME approach [Pederson and McDaniel, J. Chem. Phys. 156, 174105 (2022)]. Communication with external software packages Psi4 and OpenMM is facilitated through Python application programming interfaces (APIs). The core library contains basic utilities for running QM/MM molecular dynamics simulations, and plug-in entry-points are provided for users to implement custom energy/force calculation and integration routines, within an extensible architecture. The user interacts with PyDFT-QMMM primarily through its Python API, allowing for complex workflow development with Python scripting, for example, interfacing with PLUMED for free energy simulations. We provide benchmarks of forces and energy conservation for the QM/MM/PME and alternative QM/MM electrostatic embedding approaches. We further demonstrate a simple example use case for water solute in a water solvent system, for which radial distribution functions are computed from 100 ps QM/MM simulations; in this example, we highlight how the solvation structure is sensitive to different basis-set choices due to under- or over-polarization of the QM water molecule's electron density.

The C-terminal Region of D-DT Regulates Molecular Recognition for Protein-Ligand Complexes

Systematic analysis of molecular recognition is critical for understanding the biological function of macromolecules. For the immunomodulatory protein D-dopachrome tautomerase (D-DT), the mechanism of protein-ligand interactions is poorly understood. Here, 17 carefully designed protein variants and wild type (WT) D-DT were interrogated with an array of complementary techniques to elucidate the structural basis of ligand recognition. Utilization of a substrate and two selective inhibitors with distinct binding profiles offered previously unseen mechanistic insights into D-DT-ligand interactions. Our results demonstrate that the C-terminal region serves a key role in molecular recognition via regulation of the active site opening, protein-ligand interactions, and conformational flexibility of the pocket's environment. While our study is the first comprehensive analysis of molecular recognition for D-DT, the findings reported herein promote the understanding of protein functionality and enable the design of new structure-based drug discovery projects.

Unraveling the contributions to the spectral shape of flexible dyes in solution: insights on the absorption spectrum of an oxyluciferin analogue

We present a computational investigation of the absorption spectrum in water of 5,5-spirocyclopropyl-oxyluciferin (5,5-CprOxyLH), an analogue of the emitter compound responsible for the bioluminescence in fireflies. Several factors participate in determining the 5,5-CprOxyLH's spectral shape: (i) the contribution of the four close-energy excited states, which show significant non-adiabatic couplings, (ii) the flexible molecular structure and (iii) the specific interactions established with the surrounding environment, which strongly couple the protic solvent dynamics with the dye's spectral response. To tackle the challenge to capture and dissect the role of all these effects we preliminarily investigate the role of non-adiabatic couplings with quantum dynamics simulations and a linear vibronic coupling model in the gas phase. Then, we account for both the molecular flexibility and solvent interactions by resorting to a mixed quantum classical protocol, named Adiabatic Molecular Dynamics generalized Vertical Gradient (Ad-MD|gVG), which is built on a method recently proposed by some of us. It is rooted in the partition between stiff degrees of freedom of the dye, accounted for at the vibronic level within the harmonic approximation, and flexible degrees of freedom of the solute (and of the solvent), described classically through a sampling based on Molecular Dynamics (MD). Ad-MD|gVG avoids spurious effects arising in the excited state Hessians due to non-adiabatic couplings, and can therefore be applied to account for the contributions of the first four excited states to the 5,5-CprOxyLH absorption spectrum. The final simulated spectrum is in very good agreement with the experiment, especially when the MD is driven by a refined quantum-mechanically derived force-field. More importantly, the origin of each separate contribution to the spectral shape is appropriately accounted for, paving the way to future applications of the method to more complex systems or alternative spectroscopies, as emission or circular dichroism.

Computational Investigation of Substituent Effects on the Fluorescence Wavelengths of Oxyluciferin Analogs

Oxyluciferin, which is the light emitter for firefly bioluminescence, has been subjected to extensive chemical modifications to tune its emission wavelength and quantum yield. However, the exact mechanisms for various electron-donating and withdrawing groups to perturb the photophysical properties of oxyluciferin analogs are still not fully understood. To elucidate the substituent effects on the fluorescence wavelength of oxyluciferin analogs, we applied the absolutely localized molecular orbitals (ALMO)-based frontier orbital analysis to assess various types of interactions (i.e. permanent electrostatics/exchange repulsion, polarization, occupied-occupied orbital mixing, virtual-virtual orbital mixing, and charge-transfer) between the oxyluciferin and substituent orbitals. We suggested two distinct mechanisms that can lead to red-shifted oxyluciferin emission wavelength, a design objective that can help increase the tissue penetration of bioluminescence emission. Within the first mechanism, an electron-donating group (such as an amino or dimethylamino group) can contribute its highest occupied molecular orbital (HOMO) to an out-of-phase combination with oxyluciferin's HOMO, thus raising the HOMO energy of the substituted analog and narrowing its HOMO-LUMO gap. Alternatively, an electron-withdrawing group (such as a nitro or cyano group) can participate in an in-phase virtual-virtual orbital mixing of fragment LUMOs, thus lowering the LUMO energy of the substituted analog. Such an ALMO-based frontier orbital analysis is expected to lead to intuitive principles for designing analogs of not only the oxyluciferin molecule, but also many other functional dyes.

Spectrochemistry of Firefly Bioluminescence

The chemical reactions underlying the emission of light in fireflies and other bioluminescent beetles are some of the most thoroughly studied processes by scientists worldwide. Despite these remarkable efforts, fierce academic arguments continue around even some of the most fundamental aspects of the reaction mechanism behind the beetle bioluminescence. In an attempt to reach a consensus, we made an exhaustive search of the available literature and compiled the key discoveries on the fluorescence and chemiluminescence spectrochemistry of the emitting molecule, the firefly oxyluciferin, and its chemical analogues reported over the past 50+ years. The factors that affect the light emission, including intermolecular interactions, solvent polarity, and electronic effects, were analyzed in the context of both the reaction mechanism and the different colors of light emitted by different luciferases. The collective data points toward a combined emission of multiple coexistent forms of oxyluciferin as the most probable explanation for the variation in color of the emitted light. We also highlight realistic research directions to eventually address some of the remaining questions related to firefly bioluminescence. It is our hope that this extensive compilation of data and detailed analysis will not only consolidate the existing body of knowledge on this important phenomenon but will also aid in reaching a wider consensus on some of the mechanistic details of firefly bioluminescence.

Reaction mechanisms of C(3PJ) and C+(2PJ) with benzene in the interstellar medium from quantum mechanical molecular dynamics simulations

While spectroscopic data on small hydrocarbons in interstellar media in combination with crossed molecular beam (CMB) experiments have provided a wealth of information on astrochemically relevant species, much of the underlying mechanistic pathways of their formation remain elusive. Therefore, in this work, the chemical reaction mechanisms of C(3PJ) + C6H6 and C+(2P) + C6H6 systems using the quantum mechanical molecular dynamics (QMMD) technique at the PBE0-D3(BJ) level of theory is investigated, mimicking a CMB experiment. Both the dynamics of the reactions as well as the electronic structure for the purpose of the reaction network are evaluated. The method is validated for the first reaction by comparison to the available experimental data. The reaction scheme for the C(3PJ) + C6H6 system covers the literature data, e.g. the major products are the 1,2-didehydrocycloheptatrienyl radical (C7H5) and benzocyclopropenyl radical (C6H5-CH), and it reveals the existence of less common pathways for the first time. The chemistry of the C+(2PJ) + C6H6 system is found to be much richer, and we have found that this is because of more exothermic reactions in this system in comparison to those in the C(3PJ) + C6H6 system. Moreover, using the QMMD simulation, a number of reaction paths have been revealed that produce three distinct classes of reaction products with different ring sizes. All in all, at all the collision energies and orientations, the major product is the heptagon molecular ion for the ionic system. It is also revealed that the collision orientation has a dominant effect on the reaction products in both systems, while the collision energy mostly affects the charged system. These simulations both prove the applicability of this approach to simulate crossed molecular beams, and provide fundamental information on reactions relevant for the interstellar medium.

The role of CO2 detachment in fungal bioluminescence: thermally vs. excited state induced pathways

Different fungi lineages are known to emit light on Earth, mainly in tropical climates. Although the preparation of bioluminescent cell-free extracts allowed one to characterize the enzymatic requirements, the molecular mechanism underlying luminescence is still largely unknown and is based on the experimental putative assumption that a high-energy intermediate should be formed by reaction with O2 and formation of an endoperoxide. Here, we aim at determining, through state-of-the-art multiconfigurational quantum chemistry, the full mechanistic landscape leading from the endoperoxide to the emitting species, envisaging different possible pathways and proposing their viability. Especially, thermal CO2 detachment followed by excited-state peroxide opening (thermal-chemiluminescence) can compete with a parallel pathway, i.e., first excited-state endoperoxide opening, followed by CO2 detachment on the same excited-state (excited state-chemiluminescence). Clear differences in the energy supplies, as well as the possibility to directly populate the emitting species from the intersection seam between ground and excited states, land credence to a kinetically efficient thermal-chemiluminescent pathway, establishing for the first time a detailed description of fungal bioluminescence.

Luciferin Regeneration in Firefly Bioluminescence via Proton-Transfer-Facilitated Hydrolysis, Condensation and Chiral Inversion

Firefly bioluminescence is produced via luciferin enzymatic reactions in luciferase. Luciferin has to be unceasingly replenished to maintain bioluminescence. How is the luciferin reproduced after it has been exhausted? In the early 1970s, Okada proposed the hypothesis that the oxyluciferin produced by the previous bioluminescent reaction could be converted into new luciferin for the next bioluminescent reaction. To some extent, this hypothesis was evidenced by several detected intermediates. However, the detailed process and mechanism of luciferin regeneration remained largely unknown. For the first time, we investigated the entire process of luciferin regeneration in firefly bioluminescence by density functional theory calculations. This theoretical study suggests that luciferin regeneration consists of three sequential steps: the oxyluciferin produced from the last bioluminescent reaction generates 2-cyano-6-hydroxybenzothiazole (CHBT) in the luciferin regenerating enzyme (LRE) via a hydrolysis reaction; CHBT combines with L-cysteine in vivo to form L-luciferin via a condensation reaction; and L-luciferin inverts into D-luciferin in luciferase and thioesterase. The presently proposed mechanism not only supports the sporadic evidence from previous experiments but also clearly describes the complete process of luciferin regeneration. This work is of great significance for understanding the long-term flashing of fireflies without an in vitro energy supply.

Effect of Protein Conformation and AMP Protonation State on Fireflies' Bioluminescent Emission

The emitted color in fireflies’ bioluminescent systems depends on the beetle species the system is extracted from and on different external factors (pH, temperature…) among others. Controlling the energy of the emitted light (i.e., color) is of crucial interest for the use of such bioluminescent systems. For instance, in the biomedical field, red emitted light is desirable because of its larger tissue penetration and lower energies. In order to investigate the influence of the protein environment and the AMP protonation state on the emitted color, the emission spectra of the phenolate-keto and phenolate-enol oxyluciferin forms have been simulated by means of MD simulations and QM/MM calculations, considering: two different protein conformations (with an open or closed C-terminal domain with respect to the N-terminal) and two protonation states of AMP. The results show that the emission spectra when considering the protein characterized by a closed conformation are blue-shifted compared to the open conformation. Moreover, the complete deprotonation of AMP phosphate group (AMP²⁻) can also lead to a blue-shift of the emission spectra but only when considering the closed protein conformation (open form is not sensitive to changes of AMP protonation state). These findings can be reasoned by the different interactions (hydrogen-bonds) found between oxyluciferin and the surrounding (protein, AMP and water molecules). This study gets partial insight into the possible origin of the emitted color modulation by changes of the pH or luciferase conformations.

Amide-imide tautomerization in the glutamine side chain in enzymatic and photochemical reactions in proteins

Amide-imide tautomerization presents a pervasive class of chemical transformations in organic chemistry of natural compounds. In this Perspective, we describe two distinctively different protein systems, in which the amide-imide tautomerization in the glutamine side chain takes place in enzymatic or photochemical reactions. First, hydrolysis of guanosine triphosphate (GTP) catalyzed by the Ras-GAP protein complex suggests the occurrence of the imide tautomer of glutamine in reaction intermediates. Second, photoexcitation of flavin-binding protein domains (BLUFs) initiates a chain of reactions in the chromophore-binding pocket, including amide-imide tautomerization of glutamine. Mechanisms of these reactions at the atomic level have been revealed in quantum mechanics/molecular mechanics (QM/MM) simulations. To reinforce conclusions on the critical role of amide-imide tautomerization of glutamine in these reactions we describe results of new quantum chemistry and QM/MM calculations for relevant molecular model systems. We reexamine results of the recent IR spectroscopy studies of BLUF domains, which provide experimental evidences of Gln tautomerization in proteins. We also propose to validate the glutamine-assisted mechanism of enzymatic GTP hydrolysis by using IR spectroscopy in a proper range of wavenumbers.