Characterizing Interhelical Interactions of G-Protein Coupled Receptors with the Fragment Molecular Orbital Method.

G-protein coupled receptors (GPCRs) are the largest superfamily of membrane proteins, regulating almost every aspect of cellular activity and serving as key targets for drug discovery. We have identified an accurate and reliable computational method to characterize the strength and chemical nature of the interhelical interactions between the residues of transmembrane (TM) domains during different receptor activation states, something that cannot be characterized solely by visual inspection of structural information. Using the fragment molecular orbital (FMO) quantum mechanics method to analyze 35 crystal structures representing different branches of the class A GPCR family, we have identified 69 topologically equivalent TM residues that form a consensus network of 51 inter-TM interactions, providing novel results that are consistent with and help to rationalize experimental data. This discovery establishes a comprehensive picture of how defined molecular forces govern specific interhelical interactions which, in turn, support the structural stability, ligand binding, and activation of GPCRs.

Introduction

G-protein coupled receptors (GPCRs) have enormous physiological and biomedical importance and are involved in a wide range of diseases. It is, therefore, not surprising that 475 drugs (∼34% of all drugs approved by the U.S. Food and Drug Administration (FDA)) act on this protein family.[1] However, while the human genome contains over 800 GPCR-encoding genes, only 108 of these are targeted by currently approved therapeutics. GPCRs thus represent one of the most promising and important classes of current pharmacological targets.

The structure of a GPCR can be divided into three parts: (1) the extracellular region, consisting of the N-terminus and three extracellular loops (ECL1–ECL3); (2) the transmembrane domain, consisting of seven α-helices (TM1–TM7); and (3) the intracellular region, consisting of three intracellular loops (ICL1–ICL3), an intracellular amphipathic helix (H8), and the C-terminus. The extracellular region often modulates ligand access; the TM domain forms the structural core, binds ligands, and transduces this information to the intracellular region through conformational changes; and the intracellular region interfaces with cytosolic signaling proteins.

It is becoming increasingly clear that the structural stability, function, and ligand binding properties of GPCRs are largely driven by the strength of interactions between different transmembranes (TMs).[2,3] Although there is evidence that the thermodynamic stability of GPCRs can be manipulated via mutation of specific TM residues,[4−7] the strength and chemical nature of the molecular forces responsible for “holding” together these seven TMs of the GPCR bundle and a molecular understanding of how these forces facilitate receptor activation and ligand binding remain to be elucidated.

In 2013, a structural analysis[3] of the 20 GPCR crystal structures available at the time revealed a consensus network of 24 potential inter-TM interactions (defined as contacts) arising from the close proximity between 36 amino acids. The importance of 14 out of these 36 amino acids to the structural stability and activation of GPCRs was validated by previously published site-directed mutagenesis studies, which had shown that mutations of these residues tend to affect receptor function, resulting in either an increase or a loss of receptor activity.[8] While earlier studies[2,3] have identified potential interactions (contacts, defined based on distance criteria) arising from the close distance between TM residues, the “actual” interactions (i.e., strength in kilocalories per mole and chemical nature such as hydrophobic, electrostatic, etc.) between GPCR residues have not been identified.

The aim of this study was to identify and to characterize the size and chemical nature of inter-TM interactions of a representative set of class A GPCRs. This information will improve our knowledge of the molecular mechanisms underpinning receptor stability and function and aid GPCR structural biology and structure-based drug design (SBDD). Our understanding of the molecular mechanisms underlying the different functional properties of GPCRs is highly dependent on the availability of high-resolution structural data.[9−11] However, even with crystal structures in hand, visual inspection and the force-field-based molecular mechanics (MM) calculations often used for structural exploration cannot explain the full complexity of intramolecular interactions.[12] Recently, several notable reports have been published[12−15] that emphasize the crucial role of “underappreciated” or nonobvious intramolecular interactions involved in biomolecular recognition. These interactions include CH/π,[16,17] halogen/π,[18] cation/π,[19] and nonclassical hydrogen bonds,[20] which are often not properly parametrized in currently available force fields (FFs).[14] Furthermore, the role of hydrophobic interactions, vital for receptor stability,[21] still has no reliable predictive method for its quantification aside from quantum mechanical (QM) ones.[4,12]

Quantum mechanical methods have always been considered to be a reliable approach for the exploration of molecular interactions.[22,23] However, despite their many advantages, traditional QM approaches are generally not feasible for large biological systems such as GPCRs, due to their high computational cost.[24] We have therefore employed the fragment molecular orbital (FMO) quantum mechanical approach[8,17,23,25] in the current study. FMO offers a considerable computational speedup over traditional QM methods[26] and is an extensively validated method for the structural exploration of large biological systems.[8,17,27,28] A second key advantage of FMO is that it provides a quantitative breakdown of the interactions formed between pairs of fragments (residues), including their strength (in kilocalories per mole) and chemical nature (electrostatic or hydrophobic).[17] FMO offers an excellent solution that combines accuracy, speed, and the ability to reveal key interactions that would otherwise be hard to detect.[8]

The accuracy and speed of FMO are achieved by dividing the system into smaller pieces called fragments (Figure ). For example, each residue within a GPCR protein can be represented by a fragment. By performing QM calculations on fragments, one can make the computational cost scale almost linearly with respect to the system size. The pair interaction energy (PIE) between any two fragments calculated by FMO is a sum of four energy terms: electrostatics, exchange repulsion, charge transfer, and dispersion, and is provided by pair interaction energy decomposition analysis (PIEDA; Figure ).[29] The electrostatic and charge transfer components are important in salt bridges, hydrogen bonds, and polar interactions, while dispersion is more hydrophobic in nature. The exchange repulsion term describes the steric repulsion between electrons[24] that prevents atoms from collapsing into each other. We used FMO to calculate the pair-attraction energy (PAE; see eq and Figure ) between each residue pair (fragments i and j, see Methods) in the GPCR. The total PAE (TAE) calculated by FMO is a sum of individual PAEs and describes the overall attraction energy for each TM–TM pair.

Figure 1

Illustration of GPCR fragment generation and details of each of the four PIE components being computed using pair interaction energy decomposition analysis (PIEDA). The electrostatic term arises from the Coulomb interaction between polarized charge distributions of the fragments. The exchange repulsion term is derived from the interaction between fragments situated in close proximity and is always repulsive; it is due to Pauli repulsion and is related to the overlap of two occupied orbitals. The charge transfer term arises from the interaction between occupied orbitals of a donor and unoccupied orbitals of an acceptor. The dispersion term arises as a result of the interaction between instantaneous dipole moments of two fragments; it is hydrophobic (nonpolar) in nature and is obtained in PIEDA from the correlation energy of the electrons.

FMO is an experimentally validated method as reported in the literature,[8,17,27,28,30−36] and its usefulness has been proven in numerous drug design cases, including lead optimization of novel ITK inhibitors,[31] in binding studies of RNA-protein-translation inhibitors,[37] in the discovery of novel Hsp90 inhibitors by fragment linking,[38−40] in the discovery of novel natural products for prion disease,[41] and in many other examples.[42,43] In this work, we extended the use of FMO to interrogate molecular interactions within GPCRs and how these relate to receptor function. The principal difference between FMO and MM/FF methods arises because FMO takes into account polarization (in the self-consistent mutual polarization of fragments) and charge transfer (whereby charge is allowed to flow between fragments).[17,44] The description of electrostatics in most popular force fields is based on static charges, which neglects polarization and, in polar systems such as proteins, provides an approximation of the actual electronic state. The van der Waals forces, although perhaps reasonably well parametrized on average, fail to capture the directional nature of the dispersion terms involving halogens.[45] These theoretical considerations explain why PIE values calculated with FMO often correlate well with experimental values[8,27,28,30−32] and why FMO clearly outperformed MM methods in cases where PIE values were obtained with MMFF94x (with AM1-BCC charges on the ligand atoms)[34] and with MM/GBSA using the Amber12EHT force fields.[31] In performing detailed analyses of biological systems with the level of computational power currently available, there is no need to compromise by restricting the methodology to MM calculations when a similar and more powerful analysis can be done with FMO on a similar time frame.[32]

Methods

Test Set

In this work, we applied FMO to characterize the strength and chemical nature of the inter-TM interactions of a set of 35 class A GPCR–ligand crystal structures that represent different branches of the GPCR tree. Among this set there are six receptors with structures representing both the active and inactive states. We used the following criteria in selecting this set: (1) Selection was limited to the class A GPCR crystal structures that were available in the Protein Data Bank (PDB) at the time that this research was started. (2) We removed all structures with a resolution of >3.5 Å as not suitable for FMO. (3) We selected the highest resolution structure of each receptor as a representative (two representatives were chosen in cases where both active and inactive state structures existed). This was done to prevent biasing of the test set toward receptors with many published structures.

Residue Numbering and Structure Preparation

The position of the amino acid residues within each GPCR is described by the general numbering scheme proposed by Ballesteros and Weinstein,[46,47] a scheme for class A GPCRs whereby X.50 represents the most conserved residue (canonical residue) on helix X. The canonical residues of class A GPCRs are (with percentage of conservation): N1.50, 98%; D2.50, 90%; R3.50, 95%; W4.50, 97%; P5.50, 78%; P6.50, 99%; P7.50, 88%.[48] The remaining residues in each helix are numbered sequentially from the appropriate canonical residue; numbers decrease toward the N-terminus and increase toward the C-terminus.

During structure preparation, hydrogen atoms were added to the crystal structures at physiological pH (7.0) with the Protonate3D[49] tool implemented in MOE version 2016.08 (Chemical Computing Group), which assigns ionization states and positions hydrogens in proteins, ligands, and solvents for a given set of three-dimensional coordinates. The uncertainty associated with the position of individual atoms within a given crystal structure is dependent on the B-factor of the atom and the overall resolution of the structure. As small errors in the positions of atoms can translate to large deviations in energy terms, it is important to optimize individual crystal structures before applying any type of calculation to them.[11] In the present study, we applied a constrained minimization procedure with the semiempirical AMBER10:EHT force field[50,51] implemented in MOE version 2016.08, which allowed each atom to deviate by up to 0.5 Å from its original position in the crystal structure.

FMO Calculation Protocol

The FMO approach is a general quantum mechanical method used to understand the electronic states of the specific molecular interactions that take place within large molecules and macromolecular complexes. It is particularly useful for understanding interactions between residues within proteins and between proteins and their ligands. FMO calculations can be applied to any set of atoms within a given protein, whether the protein is soluble or membrane-bound, and can also be applied to ligands.

In FMO, the system is fragmented and the pair interaction energy (PIE; see eq ) between every fragment pair is calculated using pair interaction energy decomposition analysis (PIEDA).[17,23] For example, in proteins, each residue can be represented by a fragment. By performing QM calculations on fragments, one can achieve high computational efficiency, often resulting in linear scaling as a function of system size. The FMO method has been efficiently parallelized for CPU clusters,[26] making its calculations rapid and relatively inexpensive. One FMO calculation for a full-sized receptor took only 2 h on 340 CPU cores. A detailed description of the fragmentation strategy and the basic methodology underpinning FMO, including detailed mathematical formulations, is beyond the scope of this article, but can be found in a number of reviews.[17,23,29]

Here, the FMO method[17] was applied to GPCRs using FMO code version 5.1,[26] which is embedded in the general ab initio quantum chemistry package GAMESS (General Atomic and Molecular Electronic Structure System).[52] We used a well-established FMO protocol[16,17,34,53,54] to characterize the strength (in kilocalories per mole) and chemical nature (electrostatic or hydrophobic) of the inter-TM interactions. Our calculations were performed using the MP2 method (second order Møller–Plesset perturbation theory[55]) with the 6-31G* basis set.

The FMO calculations consisted of the following four key steps: (a) fragmentation (i.e., assigning atoms in a system to specific fragments); (b) fragment self-consistent field (SCF) calculations in the embedding polarizable potential, so that fragments mutually polarize each other in a self-consistent fashion to account for intrafragment charge transfer and other quantum effects; (c) fragment pair SCF calculations, to permit inclusion of interfragment charge transfer; and (d) total property (energy, gradient, etc.) evaluation.

As shown in eq , the PIE (Ein) between fragments i and j is a sum of four energy terms: electrostatics (Ees), exchange repulsion (Eex), charge transfer (Ect), and dispersion (Edi).

In our investigations of the interactions formed between residues of different TMs, fragment i refers to a residue of TM and j refers to a residue of TM. The pair-attraction energy (PAE or Eattr) between fragment i and fragment j is the sum of the electrostatic, charge transfer, and dispersion energy terms as shown in eq . In this work, we used PAE instead of PIE due to the fact that the GPCR crystal structures had relatively low resolutions (average 2.7 Å). Making use of the exchange repulsion term, which requires high quality structures, would have been potentially misleading.

We used eq to calculate the contribution of the electrostatic terms to the overall attraction energy. The range of f is between 0 and 1, where values approach 1 when the interaction is purely electrostatic and approach 0 when it is purely dispersive.

Results

In this work, we applied the FMO method (for more details, see Methods) to characterize the strength and chemical nature of the inter-TM interactions of a set of 35 class A GPCR–ligand crystal structures that represent different branches of the GPCR tree (Supporting Information, Table 1). The rationale for selecting this specific set of receptors is described in Methods. Among this set there are six receptors with structures representing both the active and inactive states. Topologically equivalent positions of residues were identified using the Ballesteros–Weinstein numbering scheme (BW; see section ). The BW numbering scheme was used because it allowed us to compare our current analysis with previously published reports. Based on previous reports,[17] we considered any interaction with an absolute PAE ≥ 3.0 kcal/mol to be significant.[8] The interaction is considered to be conserved if it was identified in a majority (i.e., ≥2/3 (≥65%)) of the systems studied.

Consensus Network of Inter-TM Interactions and Ligand Binding

Our FMO analysis reveals a consensus network of 51 inter-TM interactions that are mediated by 69 topologically equivalent amino acids (Figure b). On average, the level of conservation of these interactions across the 35 tested systems is 82%. According to published experimental site-directed mutagenesis (SDM) data extracted from GPCRDB[56] (Supporting Information, Figure S1), mutation of 30 of these 69 residues results in either an increase or a loss of receptor activity with ≥5-fold effect (for the remaining 39 residues, either there are no SDM data reported or the change in receptor activity was within 5-fold).

Figure 2

(a) Representative β2 adrenergic receptor (ribbons)–ligand (spheres) complex (PDB code 2RH1). The conserved inter-TM interactions are shown as white tubes. (b) Network of 51 conserved inter-TM interactions formed by 69 residues. The circles represent residues and are color-coded as follows: TM1, red; TM2, brown; TM3, yellow; TM4, gray; TM5, teal; TM6, light blue; and TM7, dark blue. Numbers denote Ballesteros–Weinstein numbering. A dashed line between a pair of circles indicates the presence of a conserved interaction. Residues previously reported[8] as involved in ligand binding in a number of different GPCRs are marked with a red triangle. (c) Schematic representation of the TM–TM interaction energies. The line between a pair of circles indicates the total TM–TM pair attraction energy (TAE, in kilocalories per mole), where the thickness of the line is proportional to the size of the TAE (only interactions < −20 kcal/mol are shown). (d–f) Three examples of conserved inter-TM interactions in a representative GPCR (the β2-adrenergic receptor). Nitrogen atoms are shown in blue, oxygen atoms are shown in red, sulfur atoms are shown in yellow, and carbon atoms are shown in green. Major contributions to residue–residue interactions are highlighted with yellow dashed lines.

As illustrated in Figure a, ligands “sit” on top of the “pile” of conserved inter-TM interactions, thereby interacting directly with neighboring residues as well as indirectly with the entire inter-TM network. This structural arrangement helps to explain why mutation of residues that are located at a significant distance from the ligand binding site can have such a strong effect on ligand binding, numerous examples of which can be observed from SDM results (Figure S1). As we previously reported,[8] the residues at positions 3.32, 3.33, 6.48, 6.51, 6.52, 7.39, and 7.43 make substantial contributions to receptor–ligand binding in >70% of all analyzed structures. We show here that six of these seven ligand-binding residues (except those at position 3.33) are also involved in the conserved network of 51 inter-TM interactions. For example, residue 6.52 simultaneously interacts with the ligand and with residue 5.47 (Figure e). SDM studies have shown that mutations in these positions frequently affect ligand binding affinity and selectivity.[56,57] This overlap between residues involved in ligand binding and those involved in inter-TM interactions can provide an explanation for how a ligand that binds at the extracellular end of the receptor is able to exert an effect on the overall structure of the GPCR.

FMO detected that some residues form more than one conserved inter-TM interaction with a neighboring residue and intriguingly this includes three canonical residues: N1.50, D2.50, and W4.50 (see section , Ballesteros–Weinstein numbering scheme for the definition of canonical residues). For example, TM2 residue D2.50 forms hydrogen bonds with the TM7 residue located at position 7.46 (Figure d); this hydrogen bond appears in 69% of the analyzed structures. Residue N1.50 interacts with both TM2 and TM5 (Figure f), and residue W4.50 frequently interacts with both TM2 and TM3 (Supporting Information, Table 2). The considerable contribution of canonical residues N1.50, D2.50, and W4.50 to the conserved inter-TM interaction network can help to rationalize why these residues are so conserved among class A GPCRs. By contrast, four other canonical residues—R3.50, P5.50, P6.50, and P7.50—form inter-TM interactions only in specific cases, with no clear segregation between the active and inactive forms of receptors.

Chemical Nature of the Conserved Inter-TM Interactions

As illustrated in Figure , out of the 51 conserved inter-TM interactions, 15 were predominantly electrostatic in nature (f > 0.5), 15 were predominantly hydrophobic (f < 0.5), and the remaining 21 had a mixed chemical nature (f ≈ 0.5, indicating equal electrostatic and hydrophobic contributions). Our findings thus emphasize the pivotal role of hydrophobic forces in TM–TM interactions, something that is quite often omitted from structure-based descriptions.[17] These findings also demonstrate that not only are these specific amino acid interactions conserved across many different GPCRs, but their strength and chemical character (hydrophobic or electrostatic) is conserved as well. This is despite the fact that the individual amino acids found at each of these 69 positions have an average sequence similarity of just 46% across the 35 receptors analyzed (Figure S2). It is an intriguing observation that the conservation of these positions during evolution appears to have taken place at the level of the interaction and not at the level of the specific amino acid residues forming the interaction.

Figure 3

Chemical character of the conserved inter-TM interactions calculated with PIEDA (Supporting Information, Table 3). Boxes are colored according to their f (chemical) factor: from dark blue (100% dispersion contribution) to yellow (100% electrostatic). The absence of a contact is represented by a white box. The bottom line (“Average”) represents the average f chemical factor of each inter-TM interaction and is color-coded using the same scheme as the matrix. The matrix is sorted by f chemical factor.

Role of Specific TM Helices

FMO identified the central role of TM3 in receptor stability and function, as it interacts strongly with five neighboring TMs (TM2, TM4, TM5, TM6, and TM7). TM3 forms 39% (20 out of 51) of all conserved inter-TM interactions (Figure b), and a specific TM3 residue, at position 3.32, is also frequently involved in receptor–ligand binding.[8]

We also observed that TM3, TM6, and TM7 form substantially weaker TM–TM pair attraction energies (TAE of −82 kcal/mol on average) with each other compared to the TAEs formed with and between other TMs (−114 kcal/mol on average). This correlates with the observation that the ligand binding site of class A GPCRs is frequently located between these three TMs.[58] Weaker TAEs between TMs suggest an area of the receptor that is more “malleable” and can adopt different conformations to accommodate ligands of different sizes and shapes. This observation also explains how GPCRs can bind very diverse ligands within the same binding site.[28] An additional intriguing point is that the TM6 residues involved in the conserved inter-TM network are located exclusively in the extracellular half of this helix (between positions 6.40 and 6.54). This is consistent with the fact that the intracellular half of TM6 is more dynamic and can “open” during the activation process.[59]

Comparing Active and Inactive Protein States

Among the 35 tested receptors, six had crystallographic coordinates for both active and inactive states (Supporting Information, Table 1). This allowed us to compare the inter-TM interactions present in the two states for these six GPCRs. It is known that during the activation process receptors undergo conformational changes, as shown in Figure a.[58] However, despite these significant structural rearrangements, 46 of the conserved 51 inter-TM interactions remained unchanged (Figure b). FMO also detected 13 state-specific, conserved interactions (six for the inactive and seven for the active state, Figure c).

Figure 4

Comparison of inter-TM interactions in inactive and active states for the six proteins that have published crystal structures for both states (PDB codes for the inactive and active structures, respectively, are rhodopsin, 1GZM and 3PQR; β1-adrenergic receptor, 4BVN and 2Y02; β2-adrenergic receptor, 2RH1 and 4LDE; M2 muscarinic receptor, 3UON and 4MQS; μ-opioid receptor, 4DKL and 5C1M; A2A adenosine receptor, 5IU4 and 4UHR). (a) Inactive (orange ribbon) and active (green ribbon) structures of the M2 muscarinic receptor are superimposed (PDB codes 3UON and 4MQS, respectively). (b) Overlap in terms of conserved inter-TM interactions between inactive and active states shown using a Venn diagram. (c) Comparison between state-specific, conserved inter-TM interactions. In the matrix, the size of the PAE between residues is shown as a heat map colored according to the gradient on the right. The absence of an interaction is shown as a gray box. (d–i) Examples of conserved changes in the inter-TM interaction network as a result of receptor activation. Nitrogen atoms are shown in blue, oxygen atoms are shown in red, sulfur atoms are shown in yellow, and carbon atoms are shown in green (active state) or in orange (inactive state). (d–f) M2 muscarinic receptor; (g–i) β2-adrenergic receptor.

It was found that the loss of an inter-TM interaction from the inactive state often coincides with the formation of a different interaction in the active state. These changes frequently involve the switch of a residue from one interacting partner to another, for example from 3.46–6.37 (inactive) to 3.46–7.53 (active) (Figure d). The changes can also happen in separate locations; for example, the loss of two hydrophobic interactions (3.43–5.54 and 3.43–6.41) in the inactive state (Figure e) coincides with the formation of a new hydrogen bond between positions 2.50 and 3.39 in the active state (Figure f). As an additional example, the loss of a nonclassical hydrogen bond 6.40–7.49 (inactive) (Figure g) coincides with the formation of a new CH−π interaction 3.40–6.44 (active) (Figure h) and a new face-to-face π-stack 5.47–6.44 (active) (Figure i). The FMO calculations support previously reported[2] SDM data showing that mutations of residues in positions 3.46, 6.37, and 7.53 resulted in a markedly reduced ability of the receptor to activate.[60] Maintaining almost the same number of interactions when moving from the inactive to the active state does not mean that this change is energetically neutral overall. The average difference between the TAEs of the inactive and active inter-TM networks is about 30 kcal/mol, with the active state being less stable. This energy difference is often at least partly compensated by the agonist binding.[8,28]

Underappreciated Interactions

A number of different types of interactions, such as classical hydrogen bonds and salt bridges, can be easily identified by visual inspection. However, there are a variety of additional interactions that also play vital roles in residue–residue binding that are not so straightforward to detect and that are not properly parametrized in many currently available force fields.[12−15] We found that approximately 50% of the 51 conserved interactions are “underappreciated” in the sense that these are not interactions that could reliably be detected by non-QM methods. Many of these underappreciated interactions are formed by backbone carbonyls and include nonclassical hydrogen bonds (see examples in Figure a,b,e,h), CH−π interactions (Figure c,d,h), and carbonyl–S interactions (Figure g). Many residues involved in underappreciated interactions are located at the junctions between TMs, allowing them to mediate the interaction networks between these TMs. For example, the residue at position 3.40 forms mediating interactions between TM5 and TM6 (see example in Figure h) by forming a CH−π interaction with F6.44 and a nonclassical hydrogen bond with the backbone carbonyl of S5.46. The role of hydrophobic interactions is also vital for biomolecular recognition, but there is still no reliable non-QM method for its quantification.[12] Therefore, we would also consider hydrophobic interactions to be underappreciated in this sense. The FMO calculations indicated that 15 of the 51 conserved inter-TM interactions were predominantly hydrophobic in nature.

Figure 5

Examples of “underappreciated” interactions. Nitrogen atoms are shown in blue, oxygen atoms are shown in red, sulfur atoms are shown in yellow, and carbon atoms are shown in green. (a–g) Active state of the β2-adrenergic receptor (PDB code 4LDE). (a) Nonclassical hydrogen bond between the side chain of V1.43 and the backbone carbonyl of G2.54. (b) Nonclassical hydrogen bond between the side chain of I1.57 and the backbone carbonyl of N2.40. These two residues also form an additional hydrophobic interaction. (c) CH−π interaction between S3.30 and F4.58. (d) CH−π interaction between S3.30 and F4.58. (e) Side chain–side chain nonclassical hydrogen bond between V2.38 and D3.49. (f) Two nonclassical hydrogen bonds formed between I3.4 and S5.46. (g) Carbonyl (backbone)–S interaction between I7.47 and C6.47. (h) Dopamine D3 receptor (PDB code 3PBL): I3.40 forms two nonclassical interactions with F6.44 (CH−π interaction) and with S5.46 (nonclassical hydrogen bond with the backbone carbonyl).

Conclusions

GPCRs regulate almost every aspect of cellular activity, making them key targets for drug discovery. It has been shown that the key properties of GPCRs are largely driven by interactions between residues on different helices of the receptor. However, the strength and chemical nature of these interactions have not previously been described. In the present study, we have used a quantum mechanical method (FMO) to characterize the molecular forces responsible for defining the positioning of the seven helices of the GPCR bundle. Our work reveals novel findings that are consistent with experimental data and provides a comprehensive picture of how interhelical interactions support the structural stability, ligand binding, and activation of GPCRs.

The FMO methodology has allowed us to discover and characterize the strength and chemical nature of 51 conserved inter-TM interactions formed by 69 residues. These are novel observations, but they are consistent with the experimental data. The strength and chemical nature of these interactions are conserved among all of the class A GPCRs analyzed. Six of these 69 residues are also frequently involved in ligand binding,[8] which helps to explain how residues involved in ligand binding can have an effect on the overall structure of the receptors. We have also found that these inter-TM interactions are conserved in their strength and chemical nature despite the fact that the participating amino acid residues are not conserved from one GPCR to another. FMO indicated that 15 of the 51 conserved inter-TM interactions were predominantly hydrophobic in nature. We also highlighted the central role of TM3 in receptor stability, function, and ligand binding, having observed that TM3 is involved in 39% of all conserved inter-TM interactions. We hypothesize that lower TAEs between TMs, such as those seen with TM3, indicate areas that are more “malleable” and able to adopt several conformations. This enables members of the GPCR superfamily to accommodate molecules of different sizes and shapes.

In 2013, the Babu’s group[3] identified 24 potential inter-TM interactions. We have determined that only 15 (Supporting Information, Table 2) of these passed the FMO cutoff (PAE ≤ −3.0 kcal/mol) with the remaining nine failing to make the cutoff. FMO also indicated that a core of 46 of the conserved 51 inter-TM interactions remained unchanged between the inactive and active states of the receptor. Where state-specific interactions were detected, the loss of one inter-TM interaction from the inactive state often coincided with the formation of a different interaction in the active state. However, moving from the inactive to the active state is not energetically neutral overall, with the active state being less stable. This energy gap is often compensated by agonist binding.[8,28]

In this study, we found that almost half of the conserved interactions are of a type that is often underappreciated and hard to detect with non-QM methods. The use of FMO therefore provides a more holistic means of identifying key molecular interactions involved in receptor structure and function. Our study also provides information on the internal energy balance required for stability, dynamics, and ligand binding of GPCR receptors. These structural insights can be applied to the design of ligands that can more efficiently interact with the inter-TM network, controlling receptor structure and flexibility and thereby affecting function. FMO is an additional resource to further our understanding of GPCR function at the atomic level, and the regular application of FMO to GPCR studies may lead to the generation of more effective GPCR-targeted drugs.[1,61] This approach can be applied to the structural exploration of other protein superfamilies and biosystems.

It would be of particular interest to apply FMO to explore the effect of water molecules on GPCR structure and activation. Recently it was reported[62] that computer simulations of diverse GPCR crystal structures revealed the presence of a conserved water molecule network. This network is maintained across the inactive and active states. As suggested by the authors, these conserved water-mediated interactions near the G-protein-coupling region, along with diverse water-mediated interactions with extracellular ligands, have direct implications for structure-based drug design and GPCR engineering. Applying FMO could provide important insights into the role of these water molecules and their interactions with the inter-TM network and with ligand binding. An interesting additional point for future exploration would be to determine how different ligands affect the inter-TM network of the same receptor and how this can lead to activation.[28]

An extension of this research to the exploration of the interactions formed between residues within the same TM (intra-TM interactions) will complete our FMO-informed understanding of key molecular interactions that influence GPCR structure and function. For example, it has been observed that ionic interaction of the conserved TM3 amino acids R3.50 and D3.49 (the E/DRY motif) maintains the receptor in its ground state.[63] This hypothesis has been confirmed by visual inspection of the rhodopsin ground-state crystal structure and by computational modeling approaches.[63] However, there are two groups of receptors within class A GPCRs that make very different uses of the E/DRY motif. In first group, nonconservative mutations of the glutamic acid/aspartic acid–arginine residues lead to an increased affinity for agonist binding, retain G protein coupling, and retain an agonist-induced response. In contrast, however, in the second group the E/DRY motif is more directly involved in governing receptor conformation and G protein coupling/recognition.[63] This example provides further evidence of residue–residue interactions that have a direct effect on ligand binding despite taking place some distance from the binding site. It is important, therefore, to look beyond the rhodopsin ground-state model of conformational activation to clarify the role of this highly conserved TM3 triplet in GPCR activation and function. The application of FMO to intrahelical interactions will be extremely useful for this purpose.

The application of FMO to the study of interhelical interactions has highlighted the utility of this methodology for interrogating features important to GPCR function, and we are now in a position to apply FMO to GPCR features such as ion binding sites. Our understanding of ion binding sites and their effects on GPCR structure and function has greatly expanded in the past few years;[64,65] however, we are only now beginning to understand the potential for applying this information to the discovery of more efficient and safer drugs with improved subtype and/or functional selectivity. The sodium site stands out as highly conserved among most class A GPCRs,[65] binding the ion in the middle of the 7TM helical bundle anchored at the most conserved aspartate residue D2.50. The analysis[65] of 45 diverse class A GPCRs has revealed the highly conserved nature of the sodium pocket, 15 residues of which are conserved with very minor variations. Moreover, those class A receptors that lack these key residues of the sodium binding site naturally, or via introduced mutations, have their ligand-induced signaling dramatically reduced or completely abolished.[65,66] It will be of especial interest to apply FMO to better understand the interplay between the sodium ion, the sodium binding residues within the GPCR, and the rest of the GPCR residues involved in the inter-TM interaction network to gain a more comprehensive appreciation of their mutual effects and the role they play in ligand binding.