Advanced Computational Methods in Molecular MedicineView this Special Issue
Molecular Modeling of the M3 Acetylcholine Muscarinic Receptor and Its Binding Site
The present study reports the results of a combined computational and site mutagenesis study designed to provide new insights into the orthosteric binding site of the human M3 muscarinic acetylcholine receptor. For this purpose a three-dimensional structure of the receptor at atomic resolution was built by homology modeling, using the crystallographic structure of bovine rhodopsin as a template. Then, the antagonist N-methylscopolamine was docked in the model and subsequently embedded in a lipid bilayer for its refinement using molecular dynamics simulations. Two different lipid bilayer compositions were studied: one component palmitoyl-oleyl phosphatidylcholine (POPC) and two-component palmitoyl-oleyl phosphatidylcholine/palmitoyl-oleyl phosphatidylserine (POPC-POPS). Analysis of the results suggested that residues F222 and T235 may contribute to the ligand-receptor recognition. Accordingly, alanine mutants at positions 222 and 235 were constructed, expressed, and their binding properties determined. The results confirmed the role of these residues in modulating the binding affinity of the ligand.
Muscarinic acetylcholine receptors (mAchRs) are integral membrane proteins that belong to the rhodopsin family of G-protein-coupled receptors (GPCRs). These proteins play a pivotal role in the regulation of many physiological functions both in the central and in the peripheral nervous systems . Like other GPCRs, these proteins exhibit a basal activity that is increased by the action of its natural agonist acetylcholine. In the central nervous system, mAChRs are known to regulate motor control, temperature homeostasis and are involved in the process of memory and learning whereas, in the peripheral system, they induce cardiovascular activity, smooth muscle contraction, gland secretion, regulation of cardiac activity, decrease blood pressure, and meiosis. The physiological actions the natural agonist acetylcholine are mediated by at least five subtypes of receptors known as M1-M5 that exhibit high sequence identity across mammalian species and exhibit different tissue distribution [2, 3].
Binding of agonists and competitive antagonists to GPCRs occurs at the orthosteric site, a hydrophobic pocket located at the extracellular side of the helix bundle that is highly conserved among the members of a subfamily . However, the binding of diffusing ligands to the orthosteric site may be modulated by residues located at the extracellular loops that flank the entrance . Interestingly, while residues defining the orthosteric binding site are well conserved within the different GPCR subfamilies, those of the extracellular loops are remarkably diverse within a subfamily, providing a good possibility to design selective allosteric antagonist . Although a wealth of information about these receptors has been accumulated in the past, there are still many questions open. Issues like the effect of the same mutation on the function of different mAchRs subtypes or the actual mapping of the binding site, as well as the mechanism of activation, are poorly understood . A deeper understanding of the structure-activity relationships of the different subtypes will help designing new selective ligands with various pharmacological profiles that may be useful for therapeutic intervention. In this sense, due to the scarcity of crystallographic structures available, molecular modeling methods can provide a deeper insight into the flexibility and the role that extracellular loops may play in ligand recognition .
In the present work, we focus on the human M3 muscarinic acetylcholine receptor (M3R), which is involved in modulation of neurotransmitter release, temperature homeostasis, and food intake in the central nervous system, as well as in the induction of smooth muscle contraction, gland secretion, indirect relaxation of vascular smooth muscle in the peripheral nervous system . M3R is arranged in the distinctive seven-transmembrane (TM) helix bundle architecture like all GPCRs, with its N-terminus located on the extracellular side of the plasma membrane and the C-terminal tail in the cytoplasmic side. Helices are connected by intracellular and extracellular loops of different lengths exhibiting preserved features of the rhodopsin family, like the two highly conserved cysteine residues C221 and C141, responsible for a disulfide bond between the second extracellular loop (ECL2) and the TM3 helix .
We undertook a combined approach involving computational and site-directed mutagenesis studies to further dissect M3R-ligand interactions with the aim of getting a deeper insight into specific structural features that can modulate the binding of the competitive antagonist N-methylscopolamine (NMS) to the orthosteric site of the receptor. For this purpose, we constructed an atomistic model of the M3R that was further refined using extended molecular dynamics (MD) simulations of the protein embedded in a lipid bilayer. The analysis of the MD trajectories provides a better understanding of the structural features that characterize the ligand-receptor interaction at the orthosteric pocket, the dynamics of the extracellular loops, as well as its putative involvement in ligand binding .
Present modeling studies suggest that residues F222 located in the ECL2 and T235 located at the edge of TM5 participate in the recognition of NMS. The former had already been shown to be important for the binding of allosteric modulators  but has not yet been shown to modulate the binding of ligands to the orthosteric site. On the other hand, the latter had already shown to be important in conferring biding specificity to the physiological ligand acetylcholine [12, 13]. Accordingly, in order to demonstrate the prediction capabilities of the model, we designed and expressed mutants F222A and T235A. Moreover, in order to have a reference, other residues of the ECL2 not expected to modulate the binding of NMS to M3 according to the atomistic model but playing key roles in the binding of allosteric modulators such as gallamine and alcuronium in the M2 receptor  were also investigated. Thus, the following mutants Y209F, K213 V, and E228N were also constructed and expressed.
Saturation binding assays were performed to determine the affinity of [3H]-NMS for these mutants. Dissociation assays on the mutants were also performed in order to get insight about the kinetics of the binding process. These results suggest that residues F222 and T235 modulate the binding of NMS  that show a slight decrease of its binding affinity in contrast with the other mutants constructed and tested.
2. Materials and Methods
2.1. Theoretical Methods
2.1.1. Construction of a M3 Receptor Model
A model of the M3 receptor was constructed by homology modeling using the crystal structure of bovine rhodopsin as a template (entry 1GZM of the Protein Data Bank)  by means of the Modeller software . The N-terminus as well as a long stretch of the C-terminus was not included in the model, due to the low sequence identity exhibited between the two sequences and that it does not presumably influence the model for the purpose of the present work. Residues 264–479, located in the third intracellular loop (ICL3), were also removed and substituted by stretch of alanine residues resulting into a shorter loop with the same length as that found in rhodopsin. The model includes the full remaining intradiscal and cytoplasmic loops, the conserved disulphide bridge between the cysteine residues located in TM3 and ECL2 (involving C141 and C221). All amino acids were modeled in the protonation state; they would have as free amino acids in water at pH 7 yielding a net charge of +17. This net charge is later compensated in the simulation box by adjusting the balance between sodium and chloride ions to give an electroneutral system, as described in detail later.
Fifty different models were generated following this procedure. These were subsequently ranked ordered using a scoring function based on structural compactness and embedded in Modeller . Structures were processed by Procheck in order to assess the stereochemical quality of the models generated . The structure finally selected exhibits the highest score with the minimum number of residues out of the preferred regions of the Ramachandran plot.
Although the constructed model provides a reasonable first approximation of the receptor, it still exhibits some limitations for its use in drug discovery . First, consideration of the rearrangement due to the bound ligand and, second, the effect of the lipid bilayer on the protein structure are disregarded in the model constructed using a crystal structure. These two effects can be accounted for in the model by means of molecular dynamics simulations.
2.1.2. MD Simulations
In order to understand the structural features of M3R, four 50 ns molecular dynamics simulations with the receptor embedded in a lipid bilayer were carried out. Specifically, four trajectories were run: with and without the antagonist NMS bound into the orthosteric binding pocket and using two different bilayer compositions: a one-component palmitoyl-oleyl phosphatidylcholine (POPC) and a two-component palmitoyl-oleyl phosphatidylcholine/palmitoyl-oleyl phosphatidylserine (POPC/POPS), where POPS lipids were added at the cytoplasmic surface of the M3R at a POPC/POPS 10 : 1 ratio. The former bilayer was considered in the present paper because it has been widely used in MD simulations; its available force field parameters are known to be carefully calibrated and it is a good representative of the membrane composition together with the fact that its gel/liquid-crystalline phase transition is above the temperature of the simulation (300 K). Selection of the latter—a two-component bilayer—was done in order to investigate the role of the inclusion of anionic lipids to neutralize the net positive charge located in the cytoplasmic region of the M3 receptor.
The M3R was placed in the center of an equilibrated POPC box of dimensions 8.5 × 8.5 × 10 nm (XYZ) containing 256 molecules, a thick layer of water molecules on each side of the bilayer, and NaCl ions at physiological concentration of 0.2 M taken from a previous study . Overlapping molecules were removed following the procedure described previously , in which all water molecules, with oxygen atoms closer than 0.40 nm to a nonhydrogen atom of the protein, as well as all lipid molecules with at least one atom closer than 0.25 nm to a nonhydrogen atom of the protein, were removed. This resulted in final system containing 189 lipids and ca. 14 600 water molecules. For convenience, the system is organized in such a way that the bilayer plane was oriented on the XY plane. For the construction of the POPC/POPS box, the 10% of the POPC molecules exhibiting the largest electrostatic potential were replaced by POPS molecules. For those simulations with N-methylscopolamine (NMS) bound into the orthosteric site, the ligand was manually docked, using information provided by site-directed mutagenesis on relevant residues involved in the binding . Figure 1 shows the system showing a cartoon of the M3 receptor with NMS bound, surrounded by POPC molecules. Water molecules have been removed for clarity.
Calculations were carried out by means of the GROMACS 3.3.2 package  using periodic boundary conditions. The all-atom OPLS force field  was used for all the molecules of the system except for the lipids, for which a specific parameterization was used . This combination has been previously employed . Systems were subjected to periodic boundary conditions in the three-coordinate directions. Calculations were done at 300 K using separate thermostats for the protein, water, ions, and lipid molecules. The time constant for the thermostats was set to 0.1 ps, except for water, for which a value of 0.01 was used. Equations of motion were integrated using the leapfrog algorithm using a time step of 2 ps.
For equilibration, the system was subjected to 0.5 ns MD simulation to allow for the removal of voids present between the protein and the lipid or water molecules. The simulation was performed with the atomic coordinates of the protein restrained to their crystallographic positions and allowing the three periodic box dimensions to change size according to a pressure of 0.1 MPa in each coordinate direction . Next, the restrains were released and the four MD trajectories containing rhodopsin were computed up to 20 ns each. As a reference, eight additional simulations were performed without protein and either with or without ions lasting 20 ns each. For all simulations, coordinates were collected every 10 ps and were stored for further analysis. In all cases, the first 10 ns were considered as equilibration period, and therefore they are not included in the analysis. All the bonds in the protein and lipid molecules were kept frozen using the LINCS algorithm. Electrostatic interactions were treated using the particle mesh Ewald summation procedure  and nonbonded interactions were computed using a cutoff of 1.0 nm.
3. Construction, Expression, and Binding Assays for Mutant M3R
3.1. Materials and Mutant Construction
1-[N-methyl-3H] scopolamine methyl chloride ([3H]-NMS; 84 Ci/mmol) was purchased from Perkin Elmer. Oligonucleotides were designed and obtained from Sigma Aldrich (for the sequence of the primers, see Figure 1(S) of Supplementary Materials available online at http://dx.doi.org/10.1155/2012/789741). The Quick Change kit was from Stratagene (La Jolla, CA, USA). Atropine was purchased from Sigma Aldrich.
Residues F222 and T235 of the human M3R were mutated to alanine using the Quick Change method. Mutant receptors, cloned into the pCD expression vector and confirmed by dideoxy sequencing, were transiently expressed in COS-7 cells by electroporation (Bio-Rad Gene Pulser, 15 μg of DNA, 260 V, 960 μF). In the case of the T235A mutant, its expression levels were improved by treatment of the cultured cells with atropine (10-6 M) for 48 h, before washing and harvesting for membrane preparation. The ability of atropine to act as a “pharmacological chaperone” has been described previously . After transfection, media was removed and cells were washed twice with phosphate buffered saline (PBS) (in the case of mutant T235A, four additional washes were needed to remove atropine), then harvesting buffer (20 mM Hepes, 10 mM EDTA) was added to the flasks. The cells were incubated for about 10 min and then scraped from the flask. The cells were pelleted and homogenized in a glass pestle (×20 plunges on ice). This homogenate was transferred to JA-17 tubes and centrifuged at 17,000 rpm at 4°C for 30 minutes. The pellet was resuspended in storage buffer (20 mM Hepes, 0.1 mM EDTA). The cells were homogenized on ice using a Polytron homogenizer setting 12,2 × 15 s, the membranes resuspended, snap frozen, and stored at −80°C. Protein concentration was determined using the Bradford assay .
3.2. Radioactive Ligand Assays
[3H]-NMS saturation binding assays were performed by binding [3H] NMS (0.01–3 nM) to membrane preparations (10 μg/mL of membrane protein). Binding was measured at 30°C in a buffer containing 20 mM Hepes, 100 mM NaCl, and 1 mM MgCl2, pH 7.5, using an assay volume of 1 mL and an incubation time of 2.5 h. Nonspecific binding was determined with 1 μM atropine. Assays were performed in quadruplicate. The binding reaction was terminated by rapid filtration on a Brandel cell harvester (Brandel Model M-245, Serial 9219). Dissociation kinetic assays were carried out to investigate the time course of [3H]-NMS dissociation from the receptor. COS7 cell membranes expressing the M3R (10 μg/mL) were equilibrated with 0.2 nM [3H]-NMS in 1-mL total volume of buffer containing 20 mM Hepes, 100 mM NaCl, and 1 mM MgCl2, pH 7.5 for 60 min at 30°C. Atropine (10 μM) was added at various time points to prevent radioligand re-association to the mAchRs. Data analysis was carried out by means of Graph Pad Prism 3.0 (Graph Pad Software, San Diego, CA, USA). Saturation binding data were fit with hyperbolae (one-site binding); Bmax and Kd values were derived from these curves. Dissociation kinetic data were normalized and fit to for a monoexponential decay curves. Experiments were repeated 4 times. Log affinities were tabulated as mean ± SEM. Statistical comparison of affinity and dissociation rate constants for mutants and WT M3R controls were carried out by one-way analysis of variance followed by Dunnett’s post hoc test.
4. Results and Discussion
As described in the Materials and Methods section, we used rhodopsin as template  for the generation of a starting structure of the M3R by homology modeling. Comparison of the transmembrane helices yields a 21.7% sequence identity and a 42.9% sequence homology. Obviously, any of the available crystal structures could have also been used for this purpose, being the choice of a template to construct the initial atomistic model of a GPCR matter of debate in the scientific community [29, 30]. Since for several years rhodopsin was the only high-resolution GPCR structure available, most of the homology modeling made use of it as template . The recently published squid rhodopsin structure  and a few ligand-activated structures the human β2 adrenergic receptor (β2AR) [33, 34], the avian β1 adrenergic receptor (β1AR) , the human A2A adenosine , and more recently the crystal structure of the CXCR4 chemokine  and the dopamine receptor  have brought the opportunity to increase the repertoire of templates that can be used for modeling purposes, opening the question of which is the most suitable structure to be used. Indeed, although the overall structure of the six proteins is similar, there are notable differences, mostly in the extracellular loops and in the ligand-binding region [39, 40]. For example, in contrast to the other structures, the human β2 adrenergic receptor as well as of the avian β1 adrenergic receptors exhibit an alpha helix in the second extracellular loop. Similarly, TM5 and TM6 in squid rhodopsin are more extended than in the other receptors, or the inverse agonist in rhodopsin is more deeply inserted inside the hydrophobic TM domain of the protein. One of the conserved structural features, the so-called ionic lock, is weaker in the case of the β1AR, β2AR and A2A adenosine receptors than in the case of rhodopsin and involves the ICL2. This latter difference may be responsible for the higher basal activity and structural instability of these receptors and may contribute to the challenges in obtaining diffraction crystals of nonrhodopsin GPCRs .
However, there are still inconclusive results about the choice of a template. Thus, at the community-wide assessment of GPCR structure modeling hold in 2008  aimed at understanding of the performance of GPCR homology modeling, 206 models of the A2A adenosine receptor submitted and assessed for their accuracy prior the release of the solved crystal structure. These models were basically assessed by the accuracy of the ligand-binding mode. The results show that among the best models there are some that were constructed using the structures of the β-adrenergic receptors and some using rhodopsin. On the other hand, in a recent report an accurate model of the transmembrane region of the β2-adrenergic receptor was reported using rhodopsin as template .
Two recently studies address the problem of template selection for GPCR homology modeling [29, 30]. The two studies conclude that the selection of the template should not exclusively be based on sequence identity but on an a priori assessment of the expected structural features of the receptor wished to model and that more high-resolution structures are still required. Moreover, the latter study concludes that the best template needs to be considered case by case, recommending the use of chimeras constructed by selecting special features from the different receptors available. Unfortunately, the six crystallographic structures available nowadays do not cover all the structural diversity necessary to have adequate templates to model a specific GPCR of the rhodopsin family leaving the question still open.
5. MD Simulations
Deviations from the initial 7TM homology-based structure were monitored through the time evolution of the root-mean-square deviation (rmsd) of the α-carbon atoms as displayed in Figure 2(a). Smaller deviations are observed when the protein is embedded in the two-component lipid bilayer (POPS-POPC) than when the protein is embedded in the one-lipid environment (POPC) (Figure 2(a)). This behavior is probably due to the presence of the anionic lipid (POPS) in the cytoplasmic area of the receptor, neutralizing the cluster of positive charges located in that region, and it is consistent with the experimental finding that this part of the protein is difficult to crystallize due to its flexibility, as previously found in MD simulations of bovine rhodopsin . The time evolution of the RMSD also permits to assess the influence of the ligand on the protein structure. The evolution of the α-carbon rmsd shows that the trajectories with the ligand bound display lower deviations in comparison to those simulations of the protein alone.
Deviations per residue were also computed for a better characterization of the contribution of each segment to the overall rmsd. Figure 2(b) shows the rmsd per residue, highlighting the lower values found in the TM regions, compared to the connecting loops. Specifically, the larger deviations correspond to residues on the ICL2 and ICL3 of the protein (deviations larger than 1 nm), whereas the smaller ones correspond to the residues located in the TM region of the protein (about 0.25 nm). This differential behavior can be correlated with the fact that receptor/G-protein coupling occurs through the intracellular loops in all mAChRs . Furthermore, the structures embedded into the two-component lipid bilayer (POPS/POPC) exhibit, in general, smaller deviations as compared to those embedded into one-component lipid bilayer with exceptions such as residues in the intracellular loop ICL1 (135–140), in ICL2 (225–233), and residues in the ICL3 (520–535).
The root mean square fluctuations (RMSF) of the α-carbon atoms, taken as reference the average structure computed from the last 25 ns were also calculated for the four trajectories. The largest fluctuations are located on the second and third intracellular loops as anticipated in the rmsd calculations (for residue rmsf see Figure 2(S) of Supplementary Materials). Fluctuations were relatively small due to the restrictions imposed by the lipid bilayer, except on the non-TM regions. In the case of TM4, TM5 and TM6 show an overall larger fluctuation in comparison to the other TM helices. These results agree with the fact that some residues within TM3-TM7 participate in the binding of the antagonist ligand NMS [26, 45, 46].
Dihedral angles were also computed using four consecutive α-carbons for the four MD simulations and disregarding the first 12 ns (Table 1), in order to identify specific conformations adopted by the extracellular loops when the NMS antagonist is bound and in the ligand-free receptor. The analysis of the M3-NMS-POPC MD trajectory shows that ECL1 adopts three different conformations, whereas it adopts only one conformation in the other three simulations. In contrast, ECL3 adopts only one conformation in the four MD trajectories. On the other hand, the time evolution of the dihedral angles defined using four consecutive Cα of ECL2 shows that it attains two different conformations along the trajectories M3R-NMS-POPC and M3R-NMS-POPC/POPS, whereas in the M3-NMS-free POPC/POPS simulation exhibits three different conformations.
This observation agrees with the putative contribution of ligands to generate microdomains on the extracellular loops that are important for the conformational changes that accompany receptor activation [47, 48]. Interestingly, we can also see that most notable changes in the variation of the dihedral angle average occur in the preceding part of the residues 218–221; therefore it seems that residues near the C221, the one that forms a disulphide bond with C141, play an important role in the conformation of this second extracellular loop. Figure 3 shows the superposition of two snapshots at different times taken from the POPC simulation, where it can clearly be seen the distinctive features the conformation adopted by the ECL2 before the first 20 ns and after 20 ns of the simulation.
In order to study different specific structural features, we focused on residue-residue hydrogen-bonding interactions, and lipid-protein interactions. For the protein-protein intramolecular interactions, we focused on hydrogen bond formation and breaking between pairs of residues of the protein along the MD trajectory (Figure 4(a)). The analysis was carried out by counting all possible atoms involved in a hydrogen/charge interaction around one residue at a time using a cutoff of 0.4 nm. Figure 7 lists the interactions between residues and their persistence along the MD trajectories. They are classified in charged-charged, charged-noncharged and noncharged-noncharged. Among the residues listed in Figure 7, those forming the so-called ionic-lock between TM3 and TM6, in subscript the Ballesteros-Weinstein numbering is used : D1653.49, R1663.50, Y1673.51 and E4866.30, R2535.60 and R1804.37 are indicated. Electrostatic interactions at these positions are known to be critical for restraining the receptor from adopting a constitutively active GPCR [50, 51]. Some of the interactions are only found in the two POPS/POPC simulations: D142-K213, D165-R166, E259-R480, E220-Y128, E220-Y530, and N548-T50. Other interactions can only be found in the POPC simulations: E257-R253, E486-R177, S121-Y534, T534-S121, T550-Q98, and T550-T554. Furthermore, some of the interactions are found in the simulations without ligand: D165-R184, E257-Y167, N104-Y544, Y149-T235. Finally, some residues are only found when the ligand is bound to the receptor: D165-R180, E220-K523, D148-S152, N104-T101, and S537-S152 (Figure 7).
On the other side, specific residues have been shown to be involved in protein-lipid interactions. These include charged like K and R, neutral like N, Q, T, and S, and aromatic like H, Y, and W . Previous analysis of hydrogen bonding among protein donors of bovine rhodopsin and lipid oxygen atoms revealed the importance of these interactions for the anchoring of the protein to the membrane. In rhodopsin, residues R2526.35, Y962.63, Y2746.57, H1524.41 and T1083.23, and to a lower extent K661.61, R691.64, Y1353.50, and W351.30 are found to participate in hydrogen bond linkages with lipid molecules, and, therefore, they can be considered to be putative hooks of the protein to the bilayer . In this work, specific hydrophilic lipid-protein interactions have been analyzed along the different trajectories using the criteria of any of the residue atoms of the receptor was within a cutoff of 0.4 nm of a lipid oxygen. Following this approach, we identified common interactions for the four trajectories studied (Figure 4(b)). Most of the lipid-protein interactions are detected at the cytoplasmic half of the receptor, similarly to what is found in the case of rhodopsin. Moreover, some of the residues in M3R, W661.30, K1001.64, T1272.63, Y1673.51, and R1844.41, correspond to homologous residues in the rhodopsin sequence that confers the visual photoreceptor the property of anchoring the protein onto the lipid bilayer .
Another point of interest arises from the comparison of the interactions using the two different bilayers. In this case, the location of POPS lipids in the cytoplasmic vicinity of the receptor would provide a larger number of these lipid-protein interactions due to the fact that carboxylate groups of POPS lipids can interact with positively charged residues such as K260, R261, K263, R480, K487, and K488 located at ICL3 of the receptor.
6. The Orthosteric Binding Pocket of M3R
NMS was docked onto the binding pocket of M3R taking into account the residues known to affect its binding from site-directed mutagenesis studies, specifically D3.32 known to interact with the ligand quaternary nitrogen , N6.52 putatively involved in a polar interaction with the ligand [53–55], and the homologous hydrophobic residues implicated in the binding of NMS at the M1 muscarinic receptor . Once docked, the complex was energy minimized and the MD simulation started. This process was done for the two different bilayer compositions studied in the present work, POPC and a 10 : 1 POPC/POPS mixture. In order to understand the stereochemical features of the binding pocket, we monitored all the residues within a cutoff of 0.35 nm of any of the atoms of the ligand along the MD trajectories. Residues involved in the binding of NMS are shown in Figure 5. Simultaneously, we monitored the different conformations the ligand attained during the simulations. The analysis reveals that the ligand actually adopts different conformations along the trajectories; however, all of them conserve the interaction between the quaternary nitrogen of NMS with the side chain of D1483.32 and the hydroxyl moiety with N5086.52 (for definition of NMS dihedral angles and their values along the MD trajectory see Figures 3(S) and 4(S) of the Supplementary Material). Two different representative conformations of the ligand along the MD trajectories are depicted in Figure 6.
In a detailed analysis of hydrogen bonding interactions between different residues of the protein (those within a cut-off of 0.35 nm) and atoms of NMS, we computed the average angle formation of each of the interactions along the dynamics, and we found that residues N5086.52, W5046.48, F222, and T235 of M3R participate in hydrogen bond formation (angle between 120°–180°) with specific atoms of the antagonist ligand NMS in the MD simulations. In agreement with previously published results literature , N6.52 is considered to be important for the binding of some muscarinic antagonists , and also W6.48 has also been shown to be important for the binding of NMS [45, 51].
We have constructed, expressed, and characterized the alanine mutants of F222 and T235 and studied them by means of saturation binding assays in order to determine the affinity of NMS for these mutants. Dissociation binding assays have also been performed to study the release of NMS as a function of time from the M3R, in 20 mM Hepes buffer. In fact, the corresponding homologous positions have shown to modulate binding of NMS to the M1 receptor [47, 58]. Thus, mutant Y177 (the corresponding position of F222) reduced the affinity of NMS about 2.2-fold in comparison to WT M1 . Furthermore, the mutant T1925.42 (the corresponding position of T235) showed a reduction of the affinity of NMS about 1.69 compared to the WT M1 .
7. Expression and Experimental Characterization of WT M3R and Mutants
The effect of the mutations studied on the expression of mAChR was measured by saturation binding assays using [3H] NMS (for representative saturation binding curves see Figure 5(S) of the Supplementary Materials). Mutant F222A showed about 55% of expression level in comparison to the WT M3R. In the case of the T235A mutant, it showed lower expression levels that were improved after atropine treatment. Alanine substitution of residues F222 and T235 somehow resulted in reduced affinities for NMS when compared to WT M3R. Mutant F222A decreased the affinity for NMS by 1.7-fold and mutant T235A by 2.6 fold with regard to WT M3R (Table 2). However, Y209F, K213 V, and E228N mutants assayed in parallel did not show any significant difference with regard to WT M3R (data not shown) which is consistent with the proposed role for these position in allosteric—but not in orthosteric—ligand binding . The measurement of dissociation time courses was carried out in order to provide information about the transition from the conformation of the receptor bound to the radioligand and the free state of the receptor. [3H] NMS dissociated from the WT M3 receptor with a monoexponential time course time course corresponding to a half-time of 10 min under this condition (20 mM Hepes buffer, 30°C) which is in good agreement with previously published values . Dissociation rates for the mutants F222A and T235A were not significantly different from the WT M3 receptor (data not shown). The effect observed in the binding affinities, together with the lack of effect in the dissociation analysis, suggests different structural environment requirements for the binding and dissociation of the ligand to and from the receptor.
The present paper addresses the effect of lipid composition and the antagonist ligand on the structure of the human M3R. The atomic resolution model of the M3R constructed by homology modeling using the crystal structure of bovine rhodopsin is compatible with experimental data from site directed mutagenesis studies and represents the inactive form of the receptor. This model was refined by MD simulations and the effect of the ligand on the conformation of the protein analyzed. Therefore, MD trajectories of different component lipid bilayers (POPC and the two-component POPC/POPS lipid bilayer), and also with or without NMS within the orthosteric binding site of M3R receptor, were analyzed.
We could determine that POPS lipid molecules provide additional stabilization to the protein structure of human M3R, especially for ICL3 because charged residues localized in this area of the receptor interact with the lipid oxygen atoms of the POPS molecule. Furthermore, our lipid-protein analysis suggests that some specific residues in the GPCR family (positions 1.30, 1.64, 2.63, 3.51, and 4.41) are responsible for anchoring the membrane protein onto the lipid bilayer (e.g., rhodopsin and M3R), as previously hypothesized [20, 60].
Our modeling identified many polar interactions that could play a role as activation microswitches [61, 62]. This includes also interactions associated to the “ionic lock” located at the cytoplasmic ends of the third and sixth transmembrane helices, respectively, (positions: 3.49, 3.50, 3.51, 6.30, and 6.34 in the GPCR family) play an important role in the receptor activation mechanism [50, 63].
In agreement with previous molecular modeling and mutagenesis for other muscarinic receptors (Avlani et al., 2007) , we found that ECL2 has a requisite of flexibility promoting the binding of GPCR ligands despite of the presence of the disulphide bridge between this loop and TM3. As seen in our average dihedral angle calculations, this loop showed more flexibility in comparison to ECL1 and ECL3. This loop fluctuates in two different conformations (most notable changes occur in the residues preceding C221—involved in the disulfide bridge—in the two MD simulations in which the ligand NMS was inside the binding pocket.
On the other hand, our experimental ligand binding and dissociation results suggest that, although F222 and T235 are not critical in the binding of NMS, these residues may participate in the conformation of a optimal binding pocket of NMS, where every small contribution in the environment confers stability to the NMS-receptor complex. Further studies, using other experimental conditions, will be needed to determine the potential role of these residues in binding of other ligands, for example, agonists or allosteric modulators.
M. Martinez Archundia wishes to acknowledge CONACYT (Mexico) for a predoctoral fellowship. This work has been supported by the Spanish Ministry of Science and Innovation through the Grants SAF2008-04943-C02-01 (to J. J. Perez) and SAF2008-04943-C02-02 (to P. Garriga).
Figure 1S. Forward and reverse primers for the mutants F222A and T235A.
Figure 2S. Root mean square fluctuations (rmsf) of Cα atoms, from an average structure computed from the last 25ns of the trajectories, were also calculated for the four trajectories.
Figure 3S. Definition of the dihedral angles of the antagonist NMS that were measured along the simulations.
Figure 4S. Values of the four dihedral angles of NMS measured in the different simulations: (a) POPC and (b) POPC/POPS.
Figure 5S. Representative Saturation Binding Curves of WT M3R and mutants F222A and T235A expressed in COS7 cells. The membranes were incubated with increasing concentrations of [3H] NMS in Hepes buffer at 30°C for about 2.5 hours. The protein concentration used was ca 10μg/mL. The specific and non-specific binding are expressed in fmol/mg protein.
M. P. Caulfield and N. J. M. Birdsall, “International union of pharmacology. XVII. Classification of muscarinic acetylcholine receptors,” Pharmacological Reviews, vol. 50, no. 2, pp. 279–290, 1998.View at: Google Scholar
V. A. Avlani, K. J. Gregory, C. J. Morton, M. W. Parker, P. M. Sexton, and A. Christopoulos, “Critical role for the second extracellular loop in the binding of both orthosteric and allosteric G protein-coupled receptor ligands,” Journal of Biological Chemistry, vol. 282, no. 35, pp. 25677–25686, 2007.View at: Publisher Site | Google Scholar
X. P. Huang and J. Ellis, “Mutational disruption of a conserved disulfide bond in muscarinic acetylcholine receptors attenuates positive homotropic cooperativity between multiple allosteric sites and has subtype-dependent effects on the affinities of muscarinic allosteric ligands,” Molecular Pharmacology, vol. 71, no. 3, pp. 759–768, 2007.View at: Publisher Site | Google Scholar
J. Wess, D. Gdula, and M. R. Brann, “Site-directed mutagenesis of the M3 muscarinic receptor: identification of a series of threonine and tyrosine residues involved in agonist but not antagonist binding,” EMBO Journal, vol. 10, no. 12, pp. 3729–3734, 1991.View at: Google Scholar
J. Wess, R. Maggio, J. R. Palmer, and Z. Vogel, “Role of conserved threonine and tyrosine residues in acetylcholine binding and muscarinic receptor activation. A study with M3 muscarinic receptor point mutants,” Journal of Biological Chemistry, vol. 267, no. 27, pp. 19313–19319, 1992.View at: Google Scholar
A. Krejčí and S. Tuček, “Changes of cooperativity between N-methylscopolamine and allosteric modulators alcuronium and gallamine induced by mutations of external loops of muscarinic M3 receptors,” Molecular Pharmacology, vol. 60, no. 4, pp. 761–767, 2001.View at: Google Scholar
Z. L. Lu, J. W. Saldanha, and E. C. Hulme, “Transmembrane domains 4 and 7 of the M1 muscarinic acetylcholine receptor are critical for ligand binding and the receptor activation switch,” Journal of Biological Chemistry, vol. 276, no. 36, pp. 34098–34104, 2001.View at: Publisher Site | Google Scholar
R. A. Laskowski, M. W. MacArthur, D. S. Moss, and J. M. Thornton, “PROCHECK—a program to check the stereochemical quality of protein structures,” Journal of Applied Crystallography, vol. 26, pp. 283–291, 1993.View at: Google Scholar
O. Berger, O. Edholm, and F. Jähnig, “Molecular dynamics simulations of a fluid bilayer of dipalmitoylphosphatidylcholine at full hydration, constant pressure, and constant temperature,” Biophysical Journal, vol. 72, no. 5, pp. 2002–2013, 1997.View at: Google Scholar
T. Darden, D. York, and L. Pedersen, “Particle mesh Ewald: an N·log(N) method for Ewald sums in large systems,” The Journal of Chemical Physics, vol. 98, no. 12, pp. 10089–10092, 1993.View at: Google Scholar
M. M. Bradford, “A rapid and sensitive method for the quantitation of microgram quantities of protein utilizing the principle of protein dye binding,” Analytical Biochemistry, vol. 72, no. 1-2, pp. 248–254, 1976.View at: Google Scholar
J. Wess, S. Nanavati, Z. Vogel, and R. Maggio, “Functional role of proline and tryptophan residues highly conserved among G protein-coupled receptors studied by mutational analysis of the M3 muscarinic receptor,” EMBO Journal, vol. 12, no. 1, pp. 331–338, 1993.View at: Google Scholar
E. Ramon, A. Cordomí, L. Bosch et al., “Critical role of electrostatic interactions of amino acids at the cytoplasmic region of helices 3 and 6 in rhodopsin conformational properties and activation,” Journal of Biological Chemistry, vol. 282, no. 19, pp. 14272–14282, 2007.View at: Publisher Site | Google Scholar
U. Ringkananont, J. Van Durme, L. Montanelli et al., “Repulsive separation of the cytoplasmic ends of transmembrane helices 3 and 6 is linked to receptor activation in a novel thyrotropin receptor mutant (M626I),” Molecular Endocrinology, vol. 20, no. 4, pp. 893–903, 2006.View at: Publisher Site | Google Scholar
E. Kurtenbach, C. A. M. Curtis, E. K. Pedder, A. Aitken, A. C. M. Harris, and E. C. Hulme, “Muscarinic acetylcholine receptors. Peptide sequencing identifies residues involved in antagonist binding and disulfide bond formation,” Journal of Biological Chemistry, vol. 265, no. 23, pp. 13702–13708, 1990.View at: Google Scholar
E. C. Hulme, M. S. Bee, and J. A. Goodwin, “Phenotypic classification of mutants: a tool for understanding ligand binding and activation of muscarinic acetylcholine receptors,” Biochemical Society Transactions, vol. 35, no. 4, pp. 742–745, 2007.View at: Google Scholar
K. Allman, K. M. Page, C. A. M. Curtis, and E. C. Hulme, “Scanning mutagenesis identifies amino acid side chains in transmembrane domain 5 of the M1 muscarinic receptor that participate in binding the acetyl methyl group of acetylcholine,” Molecular Pharmacology, vol. 58, no. 1, pp. 175–184, 2000.View at: Google Scholar
S. Lazareno and N. J. M. Birdsall, “Detection, quantitation, and verification of allosteric interactions of agents with labeled and unlabeled ligands at G protein-coupled receptors: interactions of strychnine and acetylcholine at muscarinic receptors,” Molecular Pharmacology, vol. 48, no. 2, pp. 362–378, 1995.View at: Google Scholar
P. V. Escribá, P. B. Wedegaertner, F. M. Goñi, and O. Vögler, “Lipid-protein interactions in GPCR-associated signaling,” Biochimica et Biophysica Acta, vol. 4, pp. 836–852, 2007.View at: Google Scholar
A. Cordomí, E. Ramon, P. Garriga, and J. J. Perez, “Molecular dynamics simulations of rhodopsin point mutants at the cytoplasmic side of helices 3 and 6,” Journal of Biomolecular Structure and Dynamics, vol. 25, no. 6, pp. 573–587, 2008.View at: Google Scholar