Metadynamics simulations leveraged by statistical analyses and artificial intelligence-based tools to inform the discovery of G protein-coupled receptor ligands

Abstract

G Protein-Coupled Receptors (GPCRs) are a large family of membrane proteins with pluridimensional signaling profiles. They undergo ligand-specific conformational changes, which in turn lead to the differential activation of intracellular signaling proteins and the consequent triggering of a variety of biological responses. This conformational plasticity directly impacts our understanding of GPCR signaling and therapeutic implications, as do ligand-specific kinetic differences in GPCR-induced transducer activation/coupling or GPCR-transducer complex stability. High-resolution experimental structures of ligand-bound GPCRs in the presence or absence of interacting transducers provide important, yet limited, insights into the highly dynamic process of ligand-induced activation or inhibition of these receptors. We and others have complemented these studies with computational strategies aimed at characterizing increasingly accurate metastable conformations of GPCRs using a combination of metadynamics simulations, state-of-the-art algorithms for statistical analyses of simulation data, and artificial intelligence-based tools. This minireview provides an overview of these approaches as well as lessons learned from them towards the identification of conformational states that may be difficult or even impossible to characterize experimentally and yet important to discover new GPCR ligands.

Introduction

G Protein-Coupled Receptors (GPCRs) are important drug targets consisting of seven membrane-spanning helices connected by alternating intracellular and extracellular loops, and known for transducing extracellular signals across the cell membrane. Evidence accumulated over the past decade has suggested a pluridimensional functionality of these receptors induced by ligands with different efficacies. Specifically, ligands as diverse as photons, small molecules, and peptides can stabilize different receptor conformations which, in turn, can trigger the activation of different effectors, such as several subtypes of heterotrimeric G proteins, β-arrestins, and G protein-coupled receptor kinases (1). This differential engagement of transducer subtypes with varying signal magnitudes is a phenomenon known as “functional selectivity” or “biased agonism” (2–7), and is a prerequisite for triggering therapeutic or unwanted biological effects of GPCRs via activation of different downstream signaling pathways. However, it is not the only element at play. Temporal analyses of ligand binding to GPCRs, as well as GPCR-induced G protein activation/coupling and β-arrestin recruitment evoked by ligands with different efficacies, demonstrate the existence of another dimension of functional bias by GPCRs that directly impacts our understanding of GPCR signaling and therapeutic implications, thus suggesting the importance of incorporating quantifications of ligand binding and signaling kinetics in modern drug discovery efforts (8–10).

The mechanistic and kinetic bases of GPCR-mediated functional selectivity are poorly understood, notwithstanding the past decade’s technological advances in GPCR functional and structural biology. Among them are the availability of genetically-encoded biosensors for the optical detection of signals from a variety of transduction molecules in several cell types, tissues, and whole organisms (11), as well as revolutionary methodological developments in X-ray crystallography and cryogenic electron microscopy (cryo-EM). The latter have allowed to solve high-resolution experimental structures for 140 unique GPCRs, 520 unique ligand-GPCR complexes, 95 unique GPCRs in complex with heterotrimeric G proteins, 6 unique GPCRs in complex with arrestins, and 1 unique GPCR in complex with G protein-coupled receptor kinases (data retrieved from GPCRdb (12) and the Protein Data Bank (13) on 11/4/2022). Although these structures have provided important insights into ligand-GPCR and GPCR-transducer interactions, they are heavily engineered static snapshots. Thus, inferences of long-distance conformational changes propagating from the ligand binding site in the receptor to its cellular signaling partners, and eventually translating into specific physiological cell responses, remain highly speculative.

To probe the conformational heterogeneity of GPCRs and GPCR complexes and extend information from high-resolution structural methods, researchers have resorted to biophysical techniques such as nuclear magnetic resonance (14–21), double electron-electron resonance spectroscopy (22), hydrogen/deuterium exchange mass spectroscopy (23, 24), and single-molecule fluorescence resonance energy transfer (25–27). However, current technical challenges prevent these techniques from achieving atomic-level precision for the entire GPCR alone or in complex with their natural cellular signaling partners. Although molecular dynamics (MD) simulations can provide a critical bridge between the atomic-level insight from high-resolution structural methods and molecular motions (28), standard MD algorithms limit dynamic explorations to timescales that are shorter than most biological processes notwithstanding their use of massively-parallel high-performance computing platforms. Several enhanced conformational sampling methods have been put forward to overcome these limitations (29), and our group pioneered the use of one of them, i.e., metadynamics (MetaD) [see Box 1 and (30)], for more efficient studies of the conformational plasticity of liganded or unliganded GPCRs embedded in a lipid mimetic environment. The underlying principles of MetaD are that (a) the process to be investigated can be described by a small number of reaction coordinates (collective variables, CVs) and (b) the sum of destabilizing Gaussian potentials that are added to penalize sampled conformational states for faster simulation convergence is the mirror image of the free-energy profile of said process. Although MetaD’s main problem remains that of identifying the right CVs to efficiently study the dynamic process of interest, MetaD-based strategies are increasingly utilized nowadays to study GPCR ligand binding, conformational dynamics, and kinetics (see Figure 1 for illustrations of selected strategies), thus offering an outlet where one could exploit uniquely characterized metastable states of GPCRs as targets for the discovery of functionally selective ligands.

Figure 1

This minireview provides a general overview of MetaD strategies used alone or in combination with algorithms for statistical analyses of simulation data [e.g., Markov State Models (MSMs) (32, 33), information theory-based methods (34), and transfer entropy approaches (35)] and machine learning/artificial intelligence (AI)-based tools (36) to not only characterize ligand-specific conformational states of GPCRs that are difficult to resolve experimentally but also to derive information that can help expedite the GPCR drug discovery process. This includes structural determinants of allosteric communications from the causality of correlated motions, as well as ligand-specific kinetic elements of activation [e.g., (37, 38)]. The computational strategies and applications discussed herein are not exhaustive but include representative examples that have affirmed the power of MetaD in the study of various GPCR processes, including ligand binding, receptor activation, and related kinetics. The different simulation setups and protocols of the studies reported in this minireview are summarized in Table 1.

MetaD-based strategies for the prediction of GPCR-ligand binding modes

To the best of our knowledge, the very first application of MetaD to GPCRs was published by our group more than a decade ago (31) and aimed at efficiently exploring the molecular mechanisms of GPCR ligand recognition (see Figure 1A for an example of output information). Specifically, we studied the free binding of the antagonist naloxone from the water environment to a homology model of the δ-opioid receptor (DOR) based on the X-ray crystal structure of the β₂-adrenergic receptor (β₂AR) as there were no available crystal structures of opioid receptors at the time (see Table 1 for details of system setup and simulations). Although one would have ideally used a single CV to contain the computational cost, which scales exponentially with the number of CVs (55), we could not achieve simulation convergence for this system with less than two CVs and a half microsecond of well-tempered (WT) MetaD simulations using multiple walkers (MW) (56). These simulations revealed, for the first time, unprecedented details of the free binding of a ligand to a GPCR, including the various intermediate energetic states visited by the ligand and their corresponding ligand-receptor interactions. Moreover, by restricting ligand sampling in the bulk region using an approach that allowed the ligand to only move in a conical region centered at the center of mass (COM) of the binding pocket (57), and applying the appropriate correction to the calculated free energy, we could derive ligand binding affinity estimates that were close to experimental values (58–61).

Subsequent applications by our group focused on predicting optimal ligand binding modes for atypical opioid ligands (42), as well as opioid allosteric modulators (41), which would be targeting regions on the receptor that are not conserved and highly flexible. Ligand binding in these studies was described by CVs accounting for the relative position and orientation of the ligands, as well as the number of interactions they established with the receptor. Although these simulations allowed for a thorough and efficient exploration of the dynamic process of ligand-GPCR binding at either orthosteric or allosteric sites (62, 63), they were computationally quite expensive as they required up to 7 μs of simulation time to reach convergence (42).

To ensure faster simulation convergence, a universal single CV that uses the highly conserved Trp^6.48 [superscript refers to the Ballesteros-Weinstein numbering scheme for GPCRs (64)] at the base of the orthosteric binding pocket and the orientation of the receptor in the membrane was proposed for ligand binding/unbinding to/from class A GPCRs using WT-MetaD and a funnel restraint (65) that limits the conformational sampling of the ligand in the bulk water (39) (Table 1). Funnel MetaD (FM) has been successfully applied to a number of GPCR systems (62), and its use made more easily accessible by a recently reported graphical user interface (GUI)-based protocol termed FM Advanced Protocol (FMAP) (66). Also interesting is a protocol that has recently enabled the successful prediction of preferential binding modes of different GPCR systems using conformational ensembles derived from the clustering of MW-MetaD simulations (40). These strategies, however, are still not amenable for high-throughput given their computational cost. If estimates of protein-ligand binding free energies are not required but the goal is focused on assessing the relative stability of binding mode predictions, a reasonable alternative we and others have used in applications to GPCRs [e.g., (43) for a recent example from our lab] is a combination of induced-fit docking and MetaD (67). Based on these growing examples and the ready accessibility of these strategies via user-friendly graphical interfaces, we expect that MetaD-based approaches will become a standard tool for the prediction of GPCR ligand binding in the future, and will make a real impact on drug discovery.

MetaD-based strategies to probe the activation landscape of GPCRs

We also were the first to study the activation pathway of a GPCR, specifically rhodopsin (RHO), using a combination of MetaD with adiabatic biased MD simulations (51). In particular, we used the MetaD variant known as path-sampling MetaD (PS-MetaD) to reconstruct the system’s free-energy landscape along predetermined transition trajectories between receptor inactive and active states as a function of the position along the path (s) and the distance (z) from it (Figure 1B). Our results suggested that at least four metastable macrostates containing receptor conformations with a different amplitude of the outward movement of transmembrane (TM) helix 6 are sampled by RHO in the transition from inactive to active conformations and these are connected by at least two different pathways (51). The conformations of two of these macrostates were very close to the available inactive and active experimental structures of RHO, whereas the other two macrostates contained conformations representing intermediate states. Notably, subsequent MetaD simulations we carried out on the β₂AR (46) (see section below for additional details and Table 1) confirmed the presence of one different inactive-like and one different active-like macrostates in addition to those whose conformations closely resembled the experimentally known inactive and active structures of the receptor. This finding was interesting since standard MD simulations of the β₂AR at the time had only been able to identify one intermediate state in terms of TM3-TM6 separation (68). Similar results were also obtained for the μ-opioid receptor (MOR) by carrying out adiabatic biased MD simulations and PS-MetaD simulations (38), as well as more expensive adaptive-sampling MD simulations (44). In contrast, only one intermediate state between fully inactive and active conformations was identified for the class B glucagon receptor (GCGR) (50), using parallel tempering WT-MetaD and CVs representing two linear combinations of the alpha-carbon root mean square deviation (RMSD) of TM6 to the inactive conformation of GCGR and to the active, closely related, glucagon-like peptide 1 receptor (69, 70). An interesting observation of this study based on a comparison between the computed conformational free-energy landscape associated with the activation of the receptor-agonist complex and that of the receptor-agonist-G protein complex was that the agonist stabilizes the receptor in a preactivated complex before its full activation is achieved by G protein binding (50).

MetaD was also recently used to propose the activation mechanisms of several prototypic class A (45) and class C (49) receptors. These mechanisms were class-specific and revealed an active role of the G protein in promoting conformational changes of the receptor.

MetaD-based strategies for the characterization of ligand-specific GPCR conformations

Being able to predict the receptor conformations a drug can stabilize is considered the “holy grail” in the drug discovery field as it might inform the drug’s biological outcome. With this in mind, we developed a computational strategy to study the ligand-induced modulation of the free-energy landscape of GPCRs (46) (see Figure 1C for an example of output information). Specifically, we used WT-MetaD  (56) to identify ligand-specific metastable states of the β₂AR, along pre-determined activation pathways between its high-resolution inactive and active structures using adiabatic biased MD (Table 1). The results confirmed a tendency of ligands to stabilize an inactive or active conformation of the receptor depending on their efficacy (46).

That ligands with different efficacies shift the equilibrium between inactive and active states was recently also reported for the human adenosine A2A receptor using MetaD (48). Specific conformational rearrangements of key structural elements were responsible for this shift, including rotameric changes of the conserved Trp^6.48 residue. This movement can be induced by agonists but not inverse agonists, and appears to correlate with the opening of the G protein binding site via disruption of hydrophobic packing which causes TM6 and TM5 to move away from TM3.

MetaD was also recently used to study the interplay between GPCR ligand binding and the coupling of arrestins or G proteins (47). In addition to estimate the binding free energies of ligands with different efficacies at the β₂AR in the presence or absence of Gs or β-arrestin, a combination of WT-MetaD and FM was used to study the free energy landscape of activation and transducer coupling as a function of the distance between the receptor alpha-carbons of Arg^3.50 and Leu^6.34 and the distance between Arg^3.50 and the alpha-carbons of Gαs Glu392 or β-arrestin Val71. This setup (Table 1) allowed to quantify the ligand/transducer cooperative effects on the activation of the β₂AR, thus providing a simple way to predict the functional bias of a ligand (47).

MetaD-based strategies for the prediction of ligand binding and activation kinetics

Knowledge of GPCR ligand binding kinetics, and especially drug-target residence times and involved molecular determinants, is highly desirable because it can serve as an effective guiding principle to select candidate molecules for clinical development (10). However, this collective information is difficult to obtain because on and off rates depend on the height of the highest energy barrier of the transition state between bound and unbound conformations, and this height is difficult to estimate both experimentally and computationally. MD simulations can help study ligand binding at the active site, but binding is a rare event on microscopic timescales, and as such, it is difficult to sample. Although microsecond-scale MD simulations allow to obtain multiple binding events at a GPCR embedded in an explicit lipid-water environment, from which it is possible to derive k_on estimates (71), much longer timescales would be required for the dissociation of the ligand from a GPCR orthosteric site, making it very difficult to derive k_off estimates from unbiased, standard MD simulations. To attempt more efficient predictions of dissociation rates from a GPCR, we recently designed a computational strategy that combines machine learning and infrequent MetaD (37). Specifically, the strategy used the automatic mutual information noise omission (AMINO) algorithm (72) to enable the robust and automated selection of non-redundant molecular features from unbiased MD simulations that can then be used in the reweighted autoencoded variational Bayes for enhanced sampling (RAVE) method (73, 74) to learn optimal CVs. These CVs were used in infrequent MetaD simulations and allowed to efficiently study the unbinding kinetics of two opioid drugs with very different residence times (morphine and buprenorphine) and to predict dissociation rates for these drugs from the MOR that were within one order of magnitude from experimental values (75). The simulations also gave structural information about rate-limiting transition states and metastable poses that can in principle be used in the design of drugs with a desired kinetic profile. One important lesson from these simulations was that accurate rate estimates require consideration of ligand solvation, as also concluded by others (52, 65).

Recently, MetaD was also used to investigate the dissociation of a potent drug targeting the MOR, specifically fentanyl, from the receptor (54). Not only did these simulations provide information about fentanyl’s binding kinetics and mechanism, but they suggested a role for the protonation state of the conserved residue His^6.52 in modulating fentanyl’s affinity and binding pose. While a recent study (76) reporting cryo-EM structures of the MOR bound to fentanyl showed that a His(6.52)Ala mutation reduced fentanyl-induced MOR activation of G protein and β-arrestin signaling, it did not support the predicted binding poses of fentanyl using MetaD or other computational methods (53, 54, 77).

Ligand binding kinetics is however not the only temporal aspect that is important to understand GPCR signaling. To tackle the kinetics of GPCR activation, we designed a computational strategy that uses a combination of PS-MetaD and the maximum caliber principle to not only thoroughly characterize the conformational space sampled by the receptor along its activation path, but also to derive unbiased kinetic rates from simulations with transitions accelerated by a bias potential (38) (Figure 1D). We showed that the approach could efficiently sample conformational transitions between the crystal structures of inactive and active GPCRs, specifically the MOR, at a ∼2 orders of magnitude computational cost reduction with respect to a more expensive high-throughput MD adaptive sampling protocol run on distributed computational resources. We also demonstrated that the strategy yielded thermodynamic and kinetic properties of MOR activation, as well as information about ligand-induced allosteric communications across the receptor when combined with n-body information theory (34), that were in agreement with those obtained by more computationally expensive adaptive sampling protocols.

Conclusions

Identifying the molecular determinants that underlie the different dynamic and kinetic behaviors of GPCRs is likely to provide fundamental information for drug discovery, and is therefore the focus of much research. However, obtaining this information is not a trivial undertaking, either experimentally or computationally. While strategies combining different variants of MetaD with state-of-the-art algorithms for statistical analyses of simulation data and/or AI-based tools have made strides in the study of GPCR ligand binding, activation, and related kinetics, the majority of studies published to date are limited to the receptors alone, although recent studies are focusing more and more on receptor complexes and the combined effect of ligand and transducer on the dynamic behavior of the receptor. To better understand the underlying principles of GPCR biased signaling, attention will probably need to be shifted more towards the appropriate characterization of the dynamic behavior of interacting signaling proteins, thus providing important insight into receptor-transducer complex stability and/or receptor-transducer activation kinetics. These studies will most likely require developing more sophisticated MetaD strategies or hybrid approaches and a community effort to disseminate best practices for MetaD simulations, as well as system setups, parameters, protocols, etc. to more expeditiously advance knowledge of GPCR dynamics and their relation to functional selectivity.

References

• 1

  GurevichVVGurevichEV. Molecular mechanisms of gpcr signaling: A structural perspective. Int J Mol Sci (2017) 18(12):2519. doi: 10.3390/ijms18122519

  • CrossRef
  • Google Scholar

• 2

  BohnLM. Selectivity for G protein or arrestin-mediated signaling. In: K N, editor. Functional selectivity of G protein-coupled receptor ligands. Totowa, NJ: Humana Press (2009). p. 71–85.

  • Google Scholar

• 3

  KenakinT. New concepts in drug discovery: Collateral efficacy and permissive antagonism. Nat Rev Drug Discovery (2005) 4(11):919–27. doi: 10.1038/nrd1875

  • CrossRef
  • Google Scholar

• 4

  KenakinT. Collateral efficacy in drug discovery: Taking advantage of the good (Allosteric) nature of 7tm receptors. Trends Pharmacol Sci (2007) 28(8):407–15. doi: 10.1016/j.tips.2007.06.009

  • CrossRef
  • Google Scholar

• 5

  KenakinT. Biased agonism. F1000 Biol Rep (2009) 1:87. doi: 10.3410/B1-87

  • CrossRef
  • Google Scholar

• 6

  MailmanRB. Gpcr functional selectivity has therapeutic impact. Trends Pharmacol Sci (2007) 28(8):390–6. doi: 10.1016/j.tips.2007.06.002

  • CrossRef
  • Google Scholar

• 7

  UrbanJDClarkeWPvon ZastrowMNicholsDEKobilkaBWeinsteinHet al. Functional selectivity and classical concepts of quantitative pharmacology. J Pharmacol Exp Ther (2007) 320(1):1–13. doi: 10.1124/jpet.106.104463

  • CrossRef
  • Google Scholar

• 8

  HoareSRJTewsonPHQuinnAMHughesTEBridgeLJ. Analyzing kinetic signaling data for G-Protein-Coupled receptors. Sci Rep (2020) 10(1):12263. doi: 10.1038/s41598-020-67844-3

  • CrossRef
  • Google Scholar

• 9

  HoareSRJTewsonPHSachdevSConnorMHughesTEQuinnAM. Quantifying the kinetics of signaling and arrestin recruitment by nervous system G-protein coupled receptors. Front Cell Neurosci (2021) 15:814547. doi: 10.3389/fncel.2021.814547

  • CrossRef
  • Google Scholar

• 10

  SwinneyDCHaubrichBAVan LiefdeIVauquelinG. The role of binding kinetics in gpcr drug discovery. Curr Top Med Chem (2015) 15(24):2504–22. doi: 10.2174/1568026615666150701113054

  • CrossRef
  • Google Scholar

• 11

  WrightSCBouvierM. Illuminating the complexity of gpcr pathway selectivity - advances in biosensor development. Curr Opin Struct Biol (2021) 69:142–9. doi: 10.1016/j.sbi.2021.04.006

  • CrossRef
  • Google Scholar

• 12

  KooistraAJMordalskiSPandy-SzekeresGEsguerraMMamyrbekovAMunkCet al. Gpcrdb in 2021: Integrating gpcr sequence, structure and function. Nucleic Acids Res (2021) 49(D1):D335–D43. doi: 10.1093/nar/gkaa1080

  • CrossRef
  • Google Scholar

• 13

  BermanHMWestbrookJFengZGillilandGBhatTNWeissigHet al. The protein data bank. Nucleic Acids Res (2000) 28(1):235–42. doi: 10.1093/nar/28.1.235

  • CrossRef
  • Google Scholar

• 14

  CasiraghiMPointEPozzaAMoncoqKBaneresJLCatoireLJ. Nmr analysis of gpcr conformational landscapes and dynamics. Mol Cell Endocrinol (2019) 484:69–77. doi: 10.1016/j.mce.2018.12.019

  • CrossRef
  • Google Scholar

• 15

  CongXMaurelDDemeneHVasiliauskaite-BrooksIHagelbergerJPeyssonFet al. Molecular insights into the biased signaling mechanism of the mu-opioid receptor. Mol Cell (2021) 81(20):4165–75 e6. doi: 10.1016/j.molcel.2021.07.033

  • CrossRef
  • Google Scholar

• 16

  EddyMTLeeMYGaoZGWhiteKLDidenkoTHorstRet al. Allosteric coupling of drug binding and intracellular signaling in the A2a adenosine receptor. Cell (2018) 172(1-2):68–80 e12. doi: 10.1016/j.cell.2017.12.004

  • CrossRef
  • Google Scholar

• 17

  HuangSKPandeyATranDPVillanuevaNLKitaoASunaharaRKet al. Delineating the conformational landscape of the adenosine A2a receptor during G protein coupling. Cell (2021) 184(7):1884–94 e14. doi: 10.1016/j.cell.2021.02.041

  • CrossRef
  • Google Scholar

• 18

  OkudeJUedaTKofukuYSatoMNobuyamaNKondoKet al. Identification of a conformational equilibrium that determines the efficacy and functional selectivity of the mu-opioid receptor. Angew Chem Int Ed Engl (2015) 54(52):15771–6. doi: 10.1002/anie.201508794

  • CrossRef
  • Google Scholar

• 19

  PicardLPProsserRS. Advances in the study of gpcrs by (19)F nmr. Curr Opin Struct Biol (2021) 69:169–76. doi: 10.1016/j.sbi.2021.05.001

  • CrossRef
  • Google Scholar

• 20

  SounierRMasCSteyaertJLaeremansTManglikAHuangWet al. Propagation of conformational changes during mu-opioid receptor activation. Nature (2015) 524(7565):375–8. doi: 10.1038/nature14680

  • CrossRef
  • Google Scholar

• 21

  XuJHuYKaindlJRiselPHubnerHMaedaSet al. Conformational complexity and dynamics in a muscarinic receptor revealed by nmr spectroscopy. Mol Cell (2019) 75(1):53–65 e7. doi: 10.1016/j.molcel.2019.04.028

  • CrossRef
  • Google Scholar

• 22

  ElgetiMHubbellWL. Deer analysis of gpcr conformational heterogeneity. Biomolecules (2021) 11(6):778. doi: 10.3390/biom11060778

  • CrossRef
  • Google Scholar

• 23

  DuYDucNMRasmussenSGFHilgerDKubiakXWangLet al. Assembly of a gpcr-G protein complex. Cell (2019) 177(5):1232–42 e11. doi: 10.1016/j.cell.2019.04.022

  • CrossRef
  • Google Scholar

• 24

  LiuXXuXHilgerDAschauerPTiemannJKSDuYet al. Structural insights into the process of gpcr-G protein complex formation. Cell (2019) 177(5):1243–51 e12. doi: 10.1016/j.cell.2019.04.021

  • CrossRef
  • Google Scholar

• 25

  AsherWBTerryDSGregorioGGAKahsaiAWBorgiaAXieBet al. Gpcr-mediated beta-arrestin activation deconvoluted with single-molecule precision. Cell (2022) 185(10):1661–75 e16. doi: 10.1016/j.cell.2022.03.042

  • CrossRef
  • Google Scholar

• 26

  KaukMHoffmannC. Intramolecular and intermolecular fret sensors for gpcrs - monitoring conformational changes and beyond. Trends Pharmacol Sci (2018) 39(2):123–35. doi: 10.1016/j.tips.2017.10.011

  • CrossRef
  • Google Scholar

• 27

  QuastRBMargeatE. Studying gpcr conformational dynamics by single molecule fluorescence. Mol Cell Endocrinol (2019) 493:110469. doi: 10.1016/j.mce.2019.110469

  • CrossRef
  • Google Scholar

• 28

  LatorracaNRVenkatakrishnanAJDrorRO. Gpcr dynamics: Structures in motion. Chem Rev (2017) 117(1):139–55. doi: 10.1021/acs.chemrev.6b00177

  • CrossRef
  • Google Scholar

• 29

  AbrolRSerranoESantiagoLJ. Development of enhanced conformational sampling methods to probe the activation landscape of gpcrs. Adv Protein Chem Struct Biol (2022) 128:325–59. doi: 10.1016/bs.apcsb.2021.11.001

  • CrossRef
  • Google Scholar

• 30

  LaioAParrinelloM. Escaping free-energy minima. P Natl Acad Sci USA (2002) 99(20):12562–6. doi: 10.1073/pnas.202427399

  • CrossRef
  • Google Scholar

• 31

  ProvasiDBortolatoAFilizolaM. Exploring molecular mechanisms of ligand recognition by opioid receptors with metadynamics. Biochemistry (2009) 48(42):10020–9. doi: 10.1021/bi901494n

  • CrossRef
  • Google Scholar

• 32

  Pérez-HernándezGPaulFGiorginoTDe FabritiisGNoéF. Identification of slow molecular order parameters for Markov model construction. J Chem Phys (2013) 139:15102. doi: 10.1063/1.4811489

  • CrossRef
  • Google Scholar

• 33

  SchwantesCRPandeVS. Improvements in Markov state model construction reveal many non-native interactions in the folding of Ntl9. J Chem Theory Comput (2013) 9(4):2000–9. doi: 10.1021/ct300878a

  • CrossRef
  • Google Scholar

• 34

  LeVineMVWeinsteinH. Nbit–a new information theory-based analysis of allosteric mechanisms reveals residues that underlie function in the leucine transporter leut. PloS Comput Biol (2014) 10(5):e1003603. doi: 10.1371/journal.pcbi.1003603

  • CrossRef
  • Google Scholar

• 35

  SchreiberT. Measuring information transfer. Phys Rev Lett (2000) 85:461. doi: 10.1103/PhysRevLett.85.461

  • CrossRef
  • Google Scholar

• 36

  NoeFTkatchenkoAMullerKRClementiC. Machine learning for molecular simulation. Annu Rev Phys Chem (2020) 71:361–90. doi: 10.1146/annurev-physchem-042018-052331

  • CrossRef
  • Google Scholar

• 37

  Lamim RibeiroJMProvasiDFilizolaM. A combination of machine learning and infrequent metadynamics to efficiently predict kinetic rates, transition states, and molecular determinants of drug dissociation from G protein-coupled receptors. J Chem Phys (2020) 153(12):124105. doi: 10.1063/5.0019100

  • CrossRef
  • Google Scholar

• 38

  MeralDProvasiDFilizolaM. An efficient strategy to estimate thermodynamics and kinetics of G protein-coupled receptor activation using metadynamics and maximum caliber. J Chem Phys (2018) 149(22):224101. doi: 10.1063/1.5060960

  • CrossRef
  • Google Scholar

• 39

  SalehNIbrahimPSaladinoGGervasioFLClarkT. An efficient metadynamics-based protocol to model the binding affinity and the transition state ensemble of G-Protein-Coupled receptor ligands. J Chem Inf Model (2017) 57(5):1210–7. doi: 10.1021/acs.jcim.6b00772

  • CrossRef
  • Google Scholar

• 40

  SoldnerCAHornAHCStichtH. A metadynamics-based protocol for the determination of gpcr-ligand binding modes. Int J Mol Sci (2019) 20(8):1970. doi: 10.3390/ijms20081970

  • CrossRef
  • Google Scholar

• 41

  ShangYYeatmanHRProvasiDAltAChristopoulosACanalsMet al. Proposed mode of binding and action of positive allosteric modulators at opioid receptors. ACS Chem Biol (2016) 11(5):1220–9. doi: 10.1021/acschembio.5b00712

  • CrossRef
  • Google Scholar

• 42

  CrowleyRSRileyAPSherwoodAMGroerCEShivaperumalNBiscaiaMet al. Synthetic studies of neoclerodane diterpenes from salvia divinorum: Identification of a potent and centrally acting mu opioid analgesic with reduced abuse liability. J Med Chem (2016) 59(24):11027–38. doi: 10.1021/acs.jmedchem.6b01235

  • CrossRef
  • Google Scholar

• 43

  ZhouYRamseySProvasiDEl DaibaniAAppourchauxKChakrabortySet al. Predicted mode of binding to and allosteric modulation of the mu-opioid receptor by kratom's alkaloids with reported antinociception in vivo. Biochemistry (2021) 60(18):1420–9. doi: 10.1021/acs.biochem.0c00658

  • CrossRef
  • Google Scholar

• 44

  KapoorAProvasiDFilizolaM. Atomic-level characterization of the methadone-stabilized active conformation of mu-opioid receptor. Mol Pharmacol (2020) 98(4):475–86. doi: 10.1124/mol.119.119339

  • CrossRef
  • Google Scholar

• 45

  MafiAKimSKGoddardWA3rd. The mechanism for ligand activation of the gpcr-G protein complex. Proc Natl Acad Sci U.S.A. (2022) 119(18):e2110085119. doi: 10.1073/pnas.2110085119

  • CrossRef
  • Google Scholar

• 46

  ProvasiDArtachoMCNegriAMobarecJCFilizolaM. Ligand-induced modulation of the free-energy landscape of G protein-coupled receptors explored by adaptive biasing techniques. PloS Comput Biol (2011) 7(10):e1002193. doi: 10.1371/journal.pcbi.1002193

  • CrossRef
  • Google Scholar

• 47

  SalehNSaladinoGGervasioFLClarkT. Investigating allosteric effects on the functional dynamics of Beta2-adrenergic ternary complexes with enhanced-sampling simulations. Chem Sci (2017) 8(5):4019–26. doi: 10.1039/c6sc04647a

  • CrossRef
  • Google Scholar

• 48

  LiJJonssonALBeumingTShelleyJCVothGA. Ligand-dependent activation and deactivation of the human adenosine a(2a) receptor. J Am Chem Soc (2013) 135(23):8749–59. doi: 10.1021/ja404391q

  • CrossRef
  • Google Scholar

• 49

  YangMYKimSKGoddardWA3rd. G Protein coupling and activation of the metabotropic gabab heterodimer. Nat Commun (2022) 13(1):4612. doi: 10.1038/s41467-022-32213-3

  • CrossRef
  • Google Scholar

• 50

  MattediGAcosta-GutierrezSClarkTGervasioFL. A combined activation mechanism for the glucagon receptor. Proc Natl Acad Sci U.S.A. (2020) 117(27):15414–22. doi: 10.1073/pnas.1921851117

  • CrossRef
  • Google Scholar

• 51

  ProvasiDFilizolaM. Putative active states of a prototypic G-Protein-Coupled receptor from biased molecular dynamics. Biophys J (2010) 98(10):2347–55. doi: 10.1016/j.bpj.2010.01.047

  • CrossRef
  • Google Scholar

• 52

  DeganuttiGZhukovADeflorianFFedericoSSpallutoGCookeRMet al. Impact of protein-ligand solvation and desolvation on transition state thermodynamic properties of adenosine A2a ligand binding kinetics. In Silico Pharmacol (2017) 5(1):16. doi: 10.1007/s40203-017-0037-x

  • CrossRef
  • Google Scholar

• 53

  VoQNMahinthichaichanPShenJEllisCR. How mu-opioid receptor recognizes fentanyl. Nat Commun (2021) 12(1):984. doi: 10.1038/s41467-021-21262-9

  • CrossRef
  • Google Scholar

• 54

  MahinthichaichanPVoQNEllisCRShenJ. Kinetics and mechanism of fentanyl dissociation from the mu-opioid receptor. JACS Au (2021) 1(12):2208–15. doi: 10.1021/jacsau.1c00341

  • CrossRef
  • Google Scholar

• 55

  ValssonOTiwaryPParrinelloM. Enhancing important fluctuations: Rare events and metadynamics from a conceptual viewpoint. Annu Rev Phys Chem (2016) 67:159–84. doi: 10.1146/annurev-physchem-040215-112229

  • CrossRef
  • Google Scholar

• 56

  BarducciABussiGParrinelloM. Well-tempered metadynamics: A smoothly converging and tunable free-energy method. Phys Rev Lett (2008) 100(2):20603. doi: 10.1103/PhysRevLett.100.020603

  • CrossRef
  • Google Scholar

• 57

  AllenTWAndersenOSRouxB. Energetics of ion conduction through the gramicidin channel. Proc Natl Acad Sci U.S.A. (2004) 101(1):117–22. doi: 10.1073/pnas.2635314100

  • CrossRef
  • Google Scholar

• 58

  SchlechtingenGDeHavenRNDaubertJDCasselJAChungNNSchillerPWet al. Structure-activity relationships of dynorphin a analogues modified in the address sequence. J Med Chem (2003) 46(11):2104–9. doi: 10.1021/jm020125+

  • CrossRef
  • Google Scholar

• 59

  TollLBerzetei-GurskeIPPolgarWEBrandtSRAdapaIDRodriguezLet al. Standard binding and functional assays related to medications development division testing for potential cocaine and opiate narcotic treatment medications. NIDA Res Monogr (1998) 178:440–66.

  • Google Scholar

• 60

  Tryoen-TothPDecaillotFMFilliolDBefortKLazarusLHSchillerPWet al. Inverse agonism and neutral antagonism at wild-type and constitutively active mutant delta opioid receptors. J Pharmacol Exp Ther (2005) 313(1):410–21. doi: 10.1124/jpet.104.077321

  • CrossRef
  • Google Scholar

• 61

  ValenzanoKJMillerWChenZShanSCrumleyGVictorySFet al. Dipoa ([8-(3,3-Diphenyl-Propyl)-4-Oxo-1-Phenyl-1,3,8-Triazaspiro[4.5]Dec-3-Yl]-Acetic acid), a novel, systemically available, and peripherally restricted mu opioid agonist with antihyperalgesic activity: I. In vitro pharmacological characterization and pharmacokinetic properties. J Pharmacol Exp Ther (2004) 310(2):783–92. doi: 10.1124/jpet.103.063313

  • CrossRef
  • Google Scholar

• 62

  IbrahimPClarkT. Metadynamics simulations of ligand binding to gpcrs. Curr Opin Struct Biol (2019) 55:129–37. doi: 10.1016/j.sbi.2019.04.002

  • CrossRef
  • Google Scholar

• 63

  SchneiderSProvasiDFilizolaM. The dynamic process of drug-gpcr binding at either orthosteric or allosteric sites evaluated by metadynamics. Methods Mol Biol (2015) 1335:277–94. doi: 10.1007/978-1-4939-2914-6_18

  • CrossRef
  • Google Scholar

• 64

  BallesterosJAWeinsteinH. Integrated methods for the construction of three-dimensional models and computational probing of structure-function relations in G protein-coupled receptors. Methods Neurosci (1995) 25:366–428. doi: 10.1016/S1043-9471(05)80049-7

  • CrossRef
  • Google Scholar

• 65

  LimongelliVMarinelliLCosconatiSLa MottaCSartiniSMugnainiLet al. Sampling protein motion and solvent effect during ligand binding. Proc Natl Acad Sci U.S.A. (2012) 109(5):1467–72. doi: 10.1073/pnas.1112181108

  • CrossRef
  • Google Scholar

• 66

  RanioloSLimongelliV. Ligand binding free-energy calculations with funnel metadynamics. Nat Protoc (2020) 15(9):2837–66. doi: 10.1038/s41596-020-0342-4

  • CrossRef
  • Google Scholar

• 67

  ClarkAJTiwaryPBorrelliKFengSMillerEBAbelRet al. Prediction of protein-ligand binding poses Via a combination of induced fit docking and metadynamics simulations. J Chem Theory Comput (2016) 12(6):2990–8. doi: 10.1021/acs.jctc.6b00201

  • CrossRef
  • Google Scholar

• 68

  DrorROArlowDHMaragakisPMildorfTJPanACXuHet al. Activation mechanism of the Beta2-adrenergic receptor. Proc Natl Acad Sci U.S.A. (2011) 108(46):18684–9. doi: 10.1073/pnas.1110499108

  • CrossRef
  • Google Scholar

• 69

  BussiGGervasioFLLaioAParrinelloM. Free-energy landscape for beta hairpin folding from combined parallel tempering and metadynamics. J Am Chem Soc (2006) 128(41):13435–41. doi: 10.1021/ja062463w

  • CrossRef
  • Google Scholar

• 70

  DeighanMBonomiMPfaendtnerJ. Efficient simulation of explicitly solvated proteins in the well-tempered ensemble. J Chem Theory Comput (2012) 8(7):2189–92. doi: 10.1021/ct300297t

  • CrossRef
  • Google Scholar

• 71

  DrorROPanACArlowDHBorhaniDWMaragakisPShanYet al. Pathway and mechanism of drug binding to G-Protein-Coupled receptors. Proc Natl Acad Sci U.S.A. (2011) 108(32):13118–23. doi: 10.1073/pnas.1104614108

  • CrossRef
  • Google Scholar

• 72

  RavindraPSmithZTiwaryP. Automatic mutual information noise omission (Amino): Generating order parameters for molecular systems. Mol Syst Des Eng (2020) 5:339–48. doi: 10.1039/C9ME00115H

  • CrossRef
  • Google Scholar

• 73

  WangYRibeiroJMLTiwaryP. Past-future information bottleneck for sampling molecular reaction coordinate simultaneously with thermodynamics and kinetics. Nat Commun (2019) 10(1):3573. doi: 10.1038/s41467-019-11405-4

  • CrossRef
  • Google Scholar

• 74

  WangYTiwaryP. Understanding the role of predictive time delay and biased propagator in rave. J Chem Phys (2020) 152(14):144102. doi: 10.1063/5.0004838

  • CrossRef
  • Google Scholar

• 75

  PedersenMFWrobelTMMarcher-RorstedEPedersenDSMollerTCGabrieleFet al. Biased agonism of clinically approved mu-opioid receptor agonists and Trv130 is not controlled by binding and signaling kinetics. Neuropharmacology (2020) 166:107718. doi: 10.1016/j.neuropharm.2019.107718

  • CrossRef
  • Google Scholar

• 76

  ZhuangYWangYHeBHeXZhouXEGuoSet al. Molecular recognition of morphine and fentanyl by the human M-opioid receptor. Cell (2022) 185(23):4361–75.e19. doi: 10.1016/j.cell.2022.09.041

  • CrossRef
  • Google Scholar

• 77

  RicarteADaltonJARGiraldoJ. Structural assessment of agonist efficacy in the mu-opioid receptor: Morphine and fentanyl elicit different activation patterns. J Chem Inf Model (2021) 61(3):1251–74. doi: 10.1021/acs.jcim.0c00890

  • CrossRef
  • Google Scholar