Molecular dynamics of the histamine H3 membrane receptor reveals different mechanisms of GPCR signal transduction

In this work, we studied the mechanisms of classical activation and inactivation of signal transduction by the histamine H3 receptor, a 7-helix transmembrane bundle G-Protein Coupled Receptor through long-time-scale atomistic molecular dynamics simulations of the receptor embedded in a hydrated double layer of dipalmitoyl phosphatidyl choline, a zwitterionic polysaturated ordered lipid. Three systems were prepared: the apo receptor, representing the constitutively active receptor; and two holo-receptors—the receptor coupled to the antagonist/inverse agonist ciproxifan, representing the inactive state of the receptor, and the receptor coupled to the endogenous agonist histamine and representing the active state of the receptor. An extensive analysis of the simulation showed that the three states of H3R present significant structural and dynamical differences as well as a complex behavior given that the measured properties interact in multiple and interdependent ways. In addition, the simulations described an unexpected escape of histamine from the orthosteric binding site, in agreement with the experimental modest affinities and rapid off-rates of agonists.

General physical properties of the molecular system and their definitions -The root mean-square deviation (RMSD) of the cartesian coordinates of the Cα atoms measures the structural drift of a molecule. In the absence of determining free energy landscapes for the membrane-protein system, the time dependence of the RMSD is a good indicator of the convergence of the 3D structure of the protein along the simulation to a stable state. -The 2D distance or contact maps for the average structure represent the Cα atom distance between all possible amino acid residue pairs in a 3D protein structure. They can be used to describe similarity between protein structures and to represent characteristic patterns of secondary structure. The origin, i.e., residue 1 (N-ter) for the x-and y-axes is at the upper-left corner. The last residue (C-ter) is at the lower-right corner. Red zones denote proximities between atoms, and thus between helices.
-Principal Component Analysis (PCA) represents a classic dimension reduction approach by constructing orthogonal linear combinations of many properties (in this case conformations), called principal components (PC). The greatest variance of the data lies on the first component, the second greatest variance on the second component, and so on. The PCA allows the identification of conformational communities or clusters, dividing the conformations into several populations and covering in this way a well-defined region in diversity space given by the principal S3 components. These essential degrees of freedom describe major collective modes of fluctuation that are relevant for the function of the protein. 1 PCs may be referred here to as 'meta-or super-conformations' and are very useful in investigating the molecular motions of proteins. 2,3 Through the PCA, molecular dynamics simulations (MD) allow the sampling and the grouping of similar conformational states.
For the analysis of the data generated by the MD production trajectories, we used the following programs: -EUCB 4 for the computation of the following properties, with 10 ns jumps and a switch distance of 5 Å. The residues composing the internal cavity of the representative structures of the receptor produced by CARMA were defined by finding the "wet" residues, i.e. those residues that were in contact (3 Å or less) with an internal water for more than 80% of the trajectory. The "wet" atoms of the ligands were determined in the same fashion. This approach determined the residues that had an interaction with a water molecule that had itself an interaction with the ligand. On another hand, pocket detection, ligand binding site and analysis and visualization of tunnels and channels were performed with the Fpocket, 9 SiteHound 10 and Caver software tools. 11 Those residues found in at least three out of four searches were considered cavity residues. S5 We did not deem necessary to energy-minimize the complexes at the frames we looked at during the analysis. Instead, we selected ten frames at regular intervals of the MD trajectory to obtain the interactions between the ligand and the receptor, since we compiled the longlasting neighbors of the former all along the trajectory.
When dealing with the lipid components of the membrane, we wished to find out whether there were any DPPC-specific binding sites on the receptor. For this purpose, we found out the amino acids that are in contact at least 80% of the time with at least one atom of a lipid molecule defining high-affinity, specific binding sites with the membrane protein.
We performed the computations in a Linux cluster with one master node of 8 CPUs, 10 To stocking; and 192 CPUs (16 nodes of 12 Intel XEON E5630@2.53GHz, 24 GB RAM, CPUs each).
The NAMD parallel jobs were executed through the mpiexec application. We also ran calculations in the massively parallel IBM Blue Gene of the HPC center IDRIS (http://www.idris.fr) in France.

RESULTS
-RMSD plateaus of the MD trajectories imply physically stable systems S1b that of the inverse agonist ciproxifan (CPX). Fig. S2 shows the RMSDs of the Cα atoms with respect to the initial structure. For the antagonist-H3R complex, the RMSD attains a plateau at about 150 ns, remains constant at 4.2 Å until 550 ns and then rises slightly to 4.8 Å. The agonist-H3R complex shows an RMSD that attains values of 5.0 Å at 600 ns until the end of section 1. This observation is consistent with the fact that agonist ligands tend to be less efficient in the stabilization of the structure of the protein. 12 The apo receptor shows considerable changes in the RMSD, especially in the middle of the trajectory, just before 600 ns, when it reaches almost 5 Å of RMSD, to fall afterwards to around 3.5 Å at the end of the 930 ns trajectory. All in all, the observed fluctuations suggest a good packing quality of each model.
For the structurally-conserved TM helices of the receptor, the RMSDs of each of the three systems with respect to the crystal structure of the H1R-doxepin complex at 3.1 Å resolution 13 S6 (PDB code 3RZE) are as follows: 3.1 Å for the last frame of the antagonist-H3R complex, 2.5 Å for the last frame of the agonist-H3R complex, and 2.0 Å for the apo receptor. These data indicate that the MD trajectories lead to physically stable systems.
-RMSD matrices point to high structural variations  complex shows no proximity between TM2 and TM4, TM3 and TM7, and TM4 and TM5; instead, this complex shows a proximity between TM2 and TM6. As mentioned in the Introduction, H3R possesses a high constitutive activity in the absence of agonist. Our results S7 show indeed that the overall relative positions of TM helices in the apo receptor and the agonist-H3R active complex resemble each other. In addition, ECL2 is in contact with the antagonist, in conformity with the importance of this loop in ligand binding, selectivity and (in)activation, 14 as well as stabilization of the inactive state of the receptor. 15 -Eigenvectors show intra-molecular correlated movements and inter-molecular differential movements. corresponds to cluster 1 (C1) and shows that this mode involves essentially movements of the N-ter extracellular segment of the receptor and its ECL2, which act as a flexible lid that opens and closes on top of the ligand-binding cavity. Loops ICL2 and ECL2 exert a push-n-pull motion on TM4 as a rigid body. For TM5, we observe a compression-extension movement of the fragment above Pro 5.50 and a kink movement around this residue. We also observe kink movements around Pro 6.50 and Pro 7.50 of the conserved 6.47 CWXP 6.50 and 7.49 NPXXY 7.53 motifs, respectively. All these movements are correlated. No other important movements are detected in the other helices. Cluster C2 (Fig. 1a) shows large correlated movements of ECL2, ICL2 and to a lesser extent ICL3, with slight kink motions of TM1 and TM5-TM7 around their middle-helix proline residues (Fig. S5b). In addition, TM5-TM7 move collectively in a breathing motion along the plane of the membrane. The mode corresponding to cluster C3 (Fig. 1a) shows large amplitude motions for ECL2, and a pendulum-like rigid-body motion of the lower zone of TM5-7 with their IC loops (Fig. S5c). Fig. S5d shows the collective modes of movement corresponding to cluster C4. First is an angular N-ter to C-ter motion of H8, whose C-ter gets close and far from the membrane, accompanied by a concerted motion of the N-ter, ECL3 and  -Clusters of charged residues get reconfigured from one state of the receptor to the other Electrostatic interactions are of special importance for membrane proteins because of the low dielectric environment in membranes. We thus proceeded to detect clusters of charged residues using a distance criterion of 12 Å between COMs of the side chains of amino acid residues.
For the antagonist-H3R complex, there are two charged amino acid clusters that persist 90% of the time of the production trajectory. One cluster is in the extracellular region of the receptor and the other one in the intracellular region. The two clusters are represented in Fig.   S8. The amino acid residues for the first cluster belong to TM3, ECL2, TM5-TM7, and those of the second cluster to TM1, ICL1, TM3-TM6, ECL1, ICL2, ICL3 and H8. For the agonist-H3R complex, section 1, two independent clusters persist while HSM is bound to the receptor (Fig.   S9). They involve N-ter, TM3, ECL2, TM6 and TM7 for the first network; and TM3-TM6 for the second. One residue of the first cluster, Asp 3.32, is in contact with the protonated amine (pKa ~9.4) of the aliphatic amino group of HSM, consistent with experimental data. 20 Four clusters of charged residues can be found for the apo receptor (Fig. S10). Amino acid residues from TM1, TM3 and ECL2 compose the first cluster; one residue from each TM3, 5 and 6 the second cluster; residues from the N-ter, ECL1 and TM6 the third cluster; and residues from TM1, ICL1, S9 ICL2, TM3-TM5, ICL3, TM7 and H8 the fourth cluster. Many of these residues are also involved in long-lasting H-bonds, i.e. prevailing more than 70% of the time.
-H-bonded networks between amino acid residues are a function of the state of the receptor For the antagonist-H3R complex, inter-residue H-bonds with residence times greater than 70% connect ECL2 (TM4-TM5 loop) and ICL2 (TM3-TM4 loop) to TM3; N-ter to TM7; TM1 to H8; and ECL2 to ECL2 (Supplementary Table S1). For the first segment of the production trajectory of the agonist-H3R complex, ECL2 is H-bonded to the N-ter, TM3 and TM4; and the N-ter of TM1 to H8 (Supplementary Table S1). Finally, for the apo receptor, only one nonintra-helical H-bond persists (Supplementary Table S1). Thus, it is interesting to notice that the antagonist-bound receptor shows the largest number of inter-residue H-bonds, providing the structure of this complex with an added energetic and structural stability.

-The composition of the internal cavity of the receptor abounds in Tyr and Leu
Supplementary Table S2b shows that there are 38 residues that form the orthosteric cavity of the CPX-bound receptor. These residues are contributed by TM2, TM3, ECL2, and TM5-TM7.
The environment of the ligand is rather hydrophobic with the three aromatic amino acid side chains being represented. Twenty-six residues compose the internal cavity of the agonist-H3R complex, with contributions coming this time from TM1, TM3, ECL2, and TM5-TM7. For the apo receptor, the cavity is composed of only 24 residues. As expected, this cavity is the smallest one, given the absence of ligand and the subsequent contraction of the orthosteric cavity. Secondary structures contributing to it are TM2, TM3, ECL2, and TM5-TM7. In all three systems, TM4 being eccentric, it contributes with no residues to the receptor's cavity, just like all loops (except ECL2) and, of course, H8 (the short C-ter juxta-membrane helix in the cytoplasmic side). Notice the contribution of ECL2 to the morphology of the binding cavity.
The residues that compose the internal cavity of H3R in its different states are of diverse types, with Tyr and Leu being the most abundant, followed by Phe and Ser. Residues Pro, Gly, Gln, Lys, Thr and His are absent (Supplementary Table S2b). Interestingly, the antagonist and agonist pockets show a net negative charge of 4, whereas the net charge of the pocket of the S10 apo receptor is neutral (Asp 3.32 + Arg 6.58). The cationic ligands HSM and CPX thus induce changes in the pocket that bring into play acidic residues.  Table S2a), so that water penetration to the receptor is largest upon antagonist binding, and smallest for the agonist-H3R complex, just like in the high-resolution structure of α2A-AR that reveals about 60 internal waters. 21 In addition, Yuan et al. 22 mention an increased penetration of water into the receptor cavity of μ-and κ-opioid receptors, which has been linked to the activation mechanism upon agonist binding. Therefore, the presence of water in the binding site clearly demonstrates its influence in the dynamics and conformation of the receptor.
We determined the H-bonded networks established by water molecules in the cavity with occupancy greater than 75%. In all three systems, these occupancies are characterized by very low residence times, of the order of 20% maximum, indicating a high fluidity of water molecules. Thus, even if a site may be continually hydrated, it is not so by the same water  S11). Finally, the interaction of CPX with two residues in the protein, Tyr 2.61 and Asp 3.32, is mediated by two water molecules for the former, and by one for the latter (Network3, Fig.   S11). Notice that TM1 and, of course TM4 and H8, do not participate to the H-bond network during CPX binding. The only charged amino acids in the interior of the receptor are in the binding cavity and are Asp 2.50, Asp 3.32 and Glu 5.46 (Fig. S11). The agonist-H3R complex shows two networks of H-bonds, the first of which includes HSM bonded to a water molecule and to Asp 3.32 (Fig. S11). The two networks of the apo receptor involve, each one, two internal water molecules (Fig. S11).

-Hydrophobic clusters may contain up to seven residues
We report in Supplementary Table S3 the hydrophobic clusters at 90% occupancy formed by at least three side chains for the antagonist-H3R complex, the agonist-H3R complex and the apo receptor, respectively. The cutoff used is of 6 Å between the COMs of each amino acid residue; the set of residues is composed of Ile, Leu, Val, Phe, Met, Cys, Trp, Pro and Ala.
For the antagonist-H3R complex, there are six clusters. The distribution of clusters is as follows: two clusters of three residues, two clusters of four residues, one cluster of six residues and one cluster of seven residues. One cluster, implying the N-ter and ECL1 contains 7 residues (Ala 23 from the N-ter, Phe 2.56, Cys 2.57, Leu 2.60, Trp 23.50 (ECL1), Leu 3.24 and Trp 3.28).
For the agonist-H3R complex, there are six clusters, and for the apo receptor there are five clusters. These non-covalent interactions, just like the H-bond networks, depend on the state of the receptor.
-Side-chain switches for the Met residues are not detected Kofuku et al. 23 , by monitoring the NMR signals, investigated the role of Met residue 82 (2.53) in antagonist-and partial agonist-bound states of the β2-AR, which are correlated with conformational changes of the transmembrane regions upon activation. The corresponding residue in H3R is a valine and none of its neighbors is a methionine; nevertheless, we decided to monitor the conformational states of all Met residues for the three systems. We found that S12 in all three systems, all but one of the methionine residues remain in the t or g-states throughout the trajectory. As opposed to the β2-AR then, the conformational states of the methionine residues remain unchanged between the apo receptor, and its antagonist-and agonist-bound states; there is thus no correlation of conformational changes of the methionines upon activation in H3R.

-Detection of lipid binding sites on the receptor leads to amino acid sequence motifs
We list in Supplementary Table S4 the DPPC lipids attached to the receptor in the antagonist state and the amino acid residues they interact with. The residence time of a given lipid is listed in the first line of the "occupancy" column; the next lines show the residence times of the lipid with different amino acid residues of the receptor.
For the agonist-H3R complex, the following secondary structures participate to the binding of the seven lipids: ICL1, TM2, TM5, TM6 and the TM7-H8 loop (Fig. S12, Supplementary Table   S5). Interestingly, for the apo receptor the 5-lipid binding site involves all helices except TM5 ( Fig. S12 and Supplementary Table S6).
Supplementary Table S7 shows the amino acid residues of the receptor in contact with lipids for each of the three systems grouped. First, the fraction of non-polar residues (43-56%) is the largest, followed by the aromatic residues (15-26%), and the polar and charged classes. The most frequent residue of all groups and for all three states of the receptor is leucine. In the second class, Phe is the most abundant for the antagonist and agonist (section 1) systems, and Tyr and Trp for the apo system. In the polar group, Ser and Thr for the antagonist-and agonist-H3R complexes are the most represented residues. Finally, the antagonist-H3R complex is the system showing the most contacts of charged residues with lipids -eleven Arg and six Lys. The non-polar and aromatic class residues are in contact with the larger area acyl-tails of the lipids; whereas the basic residues interact with the negatively charged lipid head groups.
Supplementary Table S4 for the antagonist complex shows 13 highest-occupancy lipid molecules associated to the receptor, of which four are in the upper-leaflet of the membrane S13 and nine in the lower-leaflet (Fig. S12). As an illustration of the interaction of DPPC lipids with the inactivated H3R, Fig. S13 shows a 2D LigPlot+ diagram of lipids 125 and 149, showing Hbonds between Arg 8.51 and one oxygen atom from the phosphate head, and hydrophobic interactions with surrounding residues for lipid 149. Both lipids share Leu 7.55 as an interacting residue. In Fig. S14, extra-and intra-cytoplasmic views of the inactivated receptor for the 490 ns frame, we can see that DPPC molecules essentially bind to opposite sides of the In order to extract the lipid binding H3R amino acid sequence motifs, we selected those amino acid residues in Supplementary Table S4 to Supplementary Table S6 in contact with lipids that show residence times ≥ 90%. The resulting binding motifs are in Supplementary Table S8.
-Leu is the most frequently membrane-exposed residue The fraction of membrane-exposed hydrophobic amino acid residues is of 71, 71 and 67% for the antagonist, agonist (section 1) and apo systems, respectively, with Leu being the most frequently exposed residue, representing 27-29% of all membrane-exposed residues. With respect to all exposed residues, the hydrophilic residues, represent 22%, 12%, and 18% for the antagonist complex, respectively. As far as the charged residues is concerned, their relative S14 exposed populations are of 8%, 17% and 15%, for the antagonist complex, for the agonist complex and for the apo receptor, respectively.
-Rotamer toggle switches are multiple and concerted The goal of this section is to determine those amino acid side chains that undergo concomitant side-chain conformational changes. We focus on the aromatic side chains of Tyr, Trp and Phe, since several rotamer toggle switches dealing with those residues have been reported in the literature. 24  Fig. 5b). Thus, in a remarkable fashion a multiple toggle rotamer switch of certain aromatic chains, along with bi-modal switches appears to be important for receptor H3R inactivation. For the agonist-H3R complex (section 1), the conformation of 1 of Phe 5.47 remains constant at g-and Phe 3.41 undergoes a g-to trans transition in the last fourth of the trajectory (not shown). Concomitant transitions S15 in 1 of several aromatic side chains take place starting at about 0.5 µs (Fig. S17). Again, a multiple toggle rotamer switch of aromatic side chains takes place during H3R activation and involves several amino acid residues. Only Trp 6.48 is common to both mechanismsactivation, inactivation. All but one of the residues are in the upper zone of the receptor and imply the ECL1 and ECL2. The number of aromatic side chains that undergo conformational transitions in the apo receptor is much reduced and includes only Phe 45.54 and Tyr 45.56 of ECL2, and Tyr 7.33 and Phe 8.54 side chains, with the first and the third toggle switches interconverting between g+ (+60°), g-(-60°) and trans states, and the second and third toggle switches common to the activation mechanism (Fig. S18). The mentioned residues are in the upper and middle zones of the receptor, except for Phe 8.54, belonging to H8. None of the residues of the DRF ( 3.49 AspArgPhe 3.51 ) motif are thus involved in H3R activation or constitutive activity implicating that the ionic lock involving Arg3.50-(Asp 6.30 and Asn 5.64) is not conserved. This is consistent with the lack of direct evidence for a role of that motif in those two situations. 26 In résumé, there are multiple concerted long-range rotamer switches for each system as the involved residues are several tens of angstroms away from each other. Finally, it is interesting to observe that aromatic side chains switches from ECL2 participate in the different mechanisms.
-Inter-residue contacts are unique to each state of the receptor; novel ionic locks are detected We have measured many selected inter-residue distances using the COMs of the side chains to detect contacts. Those pairs of residues interacting according to the criteria described in the M&M section are in Supplementary Table S9 for each of the three systems. Several of those interactions represent ionic locks.
For the antagonist-H3R complex, the conserved arginine of the ionic lock (motif DRF) among rhodopsin-like G protein-coupled receptors, Arg 3.50, interacts with Asp 3.49 and Asp 6.30 in a constant but dynamic fashion throughout the trajectory, in agreement with the behavior of the ionic lock in the MD simulations of the β2-AR receptor. 27 In addition, there is formation of a hydrophobic interaction between Met 1.39 and Trp 7.40/Trp 7.43, followed by disappearance of an initial Met 1.39-Tyr 2.61 interaction; the interaction between Met 6.55 S16 and Tyr 3.33/Phe 7.39 is sporadic. Lastly, the distribution of the Met 6.55-Tyr 6.51 distance with a major peak at 5 Å and a minor one at 9 Å is analogous to one described in the literature (Nygaard et al 2013, their Fig. 6C). 28  -Diffusion of a potassium monocation into the antagonist complex sodium allosteric site.
K + begins its course towards the interior of the receptor between TM5 and TM6 and is in contact with Leu 5.39 and His 67.00 (ECL3); it is hydrated at this time by five water molecules.
A short interval after, it interacts with Ser 5.43, Glu 5.46, Phe 5.57 and Met 6.55 from TM5 and TM6, and is coordinated by only two to three water molecules. Towards the end of its path, it establishes contacts with Asp 2.50, Cys 3.36, Ser 3.39, Glu 5.46 and Trp 6.48 (TM2, TM3, TM5 and TM6), several of which are mediated by water, like with Asp and Glu residues, to finally become rehydrated by five water molecules. From time to time, the monocation is in contact with the highly conserved Asp 3.32 of the binding pocket. In the last frame of the trajectory, the cation is coordinated by Asp 2.50, Asp 3.32 and by three structured water S17 molecules (r<3.0 Å), just like in the 2A-AR. 21 As opposed to HSM, its insertion in the pocket is irreversible in the time scale studied (Fig. S19, Fig. S20). S18    Trp7.43. c) Apo receptor. In all three systems, all but one are water-mediated interactions. S26 Many of these interactions are polyvalent. The boxes correspond to the TM residues, the red circles to water molecules and the octagons to N-ter and loop residues.