Article Outline

Article Information

PubMed

Related Articles

  • …more

Copyright © 2008 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 95, Issue 4, 1649-1657, 15 August 2008

doi:10.1529/biophysj.107.123190

Biophysical Theory and Modeling

The Interaction of Phospholipase A2 with a Phospholipid Bilayer: Coarse-Grained Molecular Dynamics Simulations

Chze Ling Wee*, 1, Kia Balali-Mood*, 1, David Gavaghan† and Mark S.P. Sansom*, Go To Corresponding Author,  

* Department of Biochemistry, University of Oxford, Oxford, United Kingdom
† Computing Laboratory, University of Oxford, Oxford, United Kingdom

Address reprint requests to Mark S. P. Sansom, Dept. of Biochemistry, University of Oxford, South Parks Road, Oxford OX1 3QU, UK. Tel.: 44-1865-275371; Fax: 44-1865-275373.

1 Chze Ling Wee and Kia Balali-Mood contributed equally to this work.

Abstract

A number of membrane-active enzymes act in a complex environment formed by the interface between a lipid bilayer and bulk water. Although x-ray diffraction studies yield structures of isolated enzyme molecules, a detailed characterization of their interactions with the interface requires a measure of how deeply such a membrane-associated protein penetrates into a lipid bilayer. Here, we apply coarse-grained (CG) molecular dynamics (MD) simulations to probe the interaction of porcine pancreatic phospholipase A2 (PLA2) with a lipid bilayer containing palmitoyl-oleoyl-phosphatidyl choline and palmitoyl-oleoyl-phosphatidyl glycerol molecules. We also used a configuration from a CG-MD trajectory to initiate two atomistic (AT) MD simulations. The results of the CG and AT simulations are evaluated by comparison with available experimental data. The membrane-binding surface of PLA2 consists of a patch of hydrophobic residues surrounded by polar and basic residues. We show this proposed footprint interacts preferentially with the anionic headgroups of the palmitoyl-oleoyl-phosphatidyl glycerol molecules. Thus, both electrostatic and hydrophobic interactions determine the location of PLA2 relative to the bilayer. From a general perspective, this study demonstrates that CG-MD simulations may be used to reveal the orientation and location of a membrane-surface-bound protein relative to a lipid bilayer, which may subsequently be refined by AT-MD simulations to probe more detailed interactions.

Introduction

The interface between a lipid bilayer and water presents a complex environment 1 that is exploited by a number of membrane surface-active enzymes. Although x-ray crystallography yields structures of the isolated enzyme molecules, a detailed model of their interactions with the interface must be approached via spectroscopic 2 and/or computational 3,4 methods. In particular, it often remains uncertain to what extent a membrane-associated protein penetrates into the hydrophobic core of a lipid bilayer.

The phospholipase A2 (PLA2) family provides a well-studied example of membrane surface-active enzymes. PLA2s hydrolyze phospholipids to produce a free fatty acid and a lysophospholipid 5. These enzymes can be broadly classified into two categories: secretory and cytosolic (sPLA2 and cPLA2) 6. sPLA2s are relatively small proteins (<200 residues) that are predominantly α-helical and may be further subdivided into two categories: venom PLA2s and nonvenom PLA2s. Whereas the primary function of venom PLA2s is to promote cell membrane lysis 7, nonvenom PLA2s (e.g., pancreatic PLA2s) act to break down phospholipid compounds in dietary fat. Thus, sPLA2s interact with the extracellular face of cell membranes 8. In contrast, cPLA2s are larger proteins (∼700 residues) that bind to the intracellular leaflet of the cell membrane via a C2 domain 9,10.

The pancreatic sPLA2 has been intensively investigated from a structural 11, biochemical/biophysical 2,3,12,13,14,15,16,17,18,19, and computational 18,20,21 perspective. It therefore provides a good test case for developing computational approaches to understand the interactions of membrane-bound enzymes with the bilayer surface.

Atomistic molecular dynamics (AT-MD) simulations are a widely used computational tool to study the conformational dynamics of integral membrane proteins 22,23 and their interactions with lipids or detergents 24,25. However, application of AT-MD to membrane-associated enzymes may be hampered by uncertainties as to the exact locations of such proteins relative to the membrane surface 4,26,27. Although AT-MD simulations can readily address a timescale in the order of ∼100ns, this may not be sufficient to enable a surface-active enzyme/bilayer system to relax to an equilibrium configuration, especially if more than one species of lipid molecule is present in the bilayer. The more approximate nature of coarse-grained (CG) models 28,29,30,31,32,33,34,35,36,37,38,39,40,41 enables substantially longer simulation times to be accessed and so enables one to optimally locate a surface-active enzyme, either via simulated self-assembly or via relaxation of an initial model. However, the application of CG-MD methods to proteins is still relatively new and has resulted in some discussion 42. It is therefore important to compare the results of CG-MD and AT-MD for well characterized membrane-interacting proteins such as PLA2.

We have used CG-MD simulations from different starting configurations to probe the position and orientation of porcine pancreatic sPLA2 with respect to an overall anionic (3:2 palmitoyl-oleoyl-phosphatidyl choline/palmitoyl-oleoyl-phosphatidyl glycerol (POPC/POPG)) phospholipid bilayer. In total, six CG-MD simulations were performed, each of duration 200ns, corresponding to a total simulation time of 1.2μs. All simulations yielded similar membrane-bound configurations of PLA2, with the enzyme positioned at the headgroup/water interface of the lipid bilayer. To enable comparison with AT-MD simulations, a snapshot from a CG-MD simulation was used as a starting point for two subsequent 20-ns AT-MD simulations. This allowed us to investigate atomic-level interactions, such as hydrogen bonding between the protein and lipid molecules.


Methods

CG simulations

The crystal structure of the porcine pancreatic PLA2 43 (Protein Data Bank code 1P2P, resolution 2.6Å) was converted into a CG representation 44,45. In the CG model, each amino acid is represented by one backbone particle, whose Cartesian coordinates correspond to that of the Cα atom, and between one and three side-chain particle(s) depending on the residue size. The CG model has been shown to be able to reproduce data derived from experiments and atomistic simulations for a number of membrane proteins and synthetic α-helical peptides 45. Transfer free energies of individual amino acids from water to cyclohexane have also been shown to correlate well with experimental data 34. The CG lipid model has also been shown to be able to reproduce structural properties of membranes derived from experiments 46. The secondary and tertiary structures of PLA2 were modeled as an elastic network 47 by imposing a harmonic restraint (force constant of 10kJmol−1Å−2) among all backbone particles that were within 7Å of each other.

Two different CG-MD simulation protocols were used (see Table 1 and Figure 1AB). Both employed an anionic (3:2 POPC/POPG) bilayer, the lipid composition of which was chosen to match experimental conditions 2. In self-assembly simulations (PLA2-SA), the bilayer was self-assembled in the presence of the protein, starting from a box of randomly positioned lipid molecules and CG water particles. In contrast, in the preformed bilayer simulations (PLA2-PF), the initial configuration consists of a PLA2 molecule with its center of mass (c.o.m.) placed ∼30Å from the surface of the lipid bilayer. Self-assembly simulations were performed to test whether the final configuration of the PLA2-PF simulations corresponded to a global or a local energy minimum.

Table 1 Summary of simulations
SystemSetupComponentsDurationFinal RMSD* (Å)
PLA2-SACG, self-assemblyPLA2+154 POPC+102 POPG+4695 water particles3×200ns2.4
PLA2-PFCG, preformed bilayerPLA2+154 POPC+102 POPG + 6417 water particles3×200ns2.3
PLA2-X-rayAT, x-ray structure of PLA2PLA2+154 POPC+102 POPG+16423 water molecules20ns2.2
PLA2-MODAT, CG-template modeled structure of PLA2PLA2+154 POPC+102 POPG+16,354 water molecules20ns2.7
* Further details of the setup of the simulations are provided in the text. The RMSDs were calculated over the final 50ns (CG) or 5ns (AT) of each simulation, relative to the starting (i.e., t=0ns) structure for that simulation, and refer to either the backbone particle (CG simulations) or the Cα atoms (AT simulations).

For the PLA2-SA simulations, the CG enzyme molecule was placed in the center of a simulation box of initial dimensions (122Å)3. A total of 256 lipids (154 POPC and 102 POPG) were added in random positions within the box. The system was solvated (with 4695 CG water particles), and Na+ counterions were added to keep the system electrically neutral. The system was then energy-minimized for 400 steps using a steepest-descent algorithm, followed by a production MD simulation of 200ns. Three repeat simulations were performed, with different initial velocities.

For the PLA2-PF simulations, we started with an equilibrated POPC bilayer with 256 lipids in the center of a simulation box of initial dimensions 93×88×170Å3 with the bilayer normal along the z axis; 102 randomly selected lipids were then converted to POPG. In the CG model, the choline and glycerol moieties are represented as one CG particle; thus, we changed the particle type (charged to nonpolar, i.e., Q0 to Nda), and the charge (+1.0 to 0.0) assigned to this CG particle (see Marrink et al. 32 and Bond and colleagues 44,45 for details of the CG model and parameter set). The system was solvated (with 6417 CG water particles), and Na+ counterions were added. The system was then energy-minimized for 200 steps using a steepest-descent algorithm, and an MD simulation of 10ns was performed to allow the lipids of this mixed bilayer to equilibrate. PLA2 was then introduced into the system, and its c.o.m. was initially positioned ∼30Å from the bilayer surface. The initial orientation of the enzyme was such that the three hydrophobic residues (W3, L19, and M20) that have been proposed to act as hydrophobic anchors 2 were oriented toward the bilayer surface (Figure 1B). A further 200 steps of energy minimization was performed, followed by an MD simulation of 5ns, during which positional restraints (force constant of 10kJmol−1Å−2) were applied to the enzyme and the lipids to allow the solvent to further equilibrate. Three repeat (with different initial velocities) production MD simulations of 200ns were then performed with all positional restraints removed.

CG-MD simulations were performed using GROMACS 3.2.1 48 (www.gromacs.org). The ionizable side chains of the enzyme were in their default charge states for pH 7 (His residues were kept uncharged). Temperature coupling used a Berendsen thermostat 49 with a weak coupling constant of 1.0 ps. The reference temperatures were 323K for PLA2-SA and 310K for PLA2-PF. Electrostatic/Coulombic interactions utilized a relative dielectric constant of 20, which was smoothly shifted to zero between 0 and 12Å. Van der Waal's interactions were smoothly shifted to zero between 9 and 12Å. In PLA2-SA, anisotropic pressure coupling with a Berendsen barostat 49 with a coupling constant of 1.0ps, a compressibility value of 1×10−5bar−1, and a reference pressure of 1 bar was used. In PLA2-PF, semiisotropic pressure coupling in x and y (the z axis of the simulation box corresponds to the bilayer normal) with a coupling constant of 1.0ps, a compressibility value of 5.0×10−6bar−1, and a reference pressure of 1 bar was used. An integration timestep of 0.04ps was used in the CG-MD simulations.


Atomistic simulations

Two AT-MD simulations (see Table 1 and Figure 1C) were performed using a configuration saved from the end of a CG-MD simulation as a starting system. The aim of these was to explore the extent to which the CG-MD simulation had generated an equilibrated system configuration that could be used as a starting point for more detailed atomistic simulations, either from the x-ray structure of the protein or from a model based on the protein structure at the end of the CG simulation. Thus, a snapshot at 200ns of the PLA2-SA CG-MD simulation was stripped of all solvent and ions and was used to define the positions of the lipid molecules and protein in the AT simulations. Two AT simulations were performed (see Table 1): one (PLA2-x-ray) in which the protein coordinates corresponded to the x-ray structure of PLA2 superimposed on the CG structure by least-squares fitting of all atoms/particles and one (PLA2-MOD) in which the CG structure was used as a template for a “remodeled” PLA2 generated using Modeller 50. In each case, the CG POPC/POPG bilayer was used to generate an equivalent atomistic bilayer by superimposing atomistic lipid structures on the corresponding CG models (see Carpenter et al. 51 for details of this procedure). The atomistic lipid molecules were taken from two libraries derived from equilibrated atomistic simulations of pure lipid bilayers, one with 1600 AT POPC molecules and one with 1100 AT POPG molecules. Each system was then solvated with SPC 52 water molecules and electrically neutralized by addition of 102 Na+ counterions. Each system was then energy-minimized for 1000 steps, using a steepest-descent algorithm, followed by an equilibration simulation of 500ps with positional restraints (force constant 10kJmol−1Å−2) applied to waters, enzyme, and lipids. The restraints were removed after 200 ps for waters, after 300ps for enzyme, and after 350ps for lipids.

AT-MD simulations were performed using GROMACS 3.2.1 with the modified GROMOS 53 force field 54. Long-range electrostatic interactions were calculated using the particle mesh Ewald method 55,56, and van der Waals interactions utilized a cutoff of 10Å. Bond lengths and angles were constrained using the LINCS algorithm 57. Semiisotropic pressure coupling with a Parrinello-Rahman barostat 58 with a coupling constant of 5.0ps, a compressibility of 4.5×10−5bar−1, and a reference pressure of 1 bar was used. Each system was temperature-coupled with a Nosé-Hoover thermostat to a reference temperature of 300K with a coupling constant of 1.0ps 59. An integration timestep of 2fs was used.



Results

Progress of the CG-MD simulations

In all three PLA2-PF simulations, the enzyme initially diffuses in the aqueous environment before anchoring into the bilayer at the headgroup/water interface (Figure 1B). The time before the onset of anchoring was ∼125, 50, and 50ns in simulations 1, 2, and 3 of PLA2-PF, respectively. In PLA2-SA (Figure 1A), a bilayer self-assembled around PLA2, and the enzyme eventually adopted an interfacial location by ∼30, 20, and 50ns in simulations 1, 2, and 3, respectively. In all six simulations, the same surface of the protein binds to the membrane. This surface (Fig. 2) consists of a cluster of hydrophobic residues (L2, W3, L19, and M20) surrounded by basic (R6, K10, K116, K121, and K122) and polar residues. W3, L19, and M20 are thought to act as hydrophobic anchors that stabilize the enzyme in a bilayer environment 2. Furthermore, R6, K10, and K116 have been proposed to form hydrogen bonds with the lipid carbonyl oxygens, and K121 and K122 are thought to form electrostatic interactions with the lipid phosphate moieties 2. Inevitably, detailed interactions such as hydrogen bonding are not represented in the CG model (although attempts were made to capture the affinity between hydrogen bonding-capable groups via parameterization of the CG force field 32). However, the CG model provides a basis for initiating detailed AT-MD simulations where such interactions can be probed (hydrogen bonding between the protein and the bilayer is discussed below). After adopting an interfacial location, PLA2 remains at the interface for the remainder of each simulation, consistent with a stable configuration. We focus our subsequent analyses on the residues highlighted in Fig. 2 because they have been proposed experimentally to be important in positioning PLA2 relative to the membrane 2.

Display large version of this figure
Figure 2
CG representation of the membrane-binding surface of PLA2. Hydrophobic, basic, acidic, and polar residues are colored green, blue, red, and white, respectively. Only the CG backbone particles are shown (i.e., side-chain particles are omitted). The broken line highlights a patch of hydrophobic residues (W3, L19, and M20), whose side chains are thought to anchor PLA2 in a bilayer. Basic residues R6, K10, K116, K121, and K122 are thought to form electrostatic interactions with the lipid carbonyl and phosphate moieties.

Depth of penetration

To assess the depth of penetration of PLA2 into the lipid bilayer, we measured the distances between the c.o.m. of the side chains of W3, L9, M20 (the hydrophobic anchor residues), and H48 (at the catalytic site) and the c.o.m. of the bilayer along the bilayer normal (z axis) for both sets of simulations once the equilibrium orientation of PLA2 had been achieved (Table 2). In both systems, the side chains of the hydrophobic anchors are positioned at an average distance of ∼20–24Å from the bilayer center, which corresponds to a location close to the headgroup/tail interface. The side chain of the catalytic site residue (H48) is located at an average distance of ∼30Å from the center, corresponding to a location just above the bilayer surface. There are no significant differences between the corresponding average distances for the PLA2-SA and PLA2-PF simulations, suggesting that the two simulation protocols result in a similar equilibrium orientation of PLA2.

Table 2 Location of key side chains relative to the bilayer center
SystemW3L19M20H48
PLA2-SA19.5±3.920.5±2.123.8±3.029.9±2.2
PLA2-PF20.4±2.122.4±2.323.6±3.029.8±2.8
PLA2-x-ray19.9±0.621.6±0.523.3±0.731.1±0.4
PLA2-MOD20.2±0.621.3±0.623.5±0.831.2±0.6
Average±SD (Å) of the distance between the center of mass of the side-chain particles of the hydrophobic anchors and active site residue of PLA2 and the center of mass of the bilayer. In the CG-MD simulations, averages are taken across all three repeats, from 50 to 200ns for each simulation of PLA2-SA and from 130 to 200ns, 60 to 200ns, and 50 to 200ns for simulations 1, 2, and 3 of PLA2-PF. In the AT-MD simulations, averages were taken over the entire 20ns.

The AT-MD simulations supported the protein location suggested by the CG-MD simulations in terms of the positions of the side chains of the hydrophobic anchor (W3, L19, and M20; Table 2). In the AT-MD simulations, we observed no significant drifts in the positions of the side chains of the hydrophobic anchors along the bilayer normal (<5%) (data not shown) from their initial positions as based on the CG-MD simulations, indicating the stability of the initial configuration, at least over 20 ns. Similarly, there was no significant drift from the initial position of the side chain of H48 in the AT-MD simulations (Table 2).

We can examine the location of PLA2 relative to the lipid bilayer in more detail via the combined time-averaged distribution (projected onto the bilayer normal) of the proposed hydrophobic anchors of PLA2 with respect to those of the components of the bilayer (Fig. 3). It can be seen that W3, L19, and M20 interact with both the lipid headgroup and tail regions and thus allow for the anchoring of the enzyme to the mixed bilayer. Tryptophan residues are known to favor an interfacial location in a bilayer environment 25,60,61,62,63,64,65. The distributions of the AT-MD simulations can be seen to be in good agreement with those from the CG-MD simulations (Fig. 3).

Display large version of this figure
Figure 3
Distribution of PLA2 (black) and its hydrophobic anchors (W3, L19, and M20; blue) with respect to the components of the bilayer system (lipid headgroups and tails; red and green, respectively). The upper panels show distributions for simulation 2 of (A) PLA2-SA (averaged over 50–200ns) and (B) PLA2-PF (averaged over 60–200ns). The distributions of the protein and hydrophobic anchors were scaled by 3× and 10×, respectively, for clarity. The lower panels show the corresponding distributions for the AT-MD simulations of (C) PLA2-x-ray and (D) PLA2-MOD.

In our CG-MD simulations, PLA2 was able to diffuse laterally in the plane of the bilayer, with the lipids dynamically repacking around the enzyme. We estimate a diffusion coefficient of 6 (±3)×10−7cm2/s (± standard error) for both PLA2-SA and PLA2-PF simulations (averaged across the equilibrated period of all three repeats for each protocol). We estimate a lower diffusion coefficient of 4 (±0.3)×10−8cm2/s (± SE) in the AT-MD simulations. This reflects the higher diffusion rate observed in CG-MD simulations. One approach to this is to normalize the timescale in CG-MD simulations 32. We have not done this because preliminary investigations 44 have suggested that different normalization factors would be needed for different systems components (e.g., lipid, protein, or water), as might be expected given the nature of the CG model.


Comparison with experiment

At first glance, our results suggest that PLA2 is positioned a little further away from the bilayer center than previously reported on the basis of depth-dependent quenching of tryptophan fluorescence by dibrominated lipids 2. However, these differences may be accounted for by the complexities of interpretation of such experiments. As has been seen in numerous lipid bilayer simulations (reviewed by, e.g., Feller 66), there is substantial thermal motion exhibited by the hydrophobic tails of a phospholipid in a bilayer environment. In particular, we analyzed the z-distributions of the sites (acyl tail sites 6, 7, 9, 10, 11, and 12) used for dibromo labels in the fluorescence quenching studies 2 in a 10-ns AT-MD simulation of a POPC bilayer (128 lipids; data not shown). As anticipated, this revealed a significant intra- and interlipid spread (both ∼15Å), which increases for sites located further down the lipid acyl chain. Furthermore, there is recent evidence that a lipid bilayer may locally deform/adapt at the protein/lipid interface 3,38,67,68. Thus, allowing for the effect of motional spread of the lipids on the Br atom positions, there may be a degree of error associated with the distances obtained from the fluorescence experiments that is sufficient to encompass the distances obtained from our CG-MD simulations. In addition, the insignificant drifts observed in the detailed AT-MD simulations from the initial configurations provide further support for the predictions made by the CG-MD simulations.


Protein/lipid interactions

Our subsequent analyses focused on a more detailed characterization of the nature of the PLA2/bilayer interactions. In Figure 4A, we show the number of contacts between PLA2 and the POPC and POPG lipids, mapped onto the surface of the protein. It can be seen that the proposed membrane-binding footprint of the enzyme interacts preferentially with the anionic POPG lipids over the zwitterionic POPC lipids, especially if one considers that overall the bilayer has three POPC lipids for every two POPGs. Residues R6, K10, K116, K121, and K122 form frequent contacts with all lipids, but preferentially with the anionic POPGs. A very similar pattern of protein-lipid interactions is seen in the AT-MD simulations (Figure 4B), confirming that no significant shifts in lipid/protein interactions occurred in the course of these (relatively short) simulations for both PLA2-x-ray or PLA2-MOD.

Display large version of this figure
Figure 4
Number of contacts between PLA2 and lipids in the (A) CG-MD and (B) AT-MD simulations. In A, a distance cutoff of 8Å between the outermost side-chain particle and any lipid particles was used to define a contact. We count the total number of contacts from 130 to 200ns, 60 to 200ns, and 50 to 200ns for simulations 1, 2, and 3 of PLA2-PF. We show the total number of contacts separately for POPC and POPG lipids in addition to the total number of contacts with all lipids. A color scale from blue (∼0–100 contacts) to white (∼200–400 contacts) to red (>600 contacts) is used. In B, a distance cutoff of 3.5Å was used to defined a contact between PLA2 and lipid. We count the total number of contacts from 0 to 20 ns of PLA2-x-ray and PLA2-MOD separately. Only the total number of contacts with all lipids is shown. A color scale from blue (∼0–150 contacts) to white (∼250–400 contacts) to red (>650 contacts) is used.

We also measured the time evolution of the combined number of contacts between the hydrophobic anchors (i.e., R3, L19, M20), the basic residues (i.e., R6, K10, K116, K121, K122), and the lipid headgroups and tails separately (Fig. 5). In both simulations, the basic residues are seen to interact preferentially with the headgroups over the tails, with an average (± SD) of 10.7±2.8 and 3.9±2.1 contacts with the headgroups and tails, respectively, in PLA2-PF, and 7.5±3.5 and 2.5±1.9 contacts, respectively, in PLA2-SA. In both PLA2-PF and PLA2-SA, we observe a slight preference for the association of the hydrophobic anchors with the lipid tails, with an average of 5.6±1.5 and 6.1±1.7 contacts with the headgroups and tails, respectively, in PLA2-PF and 4.3±1.6 and 5.4±1.9 contacts, respectively, in PLA2-SA. Overall, our results suggest that both hydrophobic and basic residues play a role binding PLA2 to the bilayer surface but that the latter predominate, explaining the observed preference for PLA2 binding to anionic membranes 2.

Display large version of this figure
Figure 5
Combined number of contacts between the hydrophobic and basic residues of PLA2 and the lipid headgroups and tails versus time. This is calculated from 50 to 200ns and 130 to 200ns of simulation 1 of PLA2-SA and PLA2-PF, respectively. The hydrophobic residues are W3, L19, and M20. The basic residues are R6, K10, K116, K121, and K122. A distance cutoff of 8Å between the outermost side-chain particle and any lipid particles was used to define a contact.

Analysis of protein/lipid hydrogen bonding over time (Fig. 6) in the AT-MD simulations suggests that a total of ∼15–20 hydrogen bonds are formed between PLA2 and the bilayer. Comparison of the two simulations suggests some drift in the number of H-bonds in PLA2-x-ray relative to that in PLA2-MOD, but the difference is not marked. A more detailed examination of the atomistic simulations reveals five to eight hydrogen bonds between R6, K10, and K116 and the lipid carbonyl oxygens (as suggested by Tatulian et al. 2) in both PLA2-x-ray and PLA2-MOD. We also observed interaction of the K121 and K122 side chains with the negatively charged PG headgroups.

Display large version of this figure
Figure 6
Hydrogen bonding between PLA2 and lipids in the AT-MD simulations of (A) PLA2-x-ray and (B) PLA2-MOD. Hydrogen bonds with POPC, POPG, and all lipids shown in red, green, and black, respectively.


Discussion and conclusion

In the current study, we have applied multiple, extended CG-MD simulations, starting from different initial configurations and followed by short (20-ns) AT-MD simulations to determine the equilibrium position of a complex monotopic protein at the surface of a mixed 3:2 POPC/POPG bilayer. Our results reveal a combination of hydrophobic and basic residues that are involved in the binding of the enzyme at the headgroup/water interface of the bilayer. This is in good agreement with recent structural analysis and homology modeling 18. In addition, we show that PLA2 interacts preferentially with the anionic POPG lipids because of the presence of basic residues on its membrane-binding surface. This is seen in both the CG-MD and AT-MD simulations.

This is encouraging in a more general sense in that we show that CG-MD simulations can reveal details of bilayer/protein selectivity (anionic lipids/basic side chains) underlying the interactions of a complex protein with the membrane surface. This is important because of the radical simplification of the treatment of electrostatic interactions in the CG model (see below). Thus, this study extends previous work on CG-MD simulations of protein/lipid interactions for integral membrane proteins 35,38,44,45 and for lipoprotein nanodisks 33,69. It also extends “mean field” approaches 70 to the prediction of the enzyme location in a bilayer that do not allow for detailed protein/lipid interactions.

How do these studies compare with previous AT-MD simulations of proteins interacting with the bilayer surface? The main advance is that CG-MD allows us to proceed in a more “assumption-free” fashion via self-assembly or preformed bilayer association simulations, thus extending previous studies where manual or semimanual docking of protein to a bilayer surface (generally guided by experimental data) was necessary. For example, AT-MD simulations combined with diffraction-derived restraints have been previously used to explore the interfacial location of membrane-active peptides, e.g., mellitin, on a di-oleoyl-phosphatidyl choline bilayer 71. Another AT-MD study used electron paramagnetic resonance data as a starting point in their simulations of the C2 domain of cytosolic PLA2 in a pure POPC bilayer 3. In the current study, we have indicated how one might combine an initial CG-MD simulation (to obtain an equilibrium location of a protein relative to a bilayer) with a subsequent AT-MD simulation (to provide more detailed information on interactions) using the outcome of the CG-MD simulation as a starting point.

What are the technical limitations of the approach we have used? The major limitations reside in the simplifications inherent in the CG force field. Arguably, the most severe limitation resides in the “truncation” of the electrostatic interactions in the CG model. However, this does not prevent the (anticipated) basic side chain/anionic headgroup interaction between PLA2 and the POPG lipids. Although it has been argued that CG simulations are too approximate to allow exploration of basic residue/lipid interactions 42, both the results of the current study (which provides good agreement with experimental data), and the good agreement between CG 68 and AT 67 simulations of the KvAP voltage sensor/lipid interactions 36 suggest this may be too pessimistic an assessment. Furthermore, a detailed comparison with experimental data on lipid/protein interactions for a number of integral membrane proteins suggests good agreement with CG simulations 35. Clearly, however, there are improvements to be made. One is in the direction of refinement of the CG force field 37. The other is in further development of multiscale simulations that allow for mixing 72 or switching between CG and AT representations.

Thus, to the best of our knowledge, this is the first such study of the interactions of a membrane interfacial enzyme with a lipid bilayer using CG-MD simulations. The dual SA/PF CG-MD approach combined with AT-MD provides additional confidence of our results, which are consistent with available experimental data. Thus, this approach can be used as a methodology for predicting the location of interfacial proteins with near atomic accuracy. Using the endpoint of the CG-MD simulations as an initial configuration for more extended (i.e., >100ns) AT-MD simulations should enable future studies to explore possible conformational transitions in PLA2 resulting from its initial encounter with a lipid bilayer. Furthermore, the computational efficiency of CG-MD simulations means it is now feasible to perform free energy calculations to define rigorously the optimum location of PLA2 in a bilayer environment, as has recently been performed for a toxin/bilayer system 73.


Acknowledgments

Thanks to Peter Bond and Henry Chew for their support and assistance.

Research in M.S.P.S.'s laboratory is supported by the Biotechnology and Biological Sciences Research Council, Engineering and Physical Sciences Research Council, and the Wellcome Trust. C.L.W. acknowledges support from the Overseas Research Students Awards Scheme.

References

1. White, S.H., and Wimley, W.C. (1999). Membrane protein folding and stability: physical principles. Annu. Rev. Biophys. Biomol. Struct. 28, 319–365. PubMed

2. Tatulian, S.A., Qin, S., Pande, A.H., and He, X. (2005). Positioning membrane proteins by novel protein engineering and biophysical approaches. J. Mol. Biol. 351, 939–947. PubMed

3. Jaud, S., Tobias, D.J., Falke, J.J., and White, S.H. (2007). Self-induced docking site of a deeply embedded peripheral membrane protein. Biophys. J. 92, 517–524. Abstract | Full Text | PDF (1056 kb) | PubMed

4. Fowler, P.W., and Coveney, P.V. (2006). A computational protocol for the integration of the monotopic protein prostaglandin H2 synthase into a phospholipid bilayer. Biophys. J. 91, 401–410. Abstract | Full Text | PDF (711 kb) | PubMed

5. Nigam, S., Shimizub, T., and Wilton, D.C. (2006). Special issue on phospholipase A2. Biochim. Biophys. Acta 1761, 1245–1400. PubMed

6. Kudo, I., and Murakami, M. (2002). Phospholipase A2 enzymes. Prostaglandins Other Lipid Mediat 68–69, 3–58. PubMed

7. Ownby, C.L., Powell, J.R., Jiang, M.S., and Fletcher, J.E. (1997). Melittin and phospholipase A2 from bee (Apis mellifera) venom cause necrosis of murine skeletal muscle in vivo. Toxicon 35, 67–80. PubMed

8. Balestrieri, B., and Arm, J.P. (2006). Group V sPLA2: classical and novel functions. Biochim. Biophys. Acta 1761, 1280–1288. PubMed

9. Cho, W., and Stahelin, R.V. (2005). Membrane-protein interactions in cell signaling and membrane trafficking. Annu. Rev. Biophys. Biomol. Struct. 34, 119–151. PubMed

10. Bhardwaj, N., Stahelin, R.V., Zhao, G., Cho, W., and Lui, H. (2007). MeTaDoR: A comprehensive resource for membrane targeting domains and their host proteins. Bioinformatics 23, 3110–3112. PubMed

11. Sekar, K. (2007). Structural biology of recombinant bovine pancreatic phospholipase A2 and its inhibitor complexes. Curr. Top. Med. Chem. 7, 779–785. PubMed

12. Verheij, H.M., Volwerk, J.J., Jansen, E.H., Puyk, W.C., Dijkstra, B.W., Drenth, J., and De Haas, G.H. (1980). Methylation of histidine-48 in pancreatic phospholipase A2. Role of histidine and calcium ion in the catalytic mechanism. Biochemistry 19, 743–750. PubMed

13. Tatulian, S.A., Biltonen, R.L., and Tamm, L.K. (1997). Structural changes in a secretory phospholipase A2 induced by membrane binding: a clue to interfacial activation?. J. Mol. Biol. 268, 809–815. PubMed

14. Meulenhoff, P.J. (1999). Interfacial action of phospholipase A2 [PhD dissertation]. (Groningen: University of Groningen), 132 p. PubMed

15. Tatulian, S.A. (2001). Toward understanding interfacial activation of secretory phospholipase A2 (PLA2): membrane surface properties and membrane-induced structural changes in the enzyme contribute synergistically to PLA2 activation. Biophys. J. 80, 789–800. Abstract | Full Text | PDF (309 kb) | PubMed

16. Tsai, Y.C., Yu, B.Z., Wang, Y.Z., Chen, J.H., and Jain, M.K. (2006). Desolvation map of the i-face of phospholipase A2. Biochim. Biophys. Acta 1758, 653–665. PubMed

17. Jain, M.K., and Berg, O.G. (2006). Coupling of the i-face and the active site of phospholipase A(2) for interfacial activation. Curr. Opin. Chem. Biol. 10, 473–479. PubMed

18. Winget, J.M., Pan, Y.H., and Bahnson, B.J. (2006). The interfacial binding surface of phospholipase A2s. Biochim. Biophys. Acta 1761, 1260–1269. PubMed

19. Yu, B.Z., Apitz-Castro, R.J., Jain, M.K., and Berg, O.G. (2007). Role of 57–72 loop in the allosteric action of bile salts on pancreatic IB phospholipase A(2): Regulation of fat and cholesterol homeostasis. Biochim. Biophys. Acta 1768, 2478–2490. PubMed

20. Jones, S.T., Ahlstrom, P., Berendsen, H.J., and Pickersgill, R.W. (1993). Molecular dynamics simulation of a phospholipase A2-substrate complex. Biochim. Biophys. Acta 1162, 135–142. PubMed

21. Zhou, F., and Schulten, K. (1996). Molecular dynamics study of phospholipase A2 on a membrane surface. Proteins Struct. Funct. Genet. 25, 12–27. PubMed

22. Ash, W.L., Zlomislic, M.R., Oloo, E.O., and Tieleman, D.P. (2004). Computer simulations of membrane proteins. Biochim. Biophys. Acta 1666, 158–189. PubMed

23. Gumbart, J., Wang, Y., Aksimentiev, A., Tajkhorshid, E., and Schulten, K. (2005). Molecular dynamics simulations of proteins in lipid bilayers. Curr. Opin. Struct. Biol. 15, 423–431. PubMed

24. Bond, P.J., and Sansom, M.S.P. (2003). Membrane protein dynamics vs. environment: simulations of OmpA in a micelle and in a bilayer. J. Mol. Biol. 329, 1035–1053. PubMed

25. Sansom, M.S.P., Bond, P.J., Deol, S.D., Grottesi, A., Haider, S., and Sands, Z.A. (2005). Molecular simulations and lipid/protein interactions: potassium channels and other membrane proteins. Biochem. Soc. Trans. 33, 916–920. PubMed

26. Nina, M., Bernèche, S., and Roux, B. (2000). Anchoring of a monotopic membrane protein: the binding of prostaglandin H2 synthase-1 to the surface of a phospholipid bilayer. Eur. Biophys. J. 29, 439–454. PubMed

27. Fowler, P.W., Balali-Mood, K., Deol, S., Coveney, P.V., and Sansom, M.S.P. (2007). Monotopic enzymes and lipid bilayers: a comparative study. Biochemistry 46, 3108–3115. PubMed

28. Shelley, J.C., Shelley, M.Y., Reeder, R.C., Bandyopadhyay, S., and Klein, M.L. (2001). A coarse grain model for phospholipid simulations. J. Phys. Chem. B 105, 4464–4470. PubMed

29. Nielsen, S.O., Lopez, C.F., Srinivas, G., and Klein, M.L. (2004). Coarse grain models and the computer simulation of soft materials. J. Phys. Condens. Matter 16, R481–R512. PubMed

30. Nielsen, S.O., Lopez, C.F., Ivanov, I., Moore, P.B., Shelley, J.C., and Klein, M.L. (2004). Transmembrane peptide-induced lipid sorting and mechanism of L-alpha-to-inverted phase transition using coarse-grain molecular dynamics. Biophys. J. 87, 2107–2115. Abstract | Full Text | PDF (544 kb) | PubMed

31. Murtola, T., Falck, E., Patra, M., Karttunen, M., and Vattulainen, I. (2004). Coarse-grained model for phospholipid/cholesterol bilayer. J. Chem. Phys. 121, 9156–9165. PubMed

32. Marrink, S.J., and H. de Vries, S Mark, A.E. (2004). Coarse grained model for semiquantitative lipid simulations. J. Phys. Chem. B 108, 750–760. PubMed

33. Shih, A.Y., Arkhipov, A., Freddolino, P.L., and Schulten, K. (2006). Coarse grained protein-lipid model with application to lipoprotein particles. J. Phys. Chem. B 110, 3674–3684. PubMed

34. Bond, P. J., C. L. Wee, and M. S. P. Sansom. 2008. Coarse-grained molecular dynamics simulations as an approach to the energetics of helix insertion into a bilayer. Biochemistry. Submitted..

35. Scott, K.A., Bond, P.J., Ivetac, A., Chetwynd, A.P., Khalid, S., and Sansom, M.S.P. (2008). Coarse-grained MD simulations of membrane protein/bilayer self assembly. Structure 16, 621–630. Abstract | Full Text | PDF (1678 kb) | PubMed

36. Sansom, M.S.P., Scott, K.A., and Bond, P.J. (2008). Coarse grained simulation: a high throughput computational approach to membrane proteins. Biochem. Soc. Trans. 36, 27–32. PubMed

37. Marrink, S.J., Risselada, J., Yefimov, S., Tieleman, D.P., and De Vries, A.H. (2007). The MARTINI forcefield: coarse grained model for biomolecular simulations. J. Phys. Chem. B 111, 7812–7824. PubMed

38. Periole, X., Huber, T., Marrink, S.J., and Sakmar, T.P. (2007). G protein-coupled receptors self-assemble in dynamics simulations of model bilayers. J. Am. Chem. Soc. 129, 10126–10132. PubMed

39. Treptow, W., Marrink, S.-J., and Tarek, M. (2008). Gating motions in voltage-gated potassium channels revealed by coarse-grained molecular dynamics simulations. J. Phys. Chem. B 112, 3277–3282. PubMed

40. Yefimov, S., van der Giessen, E., Onck, P.R., and Marrink, S.J. (2008). Mechanosensitive membrane channels in action. Biophys. J. 94, 2994–3002. Abstract | Full Text | PDF (846 kb) | PubMed

41. Monticelli, L., Kandasamy, S.K., Periole, X., Larson, R.G., Tieleman, D.P., and Marrink, S.J. (2008). The MARTINI coarse grained force field: extension to proteins. J. Chem. Theor. Comput. 4, 819–834. PubMed

42. Allen, T.W. (2007). Modeling charged protein side chains in lipid membranes. J. Gen. Physiol. 130, 237–240. PubMed

43. Dijkstra, B.W., Renetseder, R., Kalk, K.H., Hol, W.G., and Drenth, J. (1983). Structure of porcine pancreatic phospholipase A2 at 2.6 A resolution and comparison with bovine phospholipase A2. J. Mol. Biol. 168, 163–179. PubMed

44. Bond, P.J., and Sansom, M.S.P. (2006). Insertion and assembly of membrane proteins via simulation. J. Am. Chem. Soc. 128, 2697–2704. PubMed

45. Bond, P.J., Holyoake, J., Ivetac, A., Khalid, S., and Sansom, M.S.P. (2007). Coarse-grained molecular dynamics simulations of membrane proteins and peptides. J. Struct. Biol. 157, 593–605. PubMed

46. de Vries, A.H., Yefimov, S., Mark, A.E., and Marrink, S.J. (2005). Molecular structure of the lecithin ripple phase. Proc. Natl. Acad. Sci. USA 102, 5392–5396. PubMed

47. Atilgan, A.R., Durell, S.R., Jernigan, R.L., Demirel, M.C., Keskin, O., and Bahar, I. (2001). Anisotropy of fluctuation dynamics of proteins with an elastic network model. Biophys. J. 80, 505–515. Abstract | Full Text | PDF (664 kb) | PubMed

48. van der Spoel, D., Lindahl, E., Hess, B., Groenhof, G., Mark, A.E., and Berendsen, H.J. (2005). GROMACS: fast, flexible, and free. J. Comput. Chem. 26, 1701–1718. PubMed

49. Berendsen, H.J.C., Postma, J.P.M., van Gunsteren, W.F., DiNola, A., and Haak, J.R. (1984). Molecular dynamics with coupling to an external bath. J. Chem. Phys. 81, 3684–3690. PubMed

50. Sali, A., and Blundell, T.L. (1993). Comparative protein modeling by satisfaction of spatial restraints. J. Mol. Biol. 234, 779–815. PubMed

51. Carpenter, T., P. J. Bond, S. Khalid, and M. S. P. Sansom. 2008. Self-assembly of a simple membrane protein: coarse-grained molecular dynamics simulations of the influenza M2 channel. Biophys. J. Submitted..

52. Berendsen, H.J.C., Postma, J.P.M., van Gunsteren, W.F., and Hermans, J. (1981). Interactions models for water in relation to protein hydration. In Intermolecular Forces. Pullman, B., ed. (Dordrecht: Reidel), pp. 331–342. PubMed

53. Scott, W.R.P., Hunenberger, P.H., Tironi, I.G., Mark, A.E., Billeter, S.R., Fennen, J., Torda, A.E., Huber, T., Kruger, P., and van Gunsteren, W.F. (1999). The GROMOS biomolecular simulation program package. J. Phys. Chem. A 103, 3596–3607. PubMed

54. Berger, O., Edholm, O., and Jahnig, F. (1997). Molecular dynamics simulations of a fluid bilayer of dipalmitoylphosphatidycholine at full hydration, constant pressure and constant temperature. Biophys. J. 72, 2002–2013. Abstract | | PubMed

55. Darden, T., York, D., and Pedersen, L. (1993). Particle mesh Ewald—an N·log(N) method for Ewald sums in large systems. J. Chem. Phys. 98, 10089–10092. PubMed

56. Essmann, U., Perera, L., Berkowitz, M.L., Darden, T., Lee, H., and Pedersen, L.G. (1995). A smooth particle mesh Ewald method. J. Chem. Phys. 103, 8577–8593. PubMed

57. Hess, B., Bekker, H., Berendsen, H.J.C., and Fraaije, J.G.E.M. (1997). LINCS: A linear constraint solver for molecular simulations. J. Comput. Chem. 18, 1463–1472. PubMed

58. Parrinello, M., and Rahman, A. (1981). Polymorphic transitions in single crystals—a new molecular-dynamics method. J. Appl. Phys. 52, 7182–7190. PubMed

59. Hoover, W.G. (1985). Canonical dynamics: equilibrium phase-space distributions. Phys. Rev. A 31, 1695–1697. PubMed

60. Schiffer, M., Chang, C.H., and Stevens, F.J. (1992). The functions of tryptophan residues in membrane proteins. Prot. Engng. 5, 213–214. PubMed

61. Yau, W.M., Wimley, W.C., Gawrisch, K., and White, S.H. (1998). The preference of tryptophan for membrane interfaces. Biochemistry 37, 14713–14718. PubMed

62. Ridder, A., Morein, S., Stam, J.G., Kuhn, A., De Kruijff, B., and Killian, J.A. (2000). Analysis of the role of interfacial tryptophan residues in controlling the topology of membrane proteins. Biochemistry 39, 6521–6528. PubMed

63. Killian, J.A., and von Heijne, G. (2000). How proteins adapt to a membrane-water interface. Trends Biochem. Sci. 25, 429–434. Abstract | Full Text | PDF (854 kb) | PubMed

64. Grossfield, A., and Woolf, T.B. (2002). Interaction of tryptophan analogs with POPC lipid bilayers investigated by molecular dynamics calculations. Langmuir 18, 198–210. PubMed

65. Norman, K.E., and Nymeyer, H. (2006). Indole localization in lipid membranes revealed by molecular simulation. Biophys. J. 91, 2046–2054. Abstract | Full Text | PDF (508 kb) | PubMed

66. Feller, S.E. (2000). Molecular dynamics simulations of lipid bilayers. Curr. Opin. Colloid Interface Sci. 5, 217–223. PubMed

67. Sands, Z.A., and Sansom, M.S.P. (2007). How does a voltage-sensor interact with a lipid bilayer? Simulations of a potassium channel domain. Structure 15, 235–244. Abstract | Full Text | PDF (1527 kb) | PubMed

68. Bond, P.J., and Sansom, M.S.P. (2007). Bilayer deformation by the Kv channel voltage sensor domain revealed by self-assembly simulations. Proc. Natl. Acad. Sci. USA 104, 2631–2636. PubMed

69. Lu, D.Y., Aksimentiev, A., Shih, A.Y., Cruz-Chu, E., Freddolino, P.L., Arkipov, A., and Schulten, K. (2006). The role of molecular modelling in bionanotechnology. Phys. Biol. 3, S40–S53. PubMed

70. Lomize, A.L., Pogozheva, I.D., Lomize, M.A., and Mosberg, H.I. (2006). Positioning of proteins in membranes: A computational approach. Protein Sci. 15, 1318–1333. PubMed