Peer review process
Consolidated peer review report (28 November 2022)
GENERAL ASSESSMENT
The objectives of the study:
This paper aims to characterize the dynamics that drive allostery of the adenosine A1 receptor (A1R) via computational analysis of its activation free energy landscape and measurements of the appropriate geometrical parameters. This is done by focusing on the allosteric signaling pathways in different activation states, from inactive to active states via intermediate and pre-active ones, as well as the characterization of putative drug-binding pockets. The long-term objectives are to eventually be able to aid drug discovery efforts for this therapeutically important GPCR.
Key findings and major conclusions:
Conventional MD does not enable the sampling of the complete conformational landscape of receptor activation. Instead, enhanced sampling MD simulations are required to achieve this. Using metadynamics, the authors decipher the activation pathway of A1R, decode the allosteric networks and identify transient pockets. The protein energy networks computed throughout the inactive, intermediate active, pre-active and active conformational states unravel the extra and intracellular allosteric centers and the communication pathways that couple them, whereby the pathways are reinforced in the activated state. These conformations primarily differ in the dynamics of the ionic lock motif that couples TM3 to TM6 in the inactive conformation and reveal that G-proteins are required to fully stabilize the active conformation. Support for these findings comes from prior mutagenesis work on the A1R that identified key allosteric residues that in many cases map to identified communication nodes. Finally, the authors identified allosteric pockets throughout the A1R in four different conformational states that support prior experimental and MD studies on the mechanism of the positive allosteric modulator MIPS521 and which could be targeted for the design of modulators. Overall, these findings provide complementary support to a structure-based mechanism of activation and allosteric modulation of A1R, and extend the findings to incorporate dynamics across the full activation pathway.
The perceived strengths and weaknesses:
This preprint employs a combination of computational techniques to successfully reconstruct and analyze the conformational ensemble of the A1R activation. The metadynamics simulations supported the aim of the study, the results are clearly presented and the work is very well written. The authors could improve the discussion of how the protein energy network analysis could further advance rational design of specific modulators with a desired mode of action. The computational approach needs to be refined to be robust, with a focus on characterizing the convergence of the free energy landscapes. Overall, A1R is a good choice as the target for this study as there is existing structural and pharmacological data to support preliminary findings. Moreover, the framework presented herein could be adapted and scaled to other GPCRs with structural templates, which might enable comparison of allosteric pathways across families and classes.
RECOMMENDATIONS
Revisions essential for endorsement:
The paper could better demonstrate how the insights gained herein will or could lead to progress in the rational design of specific modulators with a desired effect. The authors should outline and discuss how they envision the modeling pipeline they have designed will be used towards this goal and tone-down or explain why “this information is essential to ease the design of allosteric modulators for A1R.”. A recent study on FFAR1, where the authors targeted a specific dynamic pocket could be helpful in this respect (https://www.pnas.org/doi/full/10.1073/pnas.1811066116). Specifically this might entail: How does specificity for a receptor correlate with pockets forming in a specific state? From this, how does one design an agonist vs. an antagonist vs. an inverse agonist? Does breaking a specific network select a function of the drug? How would another group follow up on this work, for example in a virtual screening campaign?
Free energy calculations:
a. A proof of convergence of the free energy calculations is missing. The authors argue that obtaining landscapes that do not change over time is proof of convergence, but this is incorrect in well-tempered metadynamics. The fact that the heights of the Gaussians decrease over time guarantees that the landscape will be stable over time, and the way to check convergence is to show that the collective variables become diffusive after convergence. In addition, to validate that the choice of collective variables (CV) is actually appropriate, they should check that CVs that were not biased are also diffusive. This would be best studied by looking at the behavior of microswitches that were not considered, such as ones describing the PIF motif, the NPxxY motif, the ligand binding pose, etc.
b. The authors should characterize the uncertainties/statistical errors on the measured free energy profiles to better evaluate the significance of change (e.g. for inspiration: https://www.plumed.org/doc-v2.7/user-doc/html/masterclass-21-2.html).
c. In the cMD trajectories, a large part of phase space is sampled, which does not appear consistent with what one would expect based on the free energy landscapes. For instance, it does not seem reasonable to cover an almost complete conformational transition in 500ns when the barrier of the system is on the order of 5-8kcal/mol. The definition of CVs may have led to an overestimation of the free energy barrier. Hence an independent validation of the free energy barrier height is needed, by e.g. changing the CV definition.
Configurations extracted from both conventional MD and wt-metadynamics are mixed in the analyses of the allosteric networks and the pockets. A more accurate way to integrate these datasets would be to modulate the weights of the configurations by their statistical weights, which can be retrieved from the metadynamics simulations.
Related to Figure S6, it is essential to compare the dynamics for all of the key class A activation motifs including the Na binding site, PIF motif, and NPxxY.
Please provide clarification on why 500 ns was chosen as the time-scale of the MD simulations and inclusion of the time course for the three independent MD simulations for each of the key structural features (e.g. TM6 torsion angles and TM3-TM6 distances).
The validation of the results in the form of previously published mutagenesis results does not appear completely convincing. Large parts of the protein are included in the allosteric network, making it likely that mutations in some of these residues will have an effect if mutated. In addition, the fact that mutations in ECL2 and ECL3 affect allostery is expected and does not constitute a good validation of the results. If no other results are included, we recommend that the language be toned down so as not to overstate the significance of the results.
What is the justification for using an energy-based scoring for network analysis, given that a conventionally correlation-based approach has been used successfully in the field? The concern with an energy-based approach is that the interaction energy calculations do not consider the dielectric effect, i.e., if water molecules interfere with two interacting residues. Since the dynamic network is one of the critical aspects of this study, we believe the authors need to explore other tools such as the one implemented in VMD (https://www.ks.uiuc.edu/Research/vmd/plugins/networkview/) and compare the results.
Provide generic residue numbers such as GPCRdb or Ballesteros Weinstein numbering for all mentioned residues in text and figures, as is standard for structural papers.
Additional suggestions for the authors to consider:
For the PEN analysis it would be useful to digest these communication networks with respect to the established structural activation motifs of class A GPCRs (Na binding site, PIF, and NPxxY) that are present at the A1R.
It is unclear why the authors chose two largely correlated CVs (See comment 2c). In addition, the choice of CV is likely contributing to the distortion of S6, as displayed in Figure 1E. It has been shown that choosing a different CV set that describes the motion between states in a more distributed way is more likely to lead to a converged conformational ensemble. We suggest repeating the metadynamics simulations with a more distributed CV set that encompasses all of the microswitches in the receptor.
To support the vision on how the analysis of activation pathway, energy networks and transient pockets could be used “to ease the design of allosteric modulators for A1R” (last sentence of the abstract), it might be necessary to show that the combination of these methods can indeed be predictive for the binding and effect of known ligands. This might provide a first step towards establishing that molecules that bind to pockets “near allosteric networks” is a promising avenue for drug discovery.
The specific TM3-TM6 residues should be specified in figures and text. Commonly used TM3-TM6 comparisons include the measured maximum distance between 2x46 to 6x37, which could be used here also (e.g. see https://docs.gpcrdb.org/structures.html#structure-descriptors).
Even though the "A1R in complex with PSB36 (PDB 5N2S)" is an inactive structure, PSB36 is an agonist. Hence, the authors should consider using the DU172 antagonist-bound structure for comparison (PPDB 5UEN)
How does adenosine and MIPS521 binding impact the different conformational states and PEN.
It would be interesting to note how the findings from this study compare/contrast to a very recently published report by Li et al, PNAS, 2022 “The full activation mechanism of the adenosine A1 receptor revealed by GaMD and Su-GaMD simulations”. Similarly with regards to the determination of allosteric binding pockets in this recent publication: “The pocketome of G-protein-coupled receptors reveals previously untargeted allosteric sites” (https://doi.org/10.1038/s41467-022-29609-6)
A major advantage of allosteric drugs is the potential to achieve higher selectivity. Expansion of this study to include other adenosine receptor subtypes or linking to other types of molecular pharmacology (e.g. biased signalling, subtype selectivity, etc.) would be a major benefit to the field.
Consider including an explanation of the physiological and pharmacological relevance of A1AR in the introduction.
Even if not entirely necessary for the results, it would be more consistent if the study would include metadynamics of the G-protein bound state.
Methods: "In other words, once the free energy surface does not change significantly during a relatively long period of time in the last part of the simulation". What is “relatively long period of time” and “change significantly”. The convergence, should be stated as a quantitative description of the observed energy differences.
The authors should strongly consider making their analysis code and simulation data publicly available (e.g. on GitHub or Zenodo) so that others can replicate and build upon this work
REVIEWING TEAM
Reviewed by:
Antonios Kolocouris, Professor, Department of Medicinal Chemistry Faculty of Pharmacy National and Kapodistrian University of Athens, Greece:
CADD and computational biophysics on adenosine receptors
David Thal, Senior Research Officer, Monash University, Australia:
structural biology and molecular pharmacology of allosteric mechanisms underlying Class A GPCRs
SciLifeLab Journal Club, Stockholm, Sweden (see Appendix for members)
Curated by:
Alexander S. Hauser, Associate Professor, University of Copenhagen, Denmark
APPENDIX
SciLifeLab Journal Club:
Feedback was generated in a meeting of the journal club involving:
Lucie Delemotte (Journal Club oversight), Associate Professor of Biophysics, KTH Royal Institute of Technology, Sweden: modeling and enhanced sampling of GPCRs and other membrane proteins.
Olivia Andén, PhD student, Stockholm University: cryo-EM and functional characterization of membrane proteins.
Cathrine Bergh, PhD student, KTH Royal Institute of Technology: enhanced sampling simulations of membrane proteins.
Koushik Choudhury, PhD student, KTH Royal Institute of Technology: membrane protein modeling, enhanced sampling simulations.
John Cowgill, postdoctoral scholar, Stockholm University: cryo-EM and simulations of membrane proteins.
Chen Fan, postdoctoral scholar, Stockholm University: cryo-EM and simulations of membrane proteins.
Nandan Haloi, postdoctoral scholar, KTH Royal Institute of Technology: membrane protein modeling, free energy calculations, structure refinement in cryo-EM maps.
Rebecca J Howard, researcher, Stockholm University: membrane protein structure-function, allosteric modulation.
Marie Lycksell, PhD student, Stockholm University: structure and simulations of membrane proteins.
Antoni Marciniak, PhD student, KTH Royal Institute of Technology: enhanced sampling simulations of GPCRs and other membrane proteins.
Darko Mitrovic, PhD student, KTH Royal Institute of Technology: membrane protein modeling, enhanced sampling, machine learning.
Alex Payne, PhD student, Memorial Sloan Kettering Center for Cancer Research: membrane protein modeling, cryo-EM structure determination, drug discovery.
Urška Rovšnik, PhD student, Stockholm University: cryo-EM and functional characterization of membrane proteins.
Akshay Sridhar, postdoctoral scholar, KTH Royal Institute of Technology: membrane protein modeling, enhanced sampling simulations.
Amanda Dyrholm Stange, PhD student, Aarhus University: membrane protein modeling, enhanced sampling simulations.
(This consolidated report is a result of peer review conducted by Biophysics Colab on version 3 of this preprint. Minor corrections and presentational issues have been omitted for brevity.)