17554
Biochemistry 1998, 37, 17554-17561
Lipid Properties and the Orientation of Aromatic Residues in OmpF, Influenza M2, and Alamethicin Systems: Molecular Dynamics Simulations† D. Peter Tieleman,*,‡ Lucy R. Forrest,§ Mark S. P. Sansom,§ and Herman J. C. Berendsen‡ BIOSON Research Institute and Laboratory of Biophysical Chemistry, UniVersity of Groningen, Nijenborgh 4, 9747 AG Groningen, The Netherlands, and Laboratory of Molecular Biophysics, Department of Biochemistry, The Rex Richards Building, UniVersity of Oxford, South Parks Road, Oxford, OX1 3QU, United Kingdom ReceiVed July 27, 1998; ReVised Manuscript ReceiVed September 17, 1998
ABSTRACT: Molecular dynamics simulations allow a direct study of the structure and dynamics of membrane proteins and lipids. We describe the behavior of aromatic residues and lipid properties in POPE and POPC bilayer models with the Escherichia coli OmpF trimer, single alamethicin and Influenza M2 helices, 4-helix M2 bundles, and two alamethicin 6-helix channel models. The total simulation time is over 24 ns, of systems containing solvent, protein, and between 104 and 318 lipids. Various types of adjustment between lipids and proteins occur, depending on the size of the protein and the degree of hydrophobic mismatch between lipid and protein. Single helices cause little measurable effect on nearby lipids whereas the 4-helix bundles, 6-helix channel models, and OmpF cause a significant lowering of order parameters in nearby lipid chains, an increased difference between odd and even chain dihedrals in the magnitude of the trans dihedral fractions and dihedral transition rates, and in most cases a decreased gauche population and a decrease in bilayer thickness. An increased tilt of the lipid chains near the proteins can account for most of the observed decrease in order parameters. The orientation of tryptophans and tyrosines on the outside of the proteins is determined by packing at the protein exterior and non-specific hydrogen bonding with lipids and solvent. The tyrosines in the broad bands that delimit the hydrophobic exterior of OmpF show little change in orientation over one nanosecond. Their rings are oriented predominantly perpendicular to the bilayer plane, with the hydroxyl group pointing toward the lipid-water interface. Phenylalanines in OmpF, alamethicin, and Influenza M2 are more mobile and assume a variety of orientations.
Lipid-protein interactions play an important role in membrane protein folding and assembly, partitioning, aggregation, and other processes, but it is difficult to obtain detailed structural information on such interactions. Experimentally, spectroscopic techniques, foremost NMR, are the most powerful tools to study the behavior of lipids near proteins and the structure of membrane-bound or inserted peptides. NMR is able to obtain information about the structure of membrane proteins and labeled lipids, although this kind of study often still presents a significant technical challenge (1-3). From the theoretical side, molecular dynamics simulations (4, 5) of lipid systems have evolved to a point where complicated systems including bacteriorhodopsin (6), gramicidin A (7), and the bacterial porin OmpF (8) have been studied in atomic detail. Simulations give detailed information about the structure and dynamics of lipids and protein side chains. This makes it possible to study important features of membrane proteins, in particular the orientation of † DPT was supported by the European Union under Contract CT940124. L.F. is an MRC student. M.S.P.S. thanks the Wellcome Trust for research support and the Oxford Supercomputing Centre for access to its facilities. * Author to whom correspondence should be addressed. E-mail:
[email protected]. Tel: +31 50 3634338. Fax: +31 50 3634800. ‡ University of Groningen. § University of Oxford.
Table 1: Overview of the Simulations system
protein
Alm 6 Alm 6H Alm 1 Alm 1H Flu 18 Flu 26 Flu 34 Flu T18 Flu T22 OmpF POPC
6 Alm 6 Alm 1 Alm 1 Alm 1 Flu M2 1 Flu M2 1 Flu M2 4 Flu M2 4 Flu M2 3 OmpF
water
ions
104 POPC 3528 104 POPC 3528 127 POPC 3822 127 POPC 3822 127 POPC 3790 127 POPC 4712 127 POPC 3803 110 POPC 4860 110 POPC 5236 318 POPE 12992 128 POPC 2460
lipid
6 Na+ 1 Na+ 1 Na+ 2 Na+ 1 Na+ 27 Na+
length (ns) ref 2 2 2 2 2 2 2 4 4 1 1.5
22 22 21 21 23 23 23 24 24 8
aromatic residues and the influence of membrane proteins on lipid behavior. Recently, several simulation studies have addressed the structure and dynamics of small proteins in bilayers, including gramicidin A (7), helices from bacteriorhodopsin (9, 10), a polyalanine model helix (11), and an amphiphilic model helix (12). However, it is difficult to obtain reasonable statistics on lipid properties as a function of the distance from a protein since relatively long simulations with a large number of lipids are necessary. Here, we analyze proteinlipid interactions in a number of simulations of membrane channels and transmembrane helices (Table 1). Combined, they possess a variety of motifs common in membrane proteins and form a unique set of data.
10.1021/bi981802y CCC: $15.00 © 1998 American Chemical Society Published on Web 11/26/1998
Lipid-Protein Interactions Studied by MD
Biochemistry, Vol. 37, No. 50, 1998 17555
Table 2: Sequences of the Peptidesa Alm Flu 18 Flu 26 Flu 34
Aib-Pro-Aib-Ala-Aib-Ala-Gln-Aib-Val-Aib-Gly-Leu-AibPro-Val-Aib-Aib-Glu-Gln-Phe Leu-Val-Ile-Ala-Ala-Ser-Ile-Ile-Gly-Ile-Leu-His-Phe-IleLeu-Trp-Ile-Leu Ser-Ser-Asp-Pro- Flu 18 -Asp-Arg-Leu-Phe Ser-Cys-Ser-Asp- Flu 26 -Phe-Lys-Cys-Ile
a Flu M2 is a 97 residues protein, only part of which is simulated. The three different lengths are different putative transmembrane segments, but that is, in the current context, of no real importance (24). Flu T22 contains 4 helices with 22 residues, with the same sequence as Flu 26 but without the first and last two residues.
Several studies have examined the behavior of aromatic residues in systems such as gramicidin A (13-15), porins (16, 17), and model peptides (18). We describe the orientation and dynamics of aromatic residues, emphasizing the differences between Phe, Trp, and Tyr. Alamethicin (Alm) contains a Phe, Influenza M2 (Flu M2) contains one or more Phe residues and a Trp, and OmpF contains a wide hydrophobic band delimited by two rings of aromatic residues. We also study how the membrane proteins influence the order parameters, chain dynamics and conformations, bilayer thickness, and lipid tilt with respect to the bilayer normal. These properties are of theoretical interest in thermodynamic models describing protein insertion and aggregation and form part of the molecular basis of the often used concept of hydrophobic mismatch (18-20). The influence of the size of the membrane protein on the properties of the neighboring lipids can be compared between the small single helices, medium-sized helix bundles, and the large OmpF trimer. The effect of different degrees of hydrophobic mismatch between proteins and lipid can be compared between Alm, the Flu M2 helices and bundles, and OmpF, which have different hydrophobic lengths. 1 METHODS
FIGURE 1: Snapshots of all proteins. Hydrophobic residues are light, hydrophilic dark, and all aromatic residues are drawn explicitly. Only one Alm helix and bundle is drawn; the other two are almost the same.
1.1 Simulations All of systems that we study here have been described in detail elsewhere (8, 21-24). An overview of the simulations is given in Table 1. In Table 2 the sequence of the peptides is given and in Figure 1 their structure. Briefly, OmpF has been simulated in a POPE bilayer with lipid force field parameters based on GROMOS (25). The other systems use POPC lipids and lipid parameters based on ref 26. They give a better density in pure lipid systems but were not available when the simulations on OmpF started. POPE has a higher melting temperature, a smaller area per lipid, and considerably more ordered chains than POPC at the same temperature. In all systems GROMOS parameters are used for the proteins (27). All systems use a group-based twin-range cutoff of 1.0 nm for the van der Waals and 1.8 nm for the Coulomb interactions. Neighborlists were updated every 15 steps (OmpF) or every 10 steps (all others). All systems were simulated at a constant pressure of 1.0 bar and a constant temperature of 315 K (OmpF) or 300 K (all others), using the weak coupling method (τp ) 1.0 ps, τT ) 0.1 ps) (28). All calculations were carried out using the GROMACS package (29, 30). Structures were saved every 100 steps (0.2
ps, OmpF) or 500 steps (1.0 ps, all others) and used for analysis. 1.2 Analyses The numbering of carbon atoms in the palmitoyl chain runs from 1 (the carbonyl carbon) to 16 (the last tail atom). The distances between lipids and protein are calculated as the minimum distance between a lipid carbon (order parameters), bond (trans fractions and transition rates), or tail carbon atom 5 (lipid tilt) and any protein CR atom. This is an arbitrary choice. A binwidth of 0.8 nm (roughly equal to the distance between two lipids) proved to be a useful compromise between sufficient statistics and good resolution. The first point in graphs with distance-dependent properties can therefore be interpreted roughly as the first layer of lipids around a protein. Order parameters, chain dihedral properties, and chain tilt are calculated for the palmitoyl chain only. Deuterium order parameters were calculated in the usual fashion (5). For the trans/gauche and dihedral transition analyses, dihedrals were divided into three bins of 60° width, around the gauche+, gauche-, and trans minima. The trans fraction
17556 Biochemistry, Vol. 37, No. 50, 1998
FIGURE 2: A graphical view of the three examples. The dashed line represents the plane of the membrane. The three thick lines with arrows represent the normal to that plane, or the z axis. The long axis in all three snapshots is dashed with an arrow, the normal on the ring thin solid with an arrow. In the left snapshot, SL ) -0.5 and SN ) 1. SN can assume all values by rotating around the Cβ-Cγ bond, which does not change SL. The angle θ is the angle between the long axis and the z axis. In the middle snapshot, the long axis is aligned with the z axis and SL has the value 1. This means that SN can only be -0.5 and is perpendicular to the z axis. Rotation around the long axis would change neither SL nor SN. In the right snapshot, θ shows the angle between the normal to the ring and the z axis. The long axis is at some angle with the z axis, and SN can assume a range of values by rotation around the CβCγ bond. However, the maximum value for SN is limited by the orientation of the long axis.
was calculated by dividing the total number of dihedrals in the trans bin by the total number of dihedrals in the three bins. If a dihedral fell outside one of the three bins it was not counted, nor was its bin changed. A transition was counted if a dihedral changed from one bin to another. The tilt of lipids was calculated as the angle between the vector from chain carbon 1 (the carbonyl carbon) to 10 and the z axis. The distance of a whole lipid from the protein was calculated as the minimum distance between carbon atom 5 and any CR atom. Membrane thickness was calculated as the distance between the carbon 2 atoms of lipids of the two monolayers, which corresponds roughly to the hydrophobic interior of the bilayer. In all analyses, standard errors were calculated assuming that each lipid is an independent sample. This is likely to underestimate the errors slightly because there is a weak correlation between the lipids, but this is compensated for by the extra information due to the long length of the simulations. The orientation of aromatic residues is described in terms of two order parameters. SN is defined as SN ) 1/2(3 cos2 θ - 1), with θ the angle between the normal vector to the plane of the phenylalanine, tryptophan, and tyrosine ring, and the axis perpendicular to the membrane. SL is defined in the same way but with θ the angle between the vector from Cγ to Cζ for tyrosine and phenylalanine and the z axis. Figure 2 gives a more detailed graphical explanation of this. A hydrogen bond is counted if the donor-acceptor distance is less than 0.35 nm and the angle between hydrogen-donor-acceptor is less than 60°. 2 RESULTS 2.1 Aromatic Residues The different systems contain several types of aromatic residues, but not all of them are in contact with the membrane. Alm has a C-terminal phenylalanine, which is located at the edge of the hydrophobic membrane interior.
Tieleman et al. The Flu peptides contain a phenylalanine, which is located in the membrane interior, and a tryptophan, which is located near the acyl chain-headgroup interface. In the 4-helix bundles this residue is buried inside the bundle. Flu 26 and Flu 34 contain additional phenylalanines, but these are located in the headgroup region and are not further analyzed. The most interesting protein is the porin trimer. Each monomer contains two tryptophans, 19 phenylalanines, and 29 tyrosines. Trp61 is buried near the trimer axis and is completely immobile in the simulation, 9 phenylalanines are buried in the barrel interior or located in an external loop, and 16 tyrosines are part of the porelining, buried in the barrel interior or located in loops. This leaves 3 tryptophans, 30 phenylalanines, and 39 tyrosines in the trimer that interact directly with the lipid bilayer. Two parameters describe the orientation of an aromatic residue in the bilayer. SN describes the angle between the z axis, which is the axis normal to the bilayer, and a normal vector on the plane of the aromatic side chain. A value of 1 means that this vector is parallel to the z axis, and thus, the ring itself is parallel to the plane of the membrane. A value of -0.5 means that this vector is perpendicular to the z axis, and thus, the ring is perpendicular to the plane of the membrane. SL describes the angle between the z axis and a vector from the Cβ to Cζ in tyrosine and phenylalanine. Three examples and Figure 2 may clarify this. (1) If a side chain has SN ) 1, the plane of the ring is parallel to the plane of the membrane and SL can only be perpendicular to the z axis, SL ) -0.5. (2) If SL ) 1, the long axis of the side chain is parallel to the z axis and SN can only be -0.5, perpendicular to the z axis. It is still possible that the plane of the ring rotates around the long axis, but this is not visible in either SN or SL. (3) If SL goes from -0.5 to close to 1, SN changes automatically, too. However, there is an additional motion, namely, rotation of the plane of the ring around the long axis, that changes SN. Thus, for most orientations of SL, SN can assume almost the entire range from -0.5 to 1 by rotation around the long axis. In Figure 3 a few representative cases are plotted. Below we examine each system in more detail. 2.1.1 Alamethicin. The phenylalanines in the four alamethicin systems have considerable freedom. Figure 3A shows that in the single helix systems many possible combinations occur. When SN is high and the phenylalanine ring is oriented parallel to the bilayer plane, the side chains are extended away from the helix, perpendicular to the z axis. When SN is lower, more conformations are accessible for the long axis. In 2 ns major changes in conformation take place. A similar picture is seen for the Alm 6 and Alm 6H bundles, with a wide distribution of orientations, including orientations with the long axis more parallel to the bundle (not shown). 2.1.2 Influenza M2. Similar to the case in alamethicin, the phenylalanine in the Flu monomers remains mostly in a conformation with the side chain extended away from the helix, almost perpendicular to the helix. Only a few short fluctuations in this orientation occur. The orientation of the ring itself fluctuates rapidly by rotation around the long axis, with fast transitions (tens of picoseconds) between a parallel and perpendicular orientation of the plane of the ring with respect to the plane of the membrane. In the short tetramer, Flu T18, the long axis is also extended and hardly changes orientation in 4 ns. The orientation of the plane of the ring
Lipid-Protein Interactions Studied by MD
FIGURE 3: The orientation of selected aromatic residues as a function of time, defined by SN (left column) and SL (right column) (see text). (A) Phe20 from Alm 1 (black) and Alm 1H (light gray), as examples of extended freely rotating phenylalanines. (B) Three examples of phenylalanine from OmpF. Phe23 (black), part of the aromatic band in OmpF, remains oriented with its ring perpendicular to the membrane plane, but its long axis is free to move. Phe145 (dark gray) is stacked between other side chains and is mostly immobile. Phe153 (light gray) is located in the middle of the hydrophobic band and is not restricted. (C) Three examples of tyrosine from OmpF. Tyr90 (black) is the most flexible tyrosine on the protein exterior. Tyr98 (dark gray) is buried in the protein exterior and completely immobilized. Tyr180 (light gray) is a typical example of the 24 tyrosines that form the aromatic band on either side of the protein.
fluctuates much more, with several sharp transitions and a wide spread in orientations between the four phenylalanines during the simulation. In Flu T22, two of the phenylalanines orient with their long axis parallel to the helix, effectively folding the side chains against the helix. This completely restricts the orientational freedom for the ring, and consequently, only very small motions are seen. The other two sidegroups are extended away from the helix, and their rings fluctuate over all available orientations. Autocorrelation times for the ring orientation can in principle be calculated, but these yield widely varying results, from picoseconds to well beyond the longest simulation time of 4 ns. In Flu 18 and Flu 26, the orientation of the Trp plane fluctuates around a mostly perpendicular orientation with respect to the plane of the membrane. The NH group in Flu 18 forms hydrogen bonds with the carbonyl oxygens of His13, Phe14, or one of two different lipids. In Flu 26 the NH group is not hydrogen-bonded, and in Flu 34 the NH group is bonded to water, lipid, or Asp28. The orientation of the Trp residues in the short bundle model is difficult to interpret, since they are pointed toward each other, at the inside and top of the bundle. In the longer bundle they are approximately stacked on top of each other on the inside of the helix bundle and are not accessible to lipids. 2.1.3 OmpF. Phe23 (Figure 3B) is part of the aromatic boundary of the hydrophobic band on the extracellular side. The plane of the ring in all three monomers remains mostly
Biochemistry, Vol. 37, No. 50, 1998 17557 perpendicular to the plane of the membrane, but the long axis has many possible orientations. Phe96 is buried in the middle of the band in the cleft on the outside between two monomers. This renders the residue nearly immobile during 1 ns, with the plane of the ring perpendicular to the plane of the membrane and the long axis mostly parallel to the z axis. Phe144 and Phe145 (Figure 3B) are stacked together and are part of the aromatic boundary on the intracellular side. Phe153 neighbors them but is part of the interior of the hydrophobic band. It is free to move and assumes a wide range of orientations (Figure 3B). Phe185 forms part of the aromatic boundary on the intracellular side and is sandwiched between Tyr180 and Tyr220. Phe265 and Phe267 are stacked and sandwiched between Tyr263 and Tyr301 at the intracellular side. Phe295 is part of the aromatic boundary on the extracellular side, and Phe303 on the intracellular side. Most of the tyrosines on the outside of OmpF are ordered in two regular bands with roughly enough space between neighboring tyrosines for another large side chain, often Phe or Trp: 313, 275, 231, and 191 on the extracellular side and 301, 263, 220, and 180 on the intracellular side. These residues are placed at regular 12-residue intervals on the β-sheets that make up the outside of the protein. If we look at the orientation of the rings and the long axis of the side chains, we see the same picture for all 8 of them. There is little motion around the long axis. In most of these residues the long axis is not moving either, with the exceptions of residues 220 and 313. These residues make brief excursions to orientations closer to perpendicular to the z axis, but return to parallel. This means that all of the tyrosines in the aromatic bands are aligned with the pore wall, parallel to the z axis, and mostly experience only small fluctuations while maintaining their primary orientation. An example of this is given in Figure 3C, Tyr180. Of the remaining 5 tyrosines on the outside, Tyr98 is deeply buried in the outside wall, with its long axis exactly parallel to the z axis, and is completely immobile over 1 ns in all three monomers (Figure 3C). The most mobile tyrosine is Tyr90, at the intracellular side in the aromatic band. This residue shows considerable variation in both the orientation of the long axis and the plane of the ring (Figure 3C). Tyr139 and Tyr157 are opposing residues at the intracellular and extracellular sides but are located more toward the middle of the hydrophobic zone than the two regular bands (not shown). The plane of their rings is more or less perpendicular to the plane of the membrane and does not change much over the simulation. Their long axes are more or less parallel to the z axis. Finally, Tyr182 (not shown) is part of a short intracellular loop but could also be seen as part of the aromatic band. In all three monomers the plane of the ring is roughly perpendicular to the plane of the membrane, but the long axes assume a whole range of orientations. All tyrosines mentioned above form hydrogen bonds, mostly with lipids and solvents. The number of hydrogen bonds per OH group for tyrosines from the outer wall with solvent and lipids is about 2. Hydrogen bonding is not very specific; tyrosines form a donor with any oxygen atom in the lipids or water and an acceptor with the NH3 of the POPE headgroups. There is one Trp on the outside of the porin, positioned with the six-ring in line with two tyrosines that are part of the regular aromatic band on the extracellular side. During
17558 Biochemistry, Vol. 37, No. 50, 1998
Tieleman et al.
FIGURE 4: The bilayer hydrophobic thickness as a function of distance from a protein. The hydrophobic thickness was estimated as the distance between two atoms at the ends of a more or less continuous hydrophobic face of the peptides.
FIGURE 5: Lipid properties as function of distance from alamethicin (Alm 1), as example of a single helix. (A) Order parameters per tail carbon. (B) Average trans fraction per tail dihedral. (C) Average time between transitions per tail dihedral.
Table 3: Estimated Hydrophobic Thickness Lh of the Peptides and Lipidsa
26), and loses much of its helical structure at each end of the peptide. The effect of larger proteins on the local thickness is more pronounced. A considerable decrease in thickness is found near Flu T18. Alm 6H and Flu T22 do not induce much change, while Alm 6, which is somewhat broader and shorter than Alm 6H, induces a slight decrease in thickness. Note that for the larger Alm bundles reliable data is only available up to 2.4 nm away from the protein. The most interesting case is the porin, for which we can calculate the bilayer thickness over 4 nm away from the protein. The first two layers of lipids are affected and show a considerable decrease in thickness, whereas the next three layers are close to the POPE bulk thickness. 2.2.2 Lipid Order and Mobility. Single Helices. Overall, the effect of single helices on the order parameters of the lipids is small. One example of a single helix is given in Figure 5. The most obvious effect is a lowering of the order parameters in the upper half the lipid chains for the first two bins. The difference between the bins is smaller for the three Flu peptides, although the overall picture is the same. The effect of the helices on the trans ratios is small. The increased alternation of trans ratios with dihedral number indicates a tilt of the lipids adjacent to the protein, because such a tilt aligns one-half of the bonds more with the z axis than the other half. Rotations around bonds parallel to the z axis occur more readily than rotations around bonds that are at an angle with the z axis. This is consistent with a direct evaluation of the tilt angle. The dihedral transition times for Alm 1 in Figure 5 are representative for all 5 single helices, both in trend and in magnitude. Adjacent lipids exhibit slightly slower transitions. Helix Bundles. The four helix bundles give better statistics on lipids close to the protein due to their larger size. As representative examples the order parameters, trans fraction, and dihedral transition times are given for Flu T22 in Figure 6. There is a significant decrease of order parameters in the first ring around the bundle, but not further away. The same is true for the trans fractions and dihedral transition times, only the first ring of lipids is affected. The trans fractions
Alm
Flu 18 Flu 22 Flu 26 Flu 34 OmpF POPC POPE
Lh 2.7 nm 2.5 nm 3.0 nm 3.2 nm 3.2 nm 2.6 nm 3.0 nm 3.2 nm a Peptide hydrophobic thickness was calculated from a single structure, and lipid hydrophobic thickness was calculated as the average distance between the carbons next to the carbonyl group in both leaflets.
the simulation, this Trp in all three monomers orients mostly with the aromatic plane perpendicular to the plane of the membrane, with short fluctuations from an angle of 90° to about 60° with the plane of the membrane. The direction of the NH-vector fluctuates slowly and is different in the three monomers. In the first the angle between the NH vector and the z axis slowly fluctuates between about 50° and close to 0°. In the second monomer this angle remains more or less constant at close to 0°, and in the third monomer a very sharp transition in a few picoseconds is seen from an orientation parallel to the z axis to an orientation perpendicular to this axis. Each tryptophan is hydrogen-bonded to 6 different water molecules during the simulation, as well as to a neighboring glutamate. 2.2 Lipid Properties 2.2.1 Local Thickness of the Membrane. How do lipids behave in the vicinity of a protein? The first and most obvious change that might occur is a change in the local thickness of the bilayer around proteins. This local thickness is plotted in Figure 4. It is clear that it is difficult to get accurate averages, particularly in the case of single helices. The shortest peptide, Flu 18, causes a small local decrease in thickness. Flu 26 causes a considerable increase in thickness. Alm 1 and Alm 1H are similar and have little effect. The behavior of these four systems can be explained by the estimated hydrophobic thickness of the peptides, summarized in Table 3. Flu 26 has a considerably longer hydrophobic segment than the other peptides. Flu 34 is even longer, but this peptide bends in the middle, has a larger overall tilt (ca. 25° with the z-axis, compared to 10 for Flu
Lipid-Protein Interactions Studied by MD
Biochemistry, Vol. 37, No. 50, 1998 17559
FIGURE 6: Lipid properties as function of distance from Flu T22, as example of a helix bundle. (A) Order parameters per tail carbon. (B) Average trans fraction per tail dihedral. (C) Average time between transitions per tail dihedral.
FIGURE 7: Lipid properties as function of distance from OmpF. (A) Order parameters per tail carbon. (B) Average trans fraction per tail dihedral. (C) Average time between transitions per tail dihedral.
are lower than in bulk lipids and the dihedral transition times longer. This picture does not differ much from the results for Alm 6 and Alm 6H, which also show slightly lower order parameters and trans fractions in the first ring around the bundle and slightly slower dihedral transition times. In Flu T18 the order parameters are much lower than in the other three bundles, but in that system too the effect is limited to the first ring of lipids. Transition rates and trans fractions are little affected. The strong alternation with dihedral number of both the trans fraction and the transition time is related to a tilt of the lipid molecules. The average tilt angle of the palmitoyl chain in bulk POPC is 32° ( 1°. In Flu T18, which shows the biggest deviation from bulk near the protein, it is 38° ( 5°. Porin. The largest protein, OmpF, has the most pronounced effect on neighboring lipids. In Figure 4 we saw a decrease of 0.4 nm in the local hydrophobic thickness near the porin. Figure 7 shows the order parameters, trans fraction, and dihedral transition times for OmpF. The first ring of lipids around the protein has much lower order parameters than bulk lipids. This effect perpetuates to the second half of the lipid tails further away from the protein. The trans fraction is correspondingly lower. The overall increased order parameters and trans fraction compared to the same properties in the other systems are due to the use of POPE. Interestingly, the average dihedral transition times are not slower for lipids closer to the protein. The average tilt angle for the palmitoyl chain less than 0.8 nm away from the porin is 30° ( 2°, whereas the bulk value is 25° ( 1°. The effect of this extra tilt on the order parameters would be a decrease to 85%, close to the observed difference.
partitioning of host-guest peptides showed that the aromatic residues are very hydrophobic at the lipid-water interface (32). In other systems, more specific roles for, in particular, tryptophans have been postulated: anti-aggregation properties (18), conductance-regulating properties (33), and as determinants of translocation (34). From our simulations, it is not clear how aromatic residues would stabilize peptides or proteins in the bilayer. In a naive picture, the rings of aromatic residues would act as floats on the rough sea which is the headgroup zone. However, the mostly perpendicular orientation of tyrosines, the random orientations of phenylalanines, and their rapid changes in orientation do not support this image. Tyrosines and tryptophans are typically hydrogen-bonded, but not in a specific way. The number of sodium ions in the simulations is too low to draw conclusions about the effect of these ions on the aromatic side chains. The orientation of the side chains appears to be determined by their general amphiphilic character and their location in the membrane protein. The anchoring and mechanical stability provided by a large scale hydrophilic-hydrophobic-hydrophilic division of the outside of the protein should outweigh the effects of individual residues. However, although the simulations give insight into the dynamics and orientation of the aromatic side chains, they do not give information about the way the side chains would be distributed without being attached to a protein, which is related to the free-energy contribution of having an aromatic side chain at the interface. Detailed simulations on the side chains by themselves, NMR, and diffraction experiments will possibly provide this information.
3 DISCUSSION
The effect of the proteins on the lipids depends on the size of the protein and the difference in hydrophobic length between the protein exterior and the membrane interior. It has been shown that this hydrophobic mismatch can give rise to drastic effects such as phase separation or a change in phase (18, 20), when the mismatch is large enough. In our simulations, we have mild cases, which will be quite common.
3.1 Aromatic Residues The abundance of aromatic residues in membrane proteins has been noticed since the first structures became available (31). Several general roles have been attributed to them, including mechanical stability (17) and special dielectric properties at the acyl chain-water interface (16). Studies of
3.2 Lipid Properties
17560 Biochemistry, Vol. 37, No. 50, 1998 The largest effects occur in OmpF. The local thickness of the bilayer near the protein is about 0.4 nm less than the bulk value, the order parameters of carbons close to the protein are about 15% lower, and the trans dihedral fraction is about 10% lower. A larger tilt angle of the lipid chains with the bilayer normal can explain the decreased order parameters, if we assume that the deuterium order parameters can be taken as a product of two contributions: the orientation of the C-D bond around the lipids molecular axis, and the orientation of the molecular axis around the normal to the bilayer. An increased tilt of the lipid long axis would decrease the contribution of the molecular axis, and thereby SCD, even if the orientation of the deuterium atoms with respect to the chain did not change. In this case, the observed increase in tilt angle from 25° to 30° can account completely for the decreased order parameters. However, this increased tilt only reduces the bilayer thickness by about 5%. The remaining decrease in thickness must be caused by an increased gauche fraction. Because in the lipid chain adjacent bonds are at an angle of about 109° with each other, the increased tilt angle of the chain as a whole, combined with the fact that the lipid molecules do not rotate isotropically about their long axes (35), aligns one-half of the bonds more with the z axis than the other half. This causes a stronger alternation in trans fractions and dihedral transition times, because rotations around bonds aligned with the z axis occur more readily than rotations around bonds that are at an angle with the z axis. A similar effect is seen in the helix bundles. Both the order parameters and trans fractions are lower near the protein, and there is a significantly stronger alternation between odd and even dihedrals for trans fractions and transition times. The largest effects are found in Flu T18, where the order parameters for carbons near the headgroup are up to 30% lower for the neighboring lipids than for the other lipids. The increase in tilt angle of lipids near the bundle is correspondingly large. This correlates with Flu T18’s short hydrophobic length and a significant mismatch with the lipids. Near the single helices, all effects are less pronounced, but in general, there is an increase in the magnitude of the alternation between even and odd dihedrals. As an overall mechanism of adjustment, it seems lipids around a protein tilt a little, causing a decrease in order parameters and a stronger alternation between even and odd dihedral trans fractions and transition times. The effect on the magnitude of transition times is small. If more adjustment is needed, an increase in gauche fractions causes a further reduction in bilayer thickness. Although our simulations did not include a clear case where the protein is much longer than the hydrophobic thickness of the bilayer, we can speculate about what would happen. In such a case, the tilt of the lipid chains might decrease, causing an increase in order parameters, and possibly the trans fraction would increase. In support of this, the order parameters close to Flu 26 and Flu 34, which have a slightly longer hydrophobic length than the bilayer, are not significantly lower than those further away from the peptides, and they are significantly higher than the order parameters close to the other proteins. It is difficult to compare our findings to experimental information. Labeled phospholipids can give information about the properties of lipids near proteins, but the long time
Tieleman et al. scale in most experiments, compared to the diffusion coefficient of lipids, makes it hard to study lipids near proteins. Some data are available on rotation rates of lipids bound to or close to proteins, studied by, for example, fluorescence (36) or ESR (37), but the conclusions are too general to be compared with our simulations. Vogt et al. studied the order parameters of a carbon chain covalently attached to gramicidin A (38). This is an elegant way to make sure the lipid chain remains close to the protein. They found that the first carbons are greatly affected by the protein, with much lower order parameters and two deuterium splittings at each of the first three carbons, indicating that the deuteriums in those positions are not equivalent. The second half of the chain showed no change in order parameters compared to bulk values. T1 measurements showed that the mobility of the first few carbons is much lower than for the rest of the chain. These results agree with ours to the extent that the effects for a chain near a small peptide are not very drastic, but the systems are too different for an exact comparison. Older deuterium NMR on bilayers with embedded proteins found little effect of the proteins on order parameters, but averages over all lipids on a long time scale in an inhomogeneous membrane (39). However, recently De Planque et al. have presented a detailed study of the effect of peptides with different lengths on the mean thickness of bilayers (40). They calculated the mean bilayer thickness from deuterium NMR order parameters in bilayers of di-C12-PC, di-C14-PC, di-C16-PC, and di-C18-PC with a family of synthetic hydrophobic polypeptides (called WALP peptides) with different lengths. In general, longer peptides in short lipids caused an increase in the mean bilayer thickness, whereas short peptides in longer lipids caused a decrease. Certain intermediate combinations showed no measurable effect on thickness. In all cases where there was a hydrophobic mismatch, the change in thickness by itself was not enough to compensate completely for the mismatch. These results are very similar to the results for our single helix simulations, but the simulations give the local thickness as a function of the distance from the helices as well. There have been several theoretical studies of proteins in lipid bilayers. In an early molecular dynamics simulation on a model membrane with simplified lipids, Edholm and Johansson found no evidence for lipids that are tightly bound to helices (41). Several effects were observed, depending on the helix side chains, including changes in order parameters, trans ratios, and dihedral dynamics, but overall the effect of the single helices on the lipids was not strong. The number of lipids used in most recent molecular dynamics simulations of membrane proteins is usually too small to calculate lipid properties as function of distance, and the focus is usually on the protein. Shen et al. found little effect of a polyalanine helix on nearby lipids (11). Gramicidin A caused a significant increase in lipid order parameters (DMPC) at a high protein concentration (7), in agreement with experimental data for this system. Sperotto and Mouritsen have studied the coherence length of the perturbation of lipids due to proteins as a function of hydrophobic mismatch and temperature (42). They use a Monte Carlo approach that does include hydrophobic mismatch and different protein sizes. Near the main phase transition temperature, the coherence length is several nanometers, but
Lipid-Protein Interactions Studied by MD well above this transition, the coherence length decreases to a few tenths of a nanometer. For larger proteins and larger hydrophobic mismatch, the coherence length increases. The results from our simulations agree qualitatively with this dependency on protein size and hydrophobic mismatch, although our resolution of 0.8 nm does not allow an accurate estimate of a coherence length from the simulations. We also did not study the effect of temperature. Future Work. Although the set of simulations described in this paper contains important differences in hydrophobic length and size, it would be interesting to use a more systematic approach to study the effect of size and hydrophobic mismatch in detail. One intriguing possibility is to combine much simplified proteins such as those used in the Monte Carlo approach of Sperotto (42) with simplified lipids in MD simulations. Potentials of mean force for such simulations can be derived from distribution functions obtained from detailed simulations. This would allow the extension of simulation methods to the study of lipid mixtures and possibly of phase segregation. REFERENCES 1. Cross, T. A., and Opella, S. J. (1994) Curr. Opin. Struct. Biol. 4, 574-581. 2. Opella, S. J. (1997) Nat. Struct. Biol. NMR supplement, 845848. 3. Thurmond, R. L., and Lindblom, G. (1997) Curr. Top. Membr. 44, 104-166. 4. Merz, K. M., Jr. (1997) Curr. Opin. Struct. Biol. 7, 511-517. 5. Tieleman, D. P., Marrink, S. J., and Berendsen, H. J. C. (1997) Biochim. Biophys. Acta 1331, 235-270. 6. Edholm, O., Berger, O., and Ja¨hnig, F. (1995) J. Mol. Biol. 250, 94-111. 7. Woolf, T. B., and Roux, B. (1996) Proteins: Struct., Funct. Genet. 24, 92-114. 8. Tieleman, D. P., and Berendsen, H. J. C. (1998) Biophys. J. 74, 2786-2801. 9. Woolf, T. B. (1997) Biophys. J. 73, 2476-2393. 10. Woolf, T. B. (1998) Biophys. J. 74, 115-131. 11. Shen, L., Bassolino, D., and Stouch, T. (1997) Biophys. J. 73, 3-21. 12. Belohorcova´, K., Davis, J. H., Woolf, T. B., and Roux, B. (1997) Biophys. J. 73, 3039-3055. 13. Killian, J. A. (1992) Biochim. Biophys. Acta 1113, 391-426. 14. Koeppe, R. E., Killian, J. A., and Greathouse, D. V. (1994) Biophys. J. 66, 14-24. 15. Ketchem, R. R., Hu, W., and Cross, T. A. (1993) Science 261, 1457-1460. 16. Cowan, S. W., Schirmer, T., Rummel, G., Steiert, M., Ghosh, R., Pauptit, R. A., Jansonius, J. N., and Rosenbusch, J. P. (1992) Nature 358, 727-733. 17. Kreusch, A., and Schulz, G. E. (1994) J. Mol. Biol. 243, 891905. 18. Killian, J. A., Salemink, I., de Planque, M. R. R., Lindblom, G., Koeppe, R. E., II, and Greathouse, D. V. (1996) Biochemistry 35, 1037-1045.
Biochemistry, Vol. 37, No. 50, 1998 17561 19. Ben-Shaul, A., Ben-Tal, N., and Honig, B. (1996) Biophys. J. 71, 130-137. 20. Mouritsen, O. G., Dammann, B., Fogedby, H. C., Ipsen, J. H., Jeppesen, C., Jorgensen, K., Risbo, J., Sabra, M. C., Sperotto, M. M., and Zuckermann, M. J. (1995) Biophys. Chem. 55, 55-68. 21. Tieleman, D. P., Berendsen, H. J. C., and Sansom, M. S. P. Biophys. J. (in press). 22. Tieleman, D. P., Sansom, M. S. P., and Berendsen, H. J. C. Biophys. J. (in press). 23. Forrest, L. R., Tieleman, D. P., and Sansom, M. S. P. Biophys. J. (submitted for publication). 24. Forrest, L. R., Tieleman, D. P., and Sansom, M. S. P. FEBS Lett. (submitted for publication). 25. Egberts, E., Marrink, S. J., and Berendsen, H. J. C. (1994) Eur. Biophys. J. 22, 423-426. 26. Berger, O., Edholm, O., and Ja¨hnig, F. (1997) Biophys. J. 72, 2002-2013. 27. van Gunsteren, W. F., Kru¨ger, P., Billeter, S. R., Mark, A. E., Eising, A. A., Scott, W. R. P., Hu¨neberger, P. H., and Tironi, I. G. (1996) Biomolecular Simulation: The GROMOS96 Manual and User Guide, Biomos/Hochschulverlag AG an der ETH Zu¨rich, Groningen/Zu¨rich. 28. Berendsen, H. J. C., Postma, J. P. M., van Gunsteren, W. F., DiNola, A., and Haak, J. R. (1984) J. Chem. Phys. 81, 36843689. 29. Berendsen, H. J. C. (1996) Science 271, 954-955. 30. van der Spoel, D., van Buuren, A. R., Apol, E., Meulenhoff, P. J., Tieleman, D. P., Sijbers, A. L. T. M., van Drunen, R., and Berendsen, H. J. C. (1996) Gromacs User Manual Version 1.2, http://rugmd0.chem.rug.nl/gmx. 31. von Heijne, G. (1997) Membrane Protein Assembly, Chapman & Hall, Austin, TX. 32. Wimley, W. C., and White, S. H. (1996) Nat. Struct. Biol. 3, 842-848. 33. Hu, W., and Cross, T. A. (1995) Biochemistry 34, 1414714155. 34. Schiffer, M., Chang, C.-H., and Stevens, F. J. (1992) Protein Eng. 5, 213-214. 35. van der Ploeg, P., and Berendsen, H. J. C. (1982) Mol. Phys. 49, 243-248. 36. Pap, E. H. W., Hanicak, A., van Hoek, A., Wirtz, K. W. A., and Visser, A. J. W. G. (1995) Biochemistry 34, 9118-9125. 37. Swamy, M. J., and Marsh, D. (1997) Biochemistry 36, 74037407. 38. Vogt, T. C. B., Killian, J. A., and De Kruijff, B. (1994) Biochemistry 33, 2063-2070. 39. Seelig, J., and Seelig, A. (1980) Q. ReV. Biophys. 13, 19-61. 40. de Planque, M. R. R., Greathouse, D. V., Koeppe, R. E., II, Scha¨fer, H., Marsh, D., and Killian, J. A. (1998) Biochemistry 37, 9333-9345. 41. Edholm, O., and Johansson, J. (1987) Eur. Biophys. J. 14, 203-209. 42. Sperotto, M. M., and Mouritsen, O. G. (1991) Biophys. J. 59, 261-270. BI981802Y
