CHARMM Force-Field Parameters for Morphine, Heroin, and Oliceridine, and Conformational Dynamics of Opioid Drugs.

Opioid drug binding to specialized G protein-coupled receptors (GPCRs) can lead to analgesia upon activation via downstream Gi protein signaling and to severe side effects via activation of the ß-arrestin signaling pathway. Knowledge of how different opioid drugs interact with receptors is essential, as it can inform and guide the design of safer therapeutics. We performed quantum and classical mechanical computations to explore the potential energy landscape of four opioid drugs: morphine and its derivatives heroin and fentanyl and for the unrelated oliceridine. From potential energy profiles for bond twists and from interactions between opioids and water, we derived a set of force-field parameters that allow a good description of structural properties and intermolecular interactions of the opioids. Potential of mean force profiles computed from molecular dynamics simulations indicate that fentanyl and oliceridine have complex energy landscapes with relatively small energy penalties, suggesting that interactions with the receptor could select different binding poses of the drugs.

Potential Energy Function for Fentanyl-Based Opioid Pain Killers.

Opioids are molecules whose binding to specialized G-Protein Coupled Receptors (GPCRs) triggers a signaling cascade that leads to the downregulation of pain pathways. Binding of an opioid to the membrane-embedded GPCR occurs when the opioid molecule is protonated, which provides a potential strategy to design nontoxic opioids that are protonated and bind to the GPCR only at the low pH of injured or inflamed tissue. Excellent model systems to study protonation-dependent binding of opioids to GPCRs are fentanyl, which is protonated and binds to the GPCR at both physiological and low pH, and the fluorinated fentanyl derivative NFEPP, which is protonated and binds to the GPCR only at low pH. The molecular mechanisms of fentanyl and NFEPP binding to the GPCR are largely unknown. To enable atomistic studies of opioid binding to GPCRs, we have carried out extensive quantum mechanical and classical mechanical computations to derive a potential energy function for fentanyl and NFEPP and present force field parameters for both opioid molecules. We find that fluorination alters the electronic ground state properties of fentanyl. As a consequence, fentanyl and NFEPP have distinct torsional and electrostatic properties likely to impact how they bind to receptors.

Loop Grafting between Similar Local Environments for Fc-Silent Antibodies.

Reduction of the affinity of the fragment crystallizable (Fc) region with immune receptors by substitution of one or a few amino acids, known as Fc-silencing, is an established approach to reduce the immune effector functions of monoclonal antibody therapeutics. This approach to Fc-silencing, however, is problematic as it can lead to instability and immunogenicity of the developed antibodies. We evaluated loop grafting as a novel approach to Fc-silencing in which the Fc loops responsible for immune receptor binding were replaced by loops of up to 20 amino acids from similar local environments in other human and mouse antibodies. Molecular dynamics simulations of the designed variants of an Fc region in a complex with the immune receptor FcγIIIa confirmed that loop grafting potentially leads to a significant reduction in the binding of the antibody variants to the receptor, while retaining their stability. In comparison, standard variants with less than eight substituted amino acids showed possible instability and a lower degree of Fc-silencing due to the occurrence of compensatory interactions. The presented approach to Fc-silencing is general and could be used to modulate undesirable side effects of other antibody therapeutics without affecting their stability or increasing their immunogenicity.

Simulations of NPC1(NTD):NPC2 Protein Complex Reveal Cholesterol Transfer Pathways.

The Niemann Pick type C (NPC) proteins, NPC1 and NPC2, are involved in the lysosomal storage disease, NPC disease. The formation of a NPC1⁻NPC2 protein⁻protein complex is believed to be necessary for the transfer of cholesterol and lipids out of the late endosomal (LE)/lysosomal (Lys) compartments. Mutations in either NPC1 or NPC2 can lead to an accumulation of cholesterol and lipids in the LE/Lys, the primary phenotype of the NPC disease. We investigated the NPC1(NTD)⁻NPC2 protein⁻protein complex computationally using two putative binding interfaces. A combination of molecular modeling and molecular dynamics simulations reveals atomic details that are responsible for interface stability. Cholesterol binding energies associated with each of the binding pockets for the two models are calculated. Analyses of the cholesterol binding in the two models support bidirectional ligand transfer when a particular interface is established. Based on the results, we propose that, depending on the location of the cholesterol ligand, a dynamical interface between the NPC2 and NPC1(NTD) proteins exists. Structural features of a particular interface can lower the energy barrier and stabilize the passage of the cholesterol substrate from NPC2 to NPC1(NTD).

Insight into Inhibitor Binding in the Eukaryotic Proteasome: Computations of the 20S CP.

A combination of molecular dynamics (MD) simulations and computational analyses uncovers structural features that may influence substrate passage and exposure to the active sites within the proteolytic chamber of the 20S proteasome core particle (CP). MD simulations of the CP reveal relaxation dynamics in which the CP slowly contracts over the 54 ns sampling period. MD simulations of the SyringolinA (SylA) inhibitor within the proteolytic B 1 ring chamber of the CP indicate that favorable van der Waals and electrostatic interactions account for the predominant association of the inhibitor with the walls of the proteolytic chamber. The time scale required for the inhibitor to travel from the center of the proteolytic chamber to the chamber wall is on the order of 4 ns, accompanied by an average energetic stabilization of approximately -20 kcal/mol.

Crystal Structures of Diaryliodonium Fluorides and Their Implications for Fluorination Mechanisms.

The radiofluorination of diaryliodonium salts is of value for producing radiotracers for positron emission tomography. We report crystal structures for two diaryliodonium fluorides. Whereas diphenyliodonium fluoride (1 a) exists as a tetramer bridged by four fluoride ions, 2-methylphenyl(phenyl)iodonium fluoride (2 a) forms a fluoride-bridged dimer that is further halogen bonded to two other monomers. We discuss the topological relationships between the two and their implications for fluorination in solution. Both radiofluorination and NMR spectroscopy show that thermolysis of 2 a gives 2-fluorotoluene and fluorobenzene in a 2 to 1 ratio that is in good agreement with the ratio observed from the radiofluorination of 2-methylphenyl(phenyl)iodonium chloride (2 b). The constancy of the product ratio affirms that the fluorinations occur via the same two rapidly interconverting transition states whose energy difference dictates chemoselectivity. From quantum chemical studies with density functional theory we attribute the "ortho-effect" to the favorable electrostatic interaction between the incoming fluoride and the o-methyl in the transition state. By utilizing the crystal structures of 1 a and 2 a, the mechanisms of fluoroarene formation from diaryliodonium fluorides in their monomeric, homodimeric, heterodimeric, and tetrameric states were also investigated. We propose that oligomerization energy dictates whether the fluorination occurs through a monomeric or an oligomeric pathway.

Efficiently computing pathway free energies: New approaches based on chain-of-replica and Non-Boltzmann Bennett reweighting schemes.

BACKGROUND: Accurately modeling condensed phase processes is one of computation's most difficult challenges. Include the possibility that conformational dynamics may be coupled to chemical reactions, where multiscale (i.e., QM/MM) methods are needed, and this task becomes even more daunting. METHODS: Free energy simulations (i.e., molecular dynamics), multiscale modeling, and reweighting schemes. RESULTS: Herein, we present two new approaches for mitigating the aforementioned challenges. The first is a new chain-of-replica method (off-path simulations, OPS) for computing potentials of mean force (PMFs) along an easily defined reaction coordinate. This development is coupled with a new distributed, highly-parallel replica framework (REPDstr) within the CHARMM package. Validation of these new schemes is carried out on two processes that undergo conformational changes. First is the simple torsional rotation of butane, while a much more challenging glycosidic rotation (in vacuo and solvated) is the second. Additionally, a new approach that greatly improves (i.e., possibly an order of magnitude) the efficiency of computing QM/MM PMFs is introduced and compared to standard schemes. Our efforts are grounded in the recently developed method for efficiently computing QM-based free energies (i.e., QM-Non-Boltzmann Bennett, QM-NBB). Again, we validate this new technique by computing the QM/MM PMF of butane's torsional rotation. CONCLUSIONS: The OPS-REPDstr method is a promising new approach that overcomes many limitations of standard pathway simulations in CHARMM. The combination of QM-NBB with pathway techniques is very promising as it offers significant advantages over current procedures. GENERAL SIGNIFICANCE: Efficiently computing potentials of mean force is a major, unresolved, area of interest. This article is part of a Special Issue entitled Recent developments of molecular dynamics.

Rapid identification of atypical tetracyclines using tandem mass spectrometric fragmentation patterns.

RATIONALE: When applying biosynthetic engineering approaches at the early stages of drug discovery, e.g. aiming to develop novel tetracycline analogues, target compounds are generally produced by engineered microorganisms in low yields. Rapid and reliable identification of metabolites with desired structural modification directly from bacterial cultures is therefore of great importance. METHODS: Structural elucidation of atypical tetracyclines was carried out by fragmentation applying electrospray ionisation tandem mass spectrometry (ESI-MS/MS) (triple quadrupole - linear ion trap; Applied Biosystems 4000 QTRAP) and a high-resolution mass spectrometer (Agilent Technologies 6224 TOF). Fragmentation patterns were obtained either with direct injection or by applying separation of target compounds with high-performance liquid chromatography (HPLC) prior to mass spectrometry. In-source and CID fragmentation were compared. Theoretical calculations of target structures using the Gaussian programme suite were carried out with the aim of strengthening experimental structural elucidation. RESULTS: Recombinant strains of Amycolatopsis sulphurea producing atypical tetracyclines chelocardin, modified chelocardin analogues (9-demethylchelocardin and 2-carboxyamido-2-deacetyl-chelocardin (CDCHD), and anhydrotetracycline (ATC) were analysed by collision-induced dissociation (CID) fragmentation with higher collision energies to yield structurally important fragments which were identified. We have demonstrated that ATC is more prone to fragmentation compared to its epimer, which was further supported by comparison of both structures calculated with ab initio calculations. CONCLUSIONS: We have demonstrated that fragmentation patterns of atypical tetracyclines in CID-MS spectra enable rapid structural elucidation of target metabolites produced by cultures of genetically engineered bacteria. This method is of significant importance for early stages of drug development considering that isolation of target metabolites produced at low concentration is challenging. Copyright © 2015 John Wiley & Sons, Ltd.

Molecular dynamics to enhance structure-based virtual screening on cathepsin B.

Molecular dynamics (MD) and molecular docking are commonly used to study molecular interactions in drug discovery. Most docking approaches consider proteins as rigid, which can decrease the accuracy of predicted docked poses. Therefore MD simulations can be used prior to docking to add flexibility to proteins. We evaluated the contribution of using MD together with docking in a docking study on human cathepsin B, a well-studied protein involved in numerous pathological processes. Using CHARMM biomolecular simulation program and AutoDock Vina molecular docking program, we found, that short MD simulations significantly improved molecular docking. Our results, expressed with the area under the receiver operating characteristic curves, show an increase in discriminatory power i.e. the ability to discriminate active from inactive compounds of molecular docking, when docking is performed to selected snapshots from MD simulations.

Structure-based function prediction of uncharacterized protein using binding sites comparison.

A challenge in structural genomics is prediction of the function of uncharacterized proteins. When proteins cannot be related to other proteins of known activity, identification of function based on sequence or structural homology is impossible and in such cases it would be useful to assess structurally conserved binding sites in connection with the protein's function. In this paper, we propose the function of a protein of unknown activity, the Tm1631 protein from Thermotoga maritima, by comparing its predicted binding site to a library containing thousands of candidate structures. The comparison revealed numerous similarities with nucleotide binding sites including specifically, a DNA-binding site of endonuclease IV. We constructed a model of this Tm1631 protein with a DNA-ligand from the newly found similar binding site using ProBiS, and validated this model by molecular dynamics. The interactions predicted by the Tm1631-DNA model corresponded to those known to be important in endonuclease IV-DNA complex model and the corresponding binding free energies, calculated from these models were in close agreement. We thus propose that Tm1631 is a DNA binding enzyme with endonuclease activity that recognizes DNA lesions in which at least two consecutive nucleotides are unpaired. Our approach is general, and can be applied to any protein of unknown function. It might also be useful to guide experimental determination of function of uncharacterized proteins.

Structural insight into the unique binding properties of pyridylethanol(phenylethyl)amine inhibitor in human CYP51.

Sterol 14α-demethylase (CYP51) is the main drug target for the treatment of fungal infections. The discovery of new efficient fungal CYP51 inhibitors requires an understanding of the structural requirements for selectivity for the fungal over the human ortholog. In this study, a binding mode of the pyridylethanol(phenylethyl)amine type CYP51 inhibitor to the human ortholog was determined at the atomic level. We isolated and purified a full-length human CYP51. The inhibitor-specific binding and its conformational and dynamic properties were evaluated using UV-visible and NMR spectroscopy. Considering the experimental data in docking calculations and molecular dynamics simulations, the location of the inhibitor moieties and their interactions with the enzyme active site were determined. The inhibitor binds to the enzyme in two diastereomeric forms, which have a common location of aromatic ring moieties, while the less bulky propyl chain can adapt to various hydrophobic regions of the enzyme active site. The halogenated phenyl ring binds in the substrate access channel making numerous contacts with the hydrophobic side chains, and its interactions with the unconserved residues are especially informative. The results reveal the unique binding properties of the investigated inhibitor in comparison to the azoles and provide novel directions for the design of selective fungal inhibitors.

Identification of Triazolopyrimidinyl Scaffold SARS-CoV-2 Papain-Like Protease (PLpro) Inhibitor.

The global impact of severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) and its companion disease, COVID-19, has reminded us of the importance of basic coronaviral research. In this study, a comprehensive approach using molecular docking, in vitro assays, and molecular dynamics simulations was applied to identify potential inhibitors for SARS-CoV-2 papain-like protease (PLpro), a key and underexplored viral enzyme target. A focused protease inhibitor library was initially created and molecular docking was performed using CmDock software (v0.2.0), resulting in the selection of hit compounds for in vitro testing on the isolated enzyme. Among them, compound 372 exhibited promising inhibitory properties against PLpro, with an IC50 value of 82 ± 34 µM. The compound also displayed a new triazolopyrimidinyl scaffold not yet represented within protease inhibitors. Molecular dynamics simulations demonstrated the favorable binding properties of compound 372. Structural analysis highlighted its key interactions with PLpro, and we stress its potential for further optimization. Moreover, besides compound 372 as a candidate for PLpro inhibitor development, this study elaborates on the PLpro binding site dynamics and provides a valuable contribution for further efforts in pan-coronaviral PLpro inhibitor development.

Double locking of an Escherichia coli promoter by two repressors prevents premature colicin expression and cell lysis.

The synthesis of Eschericha coli colicins is lethal to the producing cell and is repressed during normal growth by the LexA transcription factor, which is the master repressor of the SOS system for repair of DNA damage. Following DNA damage, LexA is inactivated and SOS repair genes are induced immediately, but colicin production is delayed and induced only in terminally damaged cells. The cause of this delay is unknown. Here we identify the global transcription repressor, IscR, as being directly responsible for the delay in colicin K expression during the SOS response, and identify the DNA target for IscR at the colicin K operon promoter. Our results suggest that, IscR stabilizes LexA at the cka promoter after DNA damage thus, preventing its cleavage and inactivation, and this cooperation ensures that suicidal colicin K production is switched on only as a last resort. A similar mechanism operates at the regulatory region of other colicins and, hence, we suggest that many promoters that control the expression of 'lethal' genes are double locked.

Ligand binding site identification by higher dimension molecular dynamics.

We propose a new molecular dynamics (MD) protocol to identify the binding site of a guest within a host. The method utilizes a four spatial (4D) dimension representation of the ligand allowing for rapid and efficient sampling within the receptor. We applied the method to two different model receptors characterized by diverse structural features of the binding site and different ligand binding affinities. The Abl kinase domain is comprised of a deep binding pocket and displays high affinity for the two chosen ligands examined here. The PDZ1 domain of PSD-95 has a shallow binding pocket that accommodates a peptide ligand involving far fewer interactions and a micromolar affinity. To ensure completely unbiased searching, the ligands were placed in the direct center of the protein receptors, away from the binding site, at the start of the 4D MD protocol. In both cases, the ligands were successfully docked into the binding site as identified in the published structures. The 4D MD protocol is able to overcome local energy barriers in locating the lowest energy binding pocket and will aid in the discovery of guest binding pockets in the absence of a priori knowledge of the site of interaction.

Common force field thermodynamics of cholesterol.

Four different force fields are examined for dynamic characteristics using cholesterol as a case study. The extent to which various types of internal degrees of freedom become thermodynamically relevant is evaluated by means of principal component analysis. More complex degrees of freedom (angle bending, dihedral rotations) show a trend towards force field independence. Moreover, charge assignments for membrane-embedded compounds are revealed to be critical with significant impact on biological reasoning.

Correlating protein hot spot surface analysis using ProBiS with simulated free energies of protein-protein interfacial residues.

A protocol was developed for the computational determination of the contribution of interfacial amino acid residues to the free energy of protein-protein binding. Thermodynamic integration, based on molecular dynamics simulation in CHARMM, was used to determine the free energy associated with single point mutations to glycine in a protein-protein interface. The hot spot amino acids found in this way were then correlated to structural similarity scores detected by the ProBiS algorithm for local structural alignment. We find that amino acids with high structural similarity scores contribute on average -3.19 kcal/mol to the free energy of protein-protein binding and are thus correlated with hot spot residues, while residues with low similarity scores contribute on average only -0.43 kcal/mol. This suggests that the local structural alignment method provides a good approximation of the contribution of a residue to the free energy of binding and is particularly useful for detection of hot spots in proteins with known structures but undetermined protein-protein complexes.

Replica exchanging self-guided Langevin dynamics for efficient and accurate conformational sampling.

This work presents a replica exchanging self-guided Langevin dynamics (RXSGLD) simulation method for efficient conformational searching and sampling. Unlike temperature-based replica exchanging simulations, which use high temperatures to accelerate conformational motion, this method uses self-guided Langevin dynamics (SGLD) to enhance conformational searching without the need to elevate temperatures. A RXSGLD simulation includes a series of SGLD simulations, with simulation conditions differing in the guiding effect and/or temperature. These simulation conditions are called stages and the base stage is one with no guiding effect. Replicas of a simulation system are simulated at the stages and are exchanged according to the replica exchanging probability derived from the SGLD partition function. Because SGLD causes less perturbation on conformational distribution than high temperatures, exchanges between SGLD stages have much higher probabilities than those between different temperatures. Therefore, RXSGLD simulations have higher conformational searching ability than temperature based replica exchange simulations. Through three example systems, we demonstrate that RXSGLD can generate target canonical ensemble distribution at the base stage and achieve accelerated conformational searching. Especially for large systems, RXSGLD has remarkable advantages in terms of replica exchange efficiency, conformational searching ability, and system size extensiveness.

Automatic Assembly and Calibration of Models of Enzymatic Reactions Based on Ordinary Differential Equations.

Enzymatic reactions have been studied for more than a 100 years. Indeed, isolation of enzymes from biological materials is no longer the main source of enzymes today, as they are now largely produced using recombinant technology, or can even be synthesized from scratch. Studies of the three-dimensional structures of enzymes can provide answers to many questions, but the kinetics of enzymatic reactions is the only method that can lead to better understanding of their function. The complexity of high-throughput analysis of progress curves of data obtained can only be achieved through numerical solutions of a suitable system of ordinary differential equations. We have developed the freely available server ENZO: a web tool for derivation and evaluation of kinetic models of enzyme-catalyzed reactions ( http://enzo.cmm.ki.si/ ). ENZO can be used for simultaneous analysis of a series of progress curves obtained under many different conditions. In this chapter, we exemplify the principles and possibilities of this type of high-throughput analysis.

Molecular Dynamics Simulations Reveal Interactions of an IgG1 Antibody With Selected Fc Receptors.

In a survey of novel interactions between an IgG1 antibody and different Fcγ receptors (FcγR), molecular dynamics simulations were performed of interactions of monoclonal antibody involved complexes with FcγRs. Free energy simulations were also performed of isolated wild-type and substituted Fc regions bound to FcγRs with the aim of assessing their relative binding affinities. Two different free energy calculation methods, Molecular Mechanical/Generalized Born Molecular Volume (MM/GBMV) and Bennett Acceptance Ratio (BAR), were used to evaluate the known effector substitution G236A that is known to selectively increase antibody dependent cellular phagocytosis. The obtained results for the MM/GBMV binding affinity between different FcγRs are in good agreement with previous experiments, and those obtained using the BAR method for the complete antibody and the Fc-FcγR simulations show increased affinity across all FcγRs when binding to the substituted antibody. The FcγRIIa, a key determinant of antibody agonistic efficacy, shows a 10-fold increase in binding affinity, which is also consistent with the published experimental results. Novel interactions between the Fab region of the antibody and the FcγRs were discovered with this in silico approach, and provide insights into the antibody-FcγR binding mechanism and show promise for future improvements of therapeutic antibodies for preclinical studies of biological drugs.

Conformational structure and energetics of 2-methylphenyl(2'-methoxyphenyl)iodonium chloride: evidence for solution clusters.

Diaryliodonium salts allow the efficient incorporation of cyclotron-produced [(18)F]fluoride ions into electron-rich and electron-deficient arenes to provide potential radiotracers for molecular imaging in vivo with positron emission tomography (PET). This process (ArI(+)Ar'+(18)F(-)→Ar(18)F+Ar'I) is still not well understood mechanistically. To better understand this and similar reactions, it would be valuable to understand the structures of diaryliodonium salts in organic media, where the reactions are typically conducted. In this endeavor, the X-ray structure of a representative iodonium salt, 2-methylphenyl(2'-methoxyphenyl)iodonium chloride (1), was determined. Our X-ray structure analysis showed 1 to have the conformational M-P dimer as the unit cell with hypervalent iodine as a stereogenic center in each conformer. With the ab initio replica path method we constructed the inversion path between the two enantiomers of 1, thereby revealing two additional pairs of enantiomers that are likely to undergo fast interconversion in solution. Also LC-MS of 1 showed the presence of dimeric and tetrameric anion-bridged clusters in weak organic solution. This observation is consistent with the energetics of 1, both as monomeric and dimeric forms in MeCN, calculated at the B3LYP/DGDZVP level. These evidences of the existence of dimeric and higher order clusters of 1 in solution are relevant to achieve a deeper general understanding of the mechanism and outcome of reactions of diaryliodonium salts in organic media with nucleophiles, such as the [(18)F]fluoride ion.