Article pubs.acs.org/jcim
Molecular Dynamics Simulations of the Adenosine A2a Receptor in POPC and POPE Lipid Bilayers: Effects of Membrane on Protein Behavior Hui Wen Ng,† Charles A. Laughton,‡ and Stephen W. Doughty†,* †
School of Pharmacy, University of Nottingham Malaysia Campus, Jalan Broga, 43500 Semenyih, Selangor, Malaysia School of Pharmacy, University of Nottingham, Nottingham, NG7 2RD, United Kingdom
‡
S Supporting Information *
ABSTRACT: Analysis of 300 ns (ns) molecular dynamics (MD) simulations of an adenosine A2a receptor (A2a AR) model, conducted in triplicate, in 1-palmitoyl-2-oleoylphosphatidylcholine (POPC) and 1-palmitoyl-2-oleoylphosphatidylethanolamine (POPE) bilayers reveals significantly different protein dynamical behavior. Principal component analysis (PCA) shows that the dissimilarities stem from interhelical rather than intrahelical motions. The difference in the hydrophobic thicknesses of these simulated lipid bilayers is potentially a significant reason for the observed difference in results. The distinct lipid headgroups might also lead to different molecular interactions and hence different protein loop motions. Overall, the A2a AR shows higher mobility and flexibility in POPC as compared to POPE.
■
budding, fission, and fusion.11 In general, sphingomyelin, phosphatidylcholine, and cholesterol are regarded as the main constituents of the membrane outer layer, while phosphatidylethanolamine, phosphatidylserine, and cholesterol are those of the membrane inner layer.12−18 To date, it remains very challenging to construct a realistic model for in silico membrane protein simulations. This is due not only to the complex and varied nature of the cell membrane composition that depends on the cell type10,19 but also to the difficulty in force field parametrization of the various membrane components. In view of this, simplified approaches have been adopted for these simulations: quite often, single species lipid bilayers are used;20−24 though in other studies somewhat more realistic bilayers containing two or more lipid species including cholesterol have been employed.25,26 Although running the risk of oversimplification, these practices have a number of advantages: (1) They greatly reduce the complications that come with force field parametrization and testing against experimental data. (2) They allow each component of the membrane protein system to be treated rigorously, and (3) they provide a well-defined system that facilitates clearer data interpretation.4 1-Palmitoyl-2-oleoylphosphatidylcholine (POPC) and 1palmitoyl-2-oleoylphosphatidylethanolamine (POPE) (see Figure 1) lipid bilayers are among the most widely used single species lipid bilayers for in silico simulations.27−30 This could be attributed to these reasons: (1) They resemble the naturally
INTRODUCTION The adenosine A2a receptor (A2a AR) belongs to the family of G-protein coupled receptors (GPCRs), a group of highly versatile receptors capable of responding to an impressive plethora of stimulus e.g. light, hormones and neurotransmitters. The A2a AR is ubiquitously expressed in the body and plays a pivotal role both in important normal response functions (e.g., inflammation1) and under pathological conditions (e.g., Parkinson’s disease2). As with all GPCRs, the A2a AR resides in the cell membrane, which has a remarkably complicated and heterogeneous architecture. The cell membrane consists of a variety of glycerophospholipids (e.g., phosphatidylcholine, phosphatidylethanolamine, and phosphatidylserine), sphingolipids (e.g., sphingomyelin, cerebrosides and gangliosides), cholesterol, and membrane proteins. These components often organize into domains of distinct compositions, physical properties, and functions.3−6 Lipid rafts, for instance, despite the fact that their existence remains controversial due to a lack of detection and isolation methods, are widely believed to be regions of cell membrane that consist of sphingolipids and cholesterols and play an important role in colocalizing membrane proteins and signaling molecules.3,6−9 The glycerophospholipids, phosphatidylcholine and phosphatidylethanolamine, are two of the major neutral lipids found in the cell membrane.10 Phosphatidylcholine, which accounts for more than 50% of eukaryotic cell membrane lipids, has a structural purpose and commonly serves as the fluid-like solvent for the lipid rafts and membrane proteins.6,10 Phosphatidylethanolamine, on the other hand, has been found to be important in biological processes such as © 2014 American Chemical Society
Received: August 5, 2013 Published: January 24, 2014 573
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
Article
achieved using the g_membed tool,42 which was available as part of GROMACS 4.5.3.43−45 The water model used was the single point charge (SPC)46 model. To remove the net charge of the system, 11 chloride ions were added as counterions. The equilibration process started with a 100 ps simulation under NVT conditions at 310 K. The simulation parameters were as follows: electrostatic interactions were treated with the smooth particle mesh Ewald (PME)47 method with a short-range cutoff of 1.2 nm; van der Waals interactions were also given a shortrange cutoff of 1.2 nm; all bonds were constrained with LINCS48 to enable a 2 fs time step to be applied; the temperature was coupled to the v-rescale49 thermostat with a time constant of 0.1 ps. Coordinates were written to the output trajectory file every 100 ps. Position restraints were applied to the heavy atoms of the A2a AR model along the x, y, and z axes as well as the phosphorus atom of POPC/POPE along the z axis only, where x−y is the membrane plane. Under these restraints, the water molecules were allowed to move in all directions while the lipids were able to do so only in the x−y plane. After NVT, three separate rounds of NPT were carried out to allow the position restraints to be gradually released. Maintaining the restraints from the NVT stage, the first round of NPT was carried out for 10 ns, with semi-isotropic pressure coupling to the Parinello−Rahman barostat50 and a time constant of 1 ps. The second round was performed for 20 ns with position restraints applied only to the protein. All the position restraints were removed in the third round, and the systems were allowed to equilibrate for 40 ns. Finally, production simulations of 300 ns were performed in triplicate under the same NPT conditions. The only difference between replicates was in the initial velocity assignments at the start of the dynamics. PCA was performed using our in-house toolkit pcazip (http://holmes.cancres.nottingham.ac.uk/pcazip), and the analysis was divided into four parts (overall, TM region, interhelical, and intrahelical motions). The six simulations (three replicates for each of the POPC and POPE systems) were combined before PCA, so the dynamics of the A2a AR model in the two different model lipid systems have been compared within a common subspace. For the top n principal components that capture 90% of the variance, the dot products and subspace overlaps between the replicate simulations were also calculated using pcazip. Further Analyses. The bilayer thickness of the POPC and POPE systems were measured using GridMAT-MD51 at the following stages of the simulations: pre-equilibration, postequilibration, and post-production. The thickness was approximated using the distance between the phosphorus (P) atoms of the upper and lower membrane layers. The compactness and the changes in the conformation of the A2a AR model in the x, y, and z axes were investigated by calculating the radius of gyration of the protein. The potential dynamic effect of the electrostatic interactions between the opposite charge groups on the lipids and protein were also studied, specifically between the positively charged quarternary nitrogens of POPC and POPE and the negatively charged aspartic and glutamic acids of the A2a AR model. This was achieved through calculating the backbone root-mean-square-fluctuations (RMSFs) of the aspartic acid and glutamic acid residues in close proximity to the lipid headgroups using GROMACS 4.5.3.43−45 All the plots and graphs were produced using gnuplot, while the pictures were generated using Chinera52 and VMD.53
Figure 1. Structures of POPC and POPE lipids.
occurring phospholipids which commonly contain a saturated sn-1 and unsaturated sn-2 acyl chains.31,32 (2) They permit the formation of secondary structures such as a bilayer.31,33 (3) They interact with other lipids such as cholesterol,31,34 and (4) their simple structure and availability of experimental data allow them to be easily modeled. In this study, we have set out to investigate the effects of using these two lipid bilayers on membrane protein behavior during molecular dynamics (MD) simulations. The apo form of the A2a AR model was inserted into both POPC and POPE lipid bilayers, and simulations of 300 ns were performed on both systems in triplicate. Principal component analysis (PCA) was used to analyze the dynamics of the A2a AR model at the level of (1) the whole protein, (2) just the transmembrane (TM) regions, (3) just between the TM helices, and (4) within helices I−VII individually. This approach has been used in our previously published studies of ligand-induced modulations into A2a AR dynamics35 and proved useful in providing a detailed analysis of the protein motions in a quantitative and unbiased way. To the best of our knowledge, this work may potentially be one of the first molecular simulation studies that compares the effects of different single species lipid bilayers on the dynamical behavior of a GPCR. This is fairly significant as a great number of GPCR membrane simulations have been conducted using single species lipid bilayers.20,35−39 Therefore, the current study is important in providing both specific insights for this system and a generalized methodology for examining how the choice of lipid bilayer can affect membrane protein simulation results, thus providing a basis for more informed decisions on the choice of simulation parameters in the future.
■
METHODS System Setup, Simulation Protocols, and PCA. The majority of the methods used in this section have been described in our previous work.35 The homology model of the complete human A2a AR was constructed on the basis of the engineered A2a AR crystal (3EML)40 with deletion of the 160residue-long T4-lysozyme portion (N1002−N1161), SO4− ions, and stearic acid molecules. Crystallographic waters were preserved, but the ligand ZM241385 was removed to generate an “apo” model. The missing residues in the crystal structure were addedi.e., M1-I3 (N-terminus), P149−H155 (extracellular loop 2, EL2), and K209-A221 (intracellular loop 3, IL3). The final structure lacked the C terminus and contained only residues 1 to 310. A brief round of steepest descent energy minimization was performed to remove any steric clashes in the models arising from the loop building. The A2a AR structure was then inserted into two pre-equilibrated and fully hydrated POPC and POPE lipid bilayers of 339 and 340 molecules, respectively (the original “Berger lipids”41 obtained from http://moose.bio.ucalgary.ca/) with the embedding process 574
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
■
Article
RESULTS Overall Protein Behavior: A2a AR Model in POPC and POPE. One of the most striking observations that can be made when the replicate simulations of the A2a AR model in POPC (A2a AR-POPC) and POPE (A2a AR-POPE) are considered together is that the replicates from these systems occupy distinctly different regions of conformational space, as defined by principal component 1 (PC1) but not by principal component 2 (PC2) (see Figure 2). The A2a AR model also
that they are directly responsible for the most fundamental difference observed here between the two systems. Elimination of the flexible loop regions from our analyses greatly reduces, but does not eliminate, the gap between the POPC and POPE simulations. A higher degree of motional flexibility remains obvious in the TM region of the POPC systems (see Figure 4) Interhelical Behavior: A2a AR in POPC and POPE. The interhelical component of the TM motions resembles strongly the overall TM region motions reported above (see Supporting Information Figure 1). This is in agreement with the findings in our previous work35 that TM flexibility is dominated by interhelical motions. Intrahelical Behavior: A2a AR in POPC and POPE. Interestingly, when analyzed at the level of individual TM helices, the protein dynamics are much less affected by the identity of the membrane lipid, though some differences remain in terms of conformational space sampled, degree of flexibility, and reproducibility in behavior between replicate simulations (see Figure 5). Helices I and VII demonstrate the highest degree of flexibility, particularly in POPC. In contrast, helix IV is the most stable. Helices V and II could also be considered stable if ignoring the anomalous behavior of replicate 1 of A2a ARPOPC (which is an acceptable approach given that the deviation is only minor and clearly not representative as shown by the tight overlap in the replicates 2 and 3 in both cases). Helix III is flexible in POPE but stable in POPC, while the opposite was true for helix VI. Eigenvector dot products and subspace overlaps are two alternative quantitative metrics that offer a different perspective on the results presented thus far. They offer the advantage of comparisons between different trajectories without the need to enforce a common subspace. The average maximum dot products for both POPC and POPE systems were in these ranges: overall protein ca. 0.50−0.60, TM region ca. 0.60−0.70, interhelical ca. 0.80−0.85, intrahelical ca. 0.80−0.95 (see Supporting Information Figure 2−6). The subspace overlaps were in these ranges: overall protein ca. 0.30−0.35, TM region ca. 0.35−0.45, interhelical ca. 0.50−0.55, intrahelical ca. 0.45− 0.70 (POPC) and 0.45−0.80 (POPE). These values suggest a continuous increase in convergence as the study progressed from the overall protein motion analysis to the intrahelical motion analysis. Moreover, the dot products and subspace overlaps for the intrahelical components recapitulate the variation in dynamical characteristics (sampling, stability, convergence) of individual helices in POPC and POPE described above. Rationalization of the Differences between A2a AR in POPC and POPE. The results from this study make it clear that the A2a AR model shows significantly different behavior when simulated in POPC and POPE lipid bilayers. One of the reasons might relate to the fact that POPC and POPE lipid bilayers possess different membrane thicknesses.27,54 We find that the thickness of the POPC bilayer, which had an initial value of 4.0 nm, remains fairly constant throughout the simulations. On the other hand, the POPE bilayer, which had a slightly greater initial thickness of 4.2 nm, expands further in the z direction to 4.3−4.5 nm. The effects of the membrane thickness on the protein overall conformation and compactness have been analyzed through the calculation of the radii of gyration of the A2a AR in POPC and POPE, in the x, y, and z directions (see Figure 6a−c). In POPC, the values for the x (red) and y (green) components
Figure 2. Projections of the A2a AR-POPE (left) and A2a AR-POPC (right) simulations onto the common subspace defined by PC1/PC2 for the whole protein. The following color scheme and symbol are used throughout this article: red, replicate 1 (POPC); green, replicate 2 (POPC); blue, replicate 3 (POPC); magenta, replicate 1 (POPE); cyan, replicate 2 (POPE); yellow, replicate 3 (POPE); and ▲, conformation at the start of production phase.
demonstrates higher flexibility and dynamical variability in POPC compared to POPE. As with our previous work,35 the protein in POPC is found to drift away from the conformation at the start of production phase as the simulations progress and shows little similarity in behavior between the replicate simulations. In POPE, however, the protein demonstrates more constrained motions whereby only a restricted area of conformational space is explored. As a result, the replicate simulations appear to be closer to convergence in POPE than in POPC. Visualization of PC1 (represented here through an overlay of the maxima of structural distortion along the principal components, see Figure 3a) reveals noticeably varied motions between the two systems in both the loop and TM regions: (1) All the loops, especially intracellular loop 3 (IL3), demonstrate major movements in the POPC membrane but more subdued motions in POPE. (2) Helices I, V, VI, and VII show fairly large motions, particularly in the extracellular portion in a POPC environment; less motions are observed in the helical bundle in POPE, although helices VI and VII of replicate 1 manifest a fairly prominent stretching, especially in the middle portion. In general, the helical motions in POPC could be described as the bending and stretching of the helices toward and away from the binding pocket leading to its expansion and contraction. Animation of PC2 (see Figure 3b) shows arguably more similar A2a AR model motions in both membrane environments, particularly in the TM regions while the loop motions remain random. In general, a concerted swinging motion of helices I and VII at the extracellular portion is observed. TM Region Behavior: A2a AR in POPC and POPE. The loop regions of a GPCR are the parts of the protein that are in most contact with the lipid headgroups, so it is not surprising 575
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
Article
Figure 3. Overlay plots representing the (a) PC1 and (b) PC2 motions of the three replicates from A2a AR-POPC (upper panel) and A2a ARPOPE (lower panel) simulations. The backbone atoms are shown for the overlapping structures, from minimum (blue) to maximum (red) value along the projection.
remain stable, fluctuating around 2.0 nm. Along the z (blue) axis, the value remains around 1.4 nm, although there is a slight increase toward the end. In POPE, the value in the x direction is larger than in POPC, but the value for the y component is comparable in both lipid systems. The value for the z component is notably greater in POPE than POPC, showing a steady climb toward 1.6 nm by the end of the 300 ns simulations. The changes in the A2a AR model conformation and compactness in response to the differing membrane thickness of POPC and POPE are shown in Figure 6d. In agreement with the radius of gyration data, this figure shows a more compact A2a AR model in POPC and an elongated A2a AR model (particularly along the z axis) in POPE. It should
however be noted that other factors might well be at play, and while the membrane thickness clearly differs, this might be an effect rather than a cause of the observed differences. Apart from membrane thickness, the difference in the POPC and POPE headgroups is expected to be a major factor that contributes toward the dynamical dissimilarities observed in A2a AR-POPC and A2a AR-POPE simulations. Comparing the −N+(CH3)3 group of POPC and −NH3+ group of POPE, the sterically hindered quaternized nature of the former reduces the potential for favorable interactions, e.g., hydrogen bonding and electrostatic interactions to take place between the lipid and the protein. 576
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
Article
between aspartic/glutamic acids and the −NH3+ group of POPE is reflected in a reduced level of dynamical flexibility of these amino acids in the POPE membrane environment. Both aspartic and glutamic acids show much lower mobility in POPE compared to POPC regardless of whether they are located in the helices or the loops. This potentially contributes to the overall reduced motion and flexibility of the A2a AR model in POPE.
■
DISCUSSION The distinct nature of the POPC and POPE lipids has resulted in a considerable difference in the A2a AR protein dynamical behavior during MD simulations. The essence of this difference is captured in the most important principal component of the PCA (i.e., PC1), indicating that it dominates the total variance in A2a AR global flexibility. The results from PCA of the whole protein, TM region, and interhelical motions collectively point toward the fact that this difference stems from a combination of the divergent loop movements (overall more mobile in POPC than POPE) and the dissimilar interhelical motions between the two lipid systems. The PCA of the intrahelical motions demonstrates that fundamentally the behavior of the individual
Figure 4. Projections of the A2a AR-POPC and A2a AR-POPE simulations onto the common subspaces defined by the PC1/PC2 for the TM-region motions.
Analysis of the interactions between charged amino acids in close proximity to the head groups suggests differences in dynamic behavior between the two systems (see Supporting Information Figure 7). The stronger attractions formed
Figure 5. Projection of the A2a AR-POPC and A2a AR-POPE simulations onto the common subspaces defined by the PC1/PC2 for the individual helices I−VII. Results are presented in three different panels for clearer visualization: combined results from both systems (main), POPC results (lower left), and POPE results (lower right). 577
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
Article
Figure 6. (a−c) Radii of gyration of the A2a AR model over 300 ns MD. Red denotes the changes occurring along the x axis, green in the y axis, and blue in the z axis. (d) The thicknesses (along the z axis) of the lipid bilayer during the simulations are shown through the black arrows. The P atoms in POPC and POPE are colored green and red, respectively.
both bilayers reveals a more compact protein in the thinner POPC and a more elongated one in the thicker POPEa possible consequence of the hydrophobic matching between the protein and the lipid bilayer.58,59 A mismatch of the protein hydrophobic belt and the lipid hydrophobic segment is a major player in membrane-protein energetics60 and the main cause of the physical adjustments that take place between the different membrane components. The adjustment of the membrane bilayer thickness could occur for stiff proteins such as those of β-barrel architecture; for softer α-helical bundles (such as GPCRs), the protein will more likely acquiesce through helical tilting or rotation of residue side chains, therefore leading to minimal membrane alteration.61 As expected, in this study, the
helices is much less sensitive to the lipid, although variations are present. This represents an interesting result hinting not only at the differences arising from simulation environments but also to providing potential understanding of GPCR−lipid interactions. A full discussion on functional selectivity and modulation is however beyond the scope of this study. The difference in interhelical structure and dynamics in the two membrane environments could potentially be a consequence of the protein’s adaptation to the different hydrophobic thickness of POPC and POPE bilayers. The range of membrane hydrophobic thicknesses observed in this investigation (3.9−4.5 nm) is in agreement with a number of other previous studies.27,55−57 The radii of gyration of the A2a AR in 578
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
■
latter appears to prevail. Anderson and Koeppe60 have discussed a number of particularly important points relevant to this study: hydrophobic matching may constrain the conformational space explored by the protein, thus stabilizing certain protein conformations in different lipid bilayers; individual α-helices tend to be rigid, and the protein deformation is more likely to be achieved through interhelical movements. Previously, the difference in the hydrogen bonding abilities of phospatidylcholine and phosphatidylethanolamine and their influence on membrane protein functions (e.g., active vs inactive) have been discussed.54,62 In addition, the electrostatic interactions of the lipid zwitterionic headgroups with the charged residues of the A2a AR model could lead to modified dynamic behavior with a subsequent impact on the rigidity and flexibility of the protein in the lipid bilayers. In essence, to capture the protein−lipid interplay more accurately, different levels of protein−lipid interactions should be considered. First, the specific interactions between protein and individual lipid molecules should be taken into account. This is demonstrated by a number of crystal structures, e.g., bacteriorhodopsin (1QHJ),63 bovine rhodopsin (1GZM),64 human β2-adrenergic receptor (3D4S),65 and human A2a AR adenosine receptor (4EIY),66 which reveal lipid molecules tightly bound to the receptors, potentially indicating the important roles of these individual lipids to the protein function. Second, the overall membrane material properties should also be carefully considered. While one of these properties, i.e., membrane thickness, has been discussed in detail, others such as the membrane curvature and elastic stress have not been covered.10,67−69 These properties have been suggested to be capable of altering (1) the free energy between conformational states, therefore changing the equilibrium of distribution and the transition kinetics, or (2) the height of the transition state, thus affecting only the transition kinetics.60
■
Article
ASSOCIATED CONTENT
S Supporting Information *
Information including the PC1/PC2 plot of the interhelical component of the TM motions, the dot product matrices, and the RMSF calculations are provided. This material is available free of charge via the Internet at http://pubs.acs.org.
■
AUTHOR INFORMATION
Corresponding Author
*Phone: +603 89248200. E-mail: Stephen.doughty@ nottingham.edu.my. Notes
The authors declare no competing financial interest.
■ ■
ACKNOWLEDGMENTS H.W.N. is grateful for the support from the University of Nottingham, Malaysia Campus for the funding for this research. ABBREVIATIONS ns, nanoseconds; MD, molecular dynamics; A2a AR, adenosine A2a receptor; PCA, principal component analysis; POPC, 1palmitoyl-2-oleoylphosphatidylcholine; POPE, 1-palmitoyl-2oleoylphosphatidylethanolamine; GPCR, G-protein coupled receptor; TM, transmembrane; IL3, intracellular loop 3; RMSF, root-mean-square fluctuations; PC1, principal component 1; PC2, principal component 2
■
REFERENCES
(1) Montesinos, M. C.; Gadangi, P.; Longaker, M.; Sung, J.; Levine, J.; Nilsen, D.; Reibman, J.; Li, M.; Jiang, C. K.; Hirschhorn, R.; Recht, P. A.; Ostad, E.; Levin, R. I.; Cronstein, B. N. Wound healing is accelerated by agonists of adenosine A2 (G alpha s-linked) receptors. J. Exp. Med. 1997, 186, 1615−1620. (2) Mally, J.; Stone, T. W. Potential of adenosine A(2A) receptor antagonists in the treatment of movement disorders. CNS Drugs 1998, 10, 311−320. (3) Lee, A. Membrane structure. Curr. Biol. 2001, 11, R811−R814. (4) Mouritsen, O. G.; Biltonen, R. L. Protein-lipid interactions and membrane heterogeneity. In Protein-Lipid Interactions, Watts, A., Ed.; Elsevier Science Publishers B.V: New York, 1993; Vol. 25. (5) Wiener, M. C.; White, S. H. Structure of a fluid dioleoylphosphatidylcholine bilayer determined by joint refinement of x-ray and neutron diffraction data. II. Distribution and packing of terminal methyl groups. Biophys. J. 1992, 61, 428−433. (6) Simons, K.; Ikonen, E. Functional rafts in cell membranes. Nature 1997, 387, 569−572. (7) Calder, P. C.; Yaqoob, P. Lipid rafts–composition, characterization, and controversies. J. Nutr. 2007, 137, 545−547. (8) Eggeling, C.; Ringemann, C.; Medda, R.; Schwarzmann, G.; Sandhoff, K.; Polyakova, S.; Belov, V. N.; Hein, B.; von Middendorff, C.; Schonle, A.; Hell, S. W. Direct observation of the nanoscale dynamics of membrane lipids in a living cell. Nature 2009, 457, 1159− 1162. (9) Leslie, M. Mysteries of the cell. Do lipid rafts exist? Science 2011, 334, 1046−1047. (10) van Meer, G.; Voelker, D. R.; Feigenson, G. W. Membrane lipids: where they are and how they behave. Nat. Rev. Mol. Cell Biol. 2008, 9, 112−124. (11) Marsh, D. Lateral pressure profile, spontaneous curvature frustration, and the incorporation and conformation of proteins in membranes. Biophys. J. 2007, 93, 3884−3899. (12) Feigenson, G. W. Phase boundaries and biological membranes. Annu. Rev. Biophys. Biomol. Struct. 2007, 36, 63−77. (13) Silvius, J. R. Fluorescence energy transfer reveals microdomain formation at physiological temperatures in lipid mixtures modeling the
CONCLUSIONS
In conclusion, this study has illustrated the differences of the A2a AR model protein dynamics in MD simulations as a result of using different lipid bilayers, i.e., POPC and POPE within in silico simulations. Clearly, the POPC and POPE bilayers bias the conformational sampling of A2a AR during simulations. While the physical adaptation of the A2a AR to membrane thickness could be one of the explanations for the differences in the interhelical motions, the molecular interactions of the lipid headgroups and the protein could help to provide a rationale for the patterns of rigidity/flexibility observed in the loops and individual helices. Finally, it is almost impossible to make a recommendation as to which lipid bilayer would be the best to use when performing GPCR membrane simulations, particularly those containing only a single species of lipid. However, it is clear from this work that it is important to be aware of the effect that a single lipid choice has on membrane protein simulations and to analyze results from such simulations in light of that knowledge. Clearly, there is a strong rationale for the development of more realistic membrane models, including components such as cholesterol,26,70 so that their potential interactions71,72 with the helices of A2a AR can also be accounted for. 579
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
Article
outer leaflet of the plasma membrane. Biophys. J. 2003, 85, 1034− 1045. (14) Wang, T. Y.; Silvius, J. R. Cholesterol does not induce segregation of liquid-ordered domains in bilayers modeling the inner leaflet of the plasma membrane. Biophys. J. 2001, 81, 2762−2773. (15) van Meer, G. Lipid traffic in animal cells. Annu. Rev. Cell. Biol. 1989, 5, 247−275. (16) Gallet, P. F.; Zachowski, A.; Julien, R.; Fellmann, P.; Devaux, P. F.; Maftah, A. Transbilayer movement and distribution of spin-labelled phospholipids in the inner mitochondrial membrane. Biochim. Biophys. Acta 1999, 1418, 61−70. (17) Gurtovenko, A. A.; Vattulainen, I. Lipid transmembrane asymmetry and intrinsic membrane potential: two sides of the same coin. J. Am. Chem. Soc. 2007, 129, 5358−5359. (18) Zachowski, A. Phospholipids in animal eukaryotic membranes: transverse asymmetry and movement. Biochem. J. 1993, 294 (Pt 1), 1− 14. (19) Pike, L. J.; Han, X.; Chung, K. N.; Gross, R. W. Lipid rafts are enriched in arachidonic acid and plasmenylethanolamine and their composition is independent of caveolin-1 expression: a quantitative electrospray ionization/mass spectrometric analysis. Biochemistry 2002, 41, 2075−2088. (20) Filizola, M.; Wang, S. X.; Weinstein, H. Dynamic models of Gprotein coupled receptor dimers: indications of asymmetry in the rhodopsin dimer from molecular dynamics simulations in a POPC bilayer. J. Comput.-Aided Mol. Des. 2006, 20, 405−416. (21) Martinez-Archundia, M.; Cordomi, A.; Garriga, P.; Perez, J. J. Molecular Modeling of the M3 Acetylcholine Muscarinic Receptor and Its Binding Site. J. Biomed. Biotechnol. 2012, 12. (22) Lai, P. C.; Crasto, C. J. Beyond modeling: all-atom olfactory receptor model simulations. Front. Genet. 2012, 3, 61. (23) Rosenbaum, D. M.; Zhang, C.; Lyons, J. A.; Holl, R.; Aragao, D.; Arlow, D. H.; Rasmussen, S. G.; Choi, H. J.; Devree, B. T.; Sunahara, R. K.; Chae, P. S.; Gellman, S. H.; Dror, R. O.; Shaw, D. E.; Weis, W. I.; Caffrey, M.; Gmeiner, P.; Kobilka, B. K. Structure and function of an irreversible agonist-beta(2) adrenoceptor complex. Nature 2011, 469, 236−240. (24) Zare, B.; Madadkar-Sobhani, A.; Dastmalchi, S.; Mahmoudian, M. Prediction of the Human EP1 Receptor Binding Site by Homology Modeling and Molecular Dynamics Simulation. Sci. Pharm. 2011, 79, 793−816. (25) Pitman, M. C.; Grossfield, A.; Suits, F.; Feller, S. E. Role of cholesterol and polyunsaturated chains in lipid-protein interactions: molecular dynamics simulation of rhodopsin in a realistic membrane environment. J. Am. Chem. Soc. 2005, 127, 4576−4577. (26) Shan, J.; Khelashvili, G.; Mondal, S.; Mehler, E. L.; Weinstein, H. Ligand-dependent conformations and dynamics of the serotonin 5HT(2A) receptor determine its activation and membrane-driven oligomerization properties. PLoS Comput. Biol. 2012, 8, e1002473. (27) Leekumjorn, S.; Sum, A. K. Molecular Characterization of Gel and Liquid-Crystalline Structures of Fully Hydrated POPC and POPE Bilayers. J. Phys. Chem. B. 2007, 111, 6026−6033. (28) Ceccarelli, M.; Marchi, M. Molecular dynamics simulation of POPC at low hydration near the liquid crystal phase transition. Biochimie 1998, 80, 415−9. (29) Marrink, S. J.; Mark, A. E. The mechanism of vesicle fusion as revealed by molecular dynamics simulations. J. Am. Chem. Soc. 2003, 125, 11144−11145. (30) Sperotto, M. M.; May, S.; Baumgaertner, A. Modelling of proteins in membranes. Chem. Phys. Lipids 2006, 141, 2−29. (31) Li, S.; Lin, H. N.; Wang, Z. Q.; Huang, C. Identification and characterization of kink motifs in 1-palmitoyl-2-oleoyl- phosphatidylcholines: a molecular mechanics study. Biophys. J. 1994, 66, 2005− 2018. (32) Barenholz, Y.; Thompson, T. E. Sphingomyelin: biophysical aspects. Chem. Phys. Lipids 1999, 102, 29−34. (33) Wang, X.; Quinn, P. J. Cubic phase is induced by cholesterol in the dispersion of 1-palmitoyl-2-oleoyl-phosphatidylethanolamine. Biochim. Biophys. Acta, Biomembr. 2002, 1564, 66−72.
(34) Ohvo-Rekila, H.; Ramstedt, B.; Leppimaki, P.; Slotte, J. P. Cholesterol interactions with phospholipids in membranes. Prog. Lipid Res. 2002, 41, 66−97. (35) Ng, H. W.; Laughton, C. A.; Doughty, S. W. Molecular Dynamics Simulations of the Adenosine A2a Receptor: Structural Stability, Sampling and Convergence. J. Chem. Inf. Model. 2013, 53 (5), 1168−1178. (36) Rodriguez, D.; Pineiro, A.; Gutierrez-de-Teran, H. Molecular Dynamics Simulations Reveal Insights into Key Structural Elements of Adenosine Receptors. Biochemistry 2011, 50, 4194−4208. (37) Schlegel, B.; Laggner, C.; Meier, R.; Langer, T.; Schnell, D.; Seifert, R.; Stark, H.; Holtje, H. D.; Sippl, W. Generation of a homology model of the human histamine H(3) receptor for ligand docking and pharmacophore-based screening. J. Comput.-Aided Mol. Des. 2007, 21, 437−453. (38) Huber, T.; Botelho, A. V.; Beyer, K.; Brown, M. F. Membrane model for the G-protein-coupled receptor rhodopsin: hydrophobic interface and dynamical structure. Biophys. J. 2004, 86, 2078−2100. (39) Goetz, A.; Lanig, H.; Gmeiner, P.; Clark, T. Molecular Dynamics Simulations of the Effect of the G-Protein and Diffusible Ligands on the I2̂ 2-Adrenergic Receptor. J. Mol. Biol. 2011, 414, 611− 623. (40) Jaakola, V.-P.; Griffith, M. T.; Hanson, M. A.; Cherezov, V.; Chien, E. Y. T.; Lane, J. R.; Ijzerman, A. P.; Stevens, R. C. The 2.6 Angstrom Crystal Structure of a Human A2A Adenosine Receptor Bound to an Antagonist. Science 2008, 322, 1211−1217. (41) Berger, O. E. O.; Jahnig, F. Molecular dynamics simulations of a fluid bilayer of dipalmitoylphosphatidylcholine at full hydration, constant pressure, and constant temperature. Biophys. J. 1997, 72, 2002−2013. (42) Wolf, M. G.; Hoefling, M.; Aponte-Santamaria, C.; Grubmuller, H.; Groenhof, G. g_membed: Efficient insertion of a membrane protein into an equilibrated lipid bilayer with minimal perturbation. J. Comput. Chem. 2010, 31, 2169−2174. (43) Berendsen, H. J. C.; van der Spoel, D.; van drunen, R. GROMACS: A message-passing parallel molecular dynamics implementation. Comput. Phys. Commun. 1995, 91, 43−56. (44) van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: Fast, flexible, and free. J. Comput. Chem. 2005, 26, 1701−1718. (45) Hess, B. K. C.; van der spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435−447. (46) Berendsen, H. P. J.; van Gunsteren, W.; Hermans, J. Interaction Models for Water in Relation to Protein Hydration. In In Intermolecular Forces; Pullman, B., Ed.; D. Reidel Publishing Company: Dordrecht, The Netherlands, 1981; pp 331−342. (47) Essmann, U.; Perera, L.; Berkowitz, M. L. A smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577−8593. (48) Hess, B. B.; Henk. LINCS: A Linear Constraint Solver for Molecular Simulations. J. Comput. Chem. 1997, 18, 1463−1472. (49) Bussi, G. D. D.; Parinello, M. Canonical sampling through velocity rescaling. J. Chem. Phys. 2007, 126, 014101. (50) Parinello, M.; Rahman, A. Polymorphic transitions in single crystals: A new molecular dynamics method. J. Appl. Phys. 1981, 52, 7182−7190. (51) Allen, W. J.; Lemkul, J. A.; Bevan, D. R. GridMAT-MD: a gridbased membrane analysis tool for use with molecular dynamics. J. Comput. Chem. 2009, 30, 1952−1958. (52) Pettersen, E. F.; Goddard, T. D.; Huang, C. C.; Couch, G. S.; Greenblatt, D. M.; Meng, E. C.; Ferrin, T. E. UCSF Chimera - A visualization system for exploratory research and analysis. J. Comput. Chem. 2004, 25, 1605−1612. (53) Humphrey, W.; Dalke, A.; Schulten, K. VMD - Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (54) Bondar, A.-N.; del Val, C.; White, S. H. Rhomboid Protease Dynamics and Lipid Interactions. Structure 2009, 17, 395−405. 580
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
Journal of Chemical Information and Modeling
Article
(55) Shityakov, S.; Dandekar, T. Molecular dynamics simulation of POPC and POPE lipid bilayers enforced by an intercalated single-wall carbon nanotube. Nano 2011, 06, 19−29. (56) Gullingsrud, J.; Schulten, K. Lipid bilayer pressure profiles and mechanosensitive channel gating. Biophys. J. 2004, 86, 3496−3509. (57) Pantano, D. A.; Klein, M. L. Characterization of membraneprotein interactions for the leucine transporter from Aquifex aeolicus by molecular dynamics calculations. J. Phys. Chem. B 2009, 113, 13715−13722. (58) Harroun, T. A.; Heller, W. T.; Weiss, T. M.; Yang, L.; Huang, H. W. Theoretical analysis of hydrophobic matching and membranemediated interactions in lipid bilayers containing gramicidin. Biophys. J. 1999, 76, 3176−3185. (59) Harroun, T. A.; Heller, W. T.; Weiss, T. M.; Yang, L.; Huang, H. W. Experimental evidence for hydrophobic matching and membranemediated interactions in lipid bilayers containing gramicidin. Biophys. J. 1999, 76, 937−945. (60) Andersen, O. S.; Koeppe, R. E., 2nd. Bilayer thickness and membrane protein function: an energetic perspective. Annu. Rev. Biophys. Biomol. Struct. 2007, 36, 107−130. (61) Lee, A. G. Lipid-protein interactions in biological membranes: a structural perspective. Biochim. Biophys. Acta 2003, 1612, 1−40. (62) Urban, S.; Wolfe, M. S. Reconstitution of intramembrane proteolysis in vitro reveals that pure rhomboid is sufficient for catalysis and specificity. Proc. Natl. Acad. Sci. U. S. A. 2005, 102, 1883−1888. (63) Belrhali, H.; Nollert, P.; Royant, A.; Menzel, C.; Rosenbusch, J. P.; Landau, E. M.; Pebay-Peyroula, E. Protein, lipid and water organization in bacteriorhodopsin crystals: a molecular view of the purple membrane at 1.9 A resolution. Structure 1999, 7, 909−917. (64) Li, J.; Edwards, P. C.; Burghammer, M.; Villa, C.; Schertler, G. F. Structure of bovine rhodopsin in a trigonal crystal form. J. Mol. Biol. 2004, 343, 1409−1438. (65) Hanson, M. A.; Cherezov, V.; Griffith, M. T.; Roth, C. B.; Jaakola, V. P.; Chien, E. Y.; Velasquez, J.; Kuhn, P.; Stevens, R. C. A specific cholesterol binding site is established by the 2.8 A structure of the human beta2-adrenergic receptor. Structure 2008, 16, 897−905. (66) Liu, W.; Chun, E.; Thompson, A. A.; Chubukov, P.; Xu, F.; Katritch, V.; Han, G. W.; Roth, C. B.; Heitman, L. H.; AP, I. J.; Cherezov, V.; Stevens, R. C. Structural basis for allosteric regulation of GPCRs by sodium ions. Science 2012, 337, 232−236. (67) Brown, M. F. Modulation of rhodopsin function by properties of the membrane bilayer. Chem. Phys. Lipids 1994, 73, 159−180. (68) Soubias, O.; Teague, W. E., Jr.; Hines, K. G.; Mitchell, D. C.; Gawrisch, K. Contribution of membrane elastic energy to rhodopsin function. Biophys. J. 2010, 99, 817−824. (69) Woolf, T. B.; Roux, B. Structure, energetics, and dynamics of lipid-protein interactions: A molecular dynamics study of the gramicidin A channel in a DMPC bilayer. Proteins 1996, 24, 92−114. (70) Paila, Y. D.; Tiwari, S.; Chattopadhyay, A. Are specific nonannular cholesterol binding sites present in G-protein coupled receptors? Biochim. Biophys. Acta 2009, 1788, 295−302. (71) Lee, J. Y.; Patel, R.; Lyman, E. Ligand-dependent cholesterol interactions with the human A2A adenosine receptor. Chem. Phys. Lipids 2012, 169, 39−45. (72) Lee, J. Y.; Lyman, E. Predictions for Cholesterol Interaction Sites on the A2A Adenosine Receptor. J. Am. Chem. Soc. 2012, 134, 16512−16515.
581
dx.doi.org/10.1021/ci400463z | J. Chem. Inf. Model. 2014, 54, 573−581
