proteins STRUCTURE O FUNCTION O BIOINFORMATICS

Coarse grain lipid–protein molecular interactions and diffusion with MsbA flippase Andrew B. Ward,1 Olgun Guvench,2 and Ronald D. Hills Jr.2* 1 Department of Molecular Biology, The Scripps Research Institute, La Jolla, California 2 Department of Pharmaceutical Sciences, University of New England, Portland, Maine

ABSTRACT Coarse-grained (CG) modeling has proven effective for simulating lipid bilayer dynamics on scales of biological interest. Modeling the dynamics of flexible membrane proteins within the bilayer, on the other hand, poses a considerable challenge due to the complexity of the folding or conformational landscape. In the present work, the multiscale coarse-graining method is applied to atomistic peptide-lipid ‘‘soup’’ simulations to develop a general set of CG protein–lipid interaction potentials. The reduced model was constructed to be compatible with recent solvent-free CG models developed for protein– protein folding and lipid–lipid model bilayer interactions. The utility of the force field was demonstrated by molecular dynamics simulation of the MsbA ABC transporter in a mixed DOPC/DOPE bilayer. An elastic network was parameterized to restrain the MsbA dimer in its open, closed and hydrolysis intermediate conformations and its impact on domain flexibility was examined. Conformational stability enabled long-time dynamics simulation of MsbA freely diffusing in a 25 nm membrane patch. Three-dimensional density analysis revealed that a shell of weakly bound ‘‘annular lipids’’ solvate the membrane accessible surface of MsbA and its internal substrate-binding chamber. The annular lipid binding modes, along with local perturbations in head group structure, are a function of the orientation of grooves formed between transmembrane helices and may influence the alternating access mechanism of substrate entry and translocation. Proteins 2012; 80:2178–2190. C 2012 Wiley Periodicals, Inc. V

Key words: bilayer thickness; membrane protein conformation; protein–phospholipid interactions; ABC transporter; multiscale coarse-graining; molecular dynamics.

INTRODUCTION Molecular dynamics simulation can yield atomic level insight into a variety of biological processes complementary to experimental data. Because of the breadth of timescales and ruggedness of the macromolecular energy landscape, traditional atomistically detailed models have limited use in problems of larger scale. Coarse-grained (CG) resolution models seek to reduce the available conformational space by eliminating unwanted molecular degrees of freedom. The removal of atomistic noise has particularly proven useful in the areas of protein folding and cell membranes. Single bead polypeptide models incorporating one interaction center at each Ca have successfully described folding mechanisms1,2 and other conformational transitions3 by encoding the network of residue contacts present in the native structure. CG representations of the head and tail groups of amphipathic phosphodiglycerides have been used extensively to model the heterogeneous distribution and diffusion of lipids in multicomponent bilayers.4–6 An outstanding question that modeling is poised to address concerns the dynami-

2178

PROTEINS

cal interplay of lipids with membrane proteins in the context of the cellular milieu.7–9 The well-characterized MARTINI force field,10 with one interaction center for every four heavy atoms on average, has enjoyed widespread use.9,11–14 Standard Lennard–Jones and Coulombic potentials represent nonbonded interactions between charged amine and phosphate head groups, two polar glycerol ester backbone sites and apolar groups for each tail. Parameters were chosen to fit the nonpolar-water phase partitioning of small chemical fragments. The MARTINI cocktail also

Abbreviations: ABC, ATP-binding cassette; CG, coarse-grained; DOPC, 1,2-dioleoyl-sn-glycero-3-phosphocholine; DOPE, 1,2-dioleoyl-sn-glycero-3-phosphoethanolamine; EL, extracellular loop; MS-CG, multiscale coarse graining; MSD, mean square displacement; NBD, nucleotide-binding domain; RMSD, root mean square deviation; TMD, transmembrane domain. Grant sponsors: University of New England; Scripps Research Institute *Correspondence to: Ronald D. Hills Jr, 716 Stevens Avenue, Portland, ME 04103. E-mail: [email protected]. Received 4 February 2012; Revised 10 April 2012; Accepted 25 April 2012 Published online 2 May 2012 in Wiley Online Library (wileyonlinelibrary.com). DOI: 10.1002/prot.24108

C 2012 WILEY PERIODICALS, INC. V

Coarse Grain Lipid–Protein Molecular Interactions

incorporates explicit water with every four water molecules mapping to a single solvent CG particle. Slightly larger antifreeze particles are included to prevent CG water freezing in standard room temperature simulations.15 After illustrating the use of the model in studying bilayer properties, Marrink et al. extended the force field to peptides.16 In contrast to backbone centric CG approaches used for protein folding,17,18 MARTINI consists of one to four amino acid sidechain sites attached to Ca beads. Sidechain centric approaches allow for a description of protein packing but require secondary structure restraints to mimic backbone hydrogen bonding and maintain structural stability.19,20 MARTINI sidechains were parameterized to reproduce amino acid water-bilayer phase partitioning.21,22 An alternative to empirical model fitting is the parameterization of CG models against atomistic simulation, termed multiscale modeling.23 The present work uses multiscale parameterization to develop a general waterfree24,25 CG model for the simulation of large membrane–protein complexes. Model parameters are obtained from detailed atomistic simulations of reference systems using the multiscale coarse-graining method (MSCG).26,27 Making no assumptions other than pairwise decomposition of interactions, MS-CG optimizes the many-body CG potential of mean force that best reproduces the atomistic ensemble of configurations through a variational procedure known as force matching. Given the complexity of conformational space, including sufficient substrates in the atomistic reference simulation is key to constructing a general model with predictive power.28 MS-CG was recently applied to unfolding and aggregation simulations of proteins and peptides to construct a generic protein model incorporating anisotropic sidechain packing and polarity.20 To extend the sidechain model of Hills et al. to membrane–protein systems, we develop a general set of lipid–protein interaction parameters. Microsecond atomistic reference simulations were performed on a ‘‘soup’’ of peptide–lipid aggregates in water and the interactions were decomposed with MSCG. The CG sites were chosen to be compatible with lipid–lipid interactions previously obtained for a model solvent-free bilayer by Lu and Voth.29 To maintain protein stability, the Ca 1 sidechain model was combined with an elastic network,30 and the impact on conformational flexibility was assessed. The CG model is illustrated by application to the ATPbinding cassette (ABC) transporter MsbA. MsbA, a bacterial lipopolysaccharide and multidrug exporter, is a 1,164-residue dimer consisting of two helical transmembrane domains (TMD) and two cytosolic nucleotidebinding domains (NBD) (Fig. 1). Crystallographic31 and spin label32 studies have characterized three conformational intermediates during the MsbA transport cycle. In the inward facing open conformation, substrate is free to diffuse along the inner leaflet and enter a large internal

Figure 1 Topology of the MsbA dimer. One monomer in the inward facing open crystal conformation is colored, including its cytosolic NBD extension. Transmembrane helices (TM1-TM6) and extracellular loops (EL1-EL3) comprise each TMD. TM4/TM5 move as a rigid unit upon closing of the flexible EL2/EL3 hinge to form the closed conformation. During ATP hydrolysis, TM1/TM2 separate from TM3/TM6 to form the outward facing conformation. [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

binding chamber between opposing TMDs. Bending of a flexible hinge between the TMDs brings the NBDs into partial contact, forming the closed conformation. Binding of the canonical ATP dimer sandwich then forms the complete NBD–NBD interface, whose structure has been trapped with the bound hydrolysis transition state analog AMP-PNP. This structure, which we refer to as the hydrolysis intermediate, has an altered (outward facing) orientation of TM helices responsible for releasing lipid substrate to the extracellular leaflet. Long time molecular dynamics simulations were performed for the three distinct ATP-dependent MsbA conformations freely diffusing in a 25 nm bilayer patch. While ATP-hydrolysis induced conformational changes have been modeled elsewhere for MsbA33,34 and transporter homolog Sav1866,35–39 the CG model enabled examination of the role of molecular level lipid–protein interactions and diffusion at the membrane accessible protein surface. Complementary to growing biophysical evidence,7,40,41 a shell of ‘‘annular lipids’’ was observed to become ordered around the protein in a conformation-dependent manner. Rather than a gradual change in the hydrophobic thickness, local pockets of lipids exhibited perturbations in bilayer structure near the membrane accessible protein surface. PROTEINS

2179

A.B. Ward et al.

MATERIALS AND METHODS Model parameterization

The membrane-protein CG force field was constructed as follows. Bonded and nonbonded protein–protein interactions were treated using the sidechain centric model of Hills et al. developed for the folding of aqueous proteins, in which five CG site types represent backbone alpha carbons and apolar, polar, positive and negative sidechain functional groups consisting of two or three heavy atoms.20 Backbone and sidechain groups are fully flexible in the model, being parameterized only to enforce chirality, Ramachandran preferences and predominant bond angle distributions from unfolded peptides. Bonded and nonbonded lipid–lipid interactions employed the MS-CG model developed by Lu and Voth29 for a 1:1 mixed DOPC/DOPE bilayer, consisting of nine unique CG site types including five acyl chain types to incorporate depth-dependent hydrophobic interactions. Flexible bonded terms were developed from Boltzmann inversion of atomistic bilayer simulations. In both MS-CG models, solvent degrees of freedom were integrated out resulting in a water-free force field allowing fast diffusion at the expense of rigorous timescales. The zwitterionic DOPC/DOPE model bilayer was chosen for its stability in the absence of solvent with a hydrophobic thickness that matches MsbA. Protein–lipid interactions were parameterized in the present work. To aid convergence of the pair potentials, only two site types were considered for the acyl chains (Fig. 2). The location and number of lipid sites were unaltered since they comprise atom groupings comparable to the Hills et al. protein model. Atomistic reference simulations were performed in Gromacs42 using all-atom OPLS43 protein, SPC solvent, and united atom lipid parameters.44,45 The default Gromacs parameters were used with particle mesh Ewald, 1.2-nm cutoff for grid-based neighbor searching every 10 steps, 1 nm short-range cutoff, Nose´-Hoover thermostat, and 2 fs timestep as in previous work.20 Two molecules of the 12-residue Trpzip b-hairpin (1LE1.pdb) were placed in a (4 nm)3 periodic water box with three DOPC molecules, bringing the solute/water concentration to 18%. Fifty NVT simulations, each 50 ns in length, were performed starting from random configurations for a total of 2.5 ls. Atom coordinates and forces were recorded every 5 ps. Inspired by the weak temperature dependence of interaction potentials observed previously,20 the concentrated peptide–lipid soup was simulated at 498 K to rapidly explore conformational space. CG protein–lipid interactions were parameterized from atomistic trajectories using the MSCGFM program.46 Pairwise force curves were optimized using 0.005 nm linear splines out to a distance of 2 nm. Block averaging was performed on the least squares solutions of sets of 6000 randomly shuffled frames. The process described

2180

PROTEINS

Figure 2 Mapping of CG site types for parameterizing protein–lipid interactions. The lipid CG sites of Lu and Voth29 are used with the modification that protein–lipid interactions are defined for six rather than nine lipid site types. DOPC and DOPE are represented by choline (blue), ethanolamine (pink), phosphatidyl (red), glyceryl (green), ester (light blue), and tail (gray) functional groups consisting of 5, 2, 6, 3, 5, and 3 heavy atoms, respectively. CG sites are connected by virtual bonds (red lines). [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

above for DOPC was subsequently repeated with a 2Trpzip 1 3DOPE system. Force curves were integrated and smoothed into tabulated potentials for input in Gromacs 4. Excepting the two cationic head groups, the best-behaved (smooth, decays to zero) potential of the DOPC and DOPE systems was chosen for each protein– lipid site pairing or, if nearly identical, the two were averaged. Long-range potentials were truncated or switched off linearly from 1 to 1.2 nm resulting in integration errors below 1 kJ mol21.

CG dynamics simulations

CG simulations were performed in Gromacs 4 using Langevin dynamics with a 2 ps inverse friction constant and 2 fs timestep to maintain high frequency vibrations, consistent with previous work.20 All nonbonded interaction pairs separated by three or more bonds and <1.2 nm were updated every step with a grid neighbor search.

Coarse Grain Lipid–Protein Molecular Interactions

Replica exchange molecular dynamics47 was performed on CG representations of the 2Trpzip 1 3DOPC and 2Trpzip 1 3DOPE systems for 400 ns. Forty-eight temperature windows were exponentially spaced from 300 to 700 K, excepting seven temperatures centered at 498 K to aid comparison with atomistic simulations. Exchanges were attempted every 200 ps with >50% acceptance ratios and snapshots were recorded at 2 ps intervals. Protein–lipid interactions were monitored from the pairwise radial distribution functions of intermolecular site–site pairings. Stable CG DOPC/DOPE bilayer simulations by Lu and Voth with 512 lipids were previously performed at constant volume with an area per head group of 0.58 nm2, closely matching the average value of 0.568 nm2 from atomistic NPT simulations.29,48 While a sophisticated implementation of MS-CG has been developed for NPT simulations of liquids,49 the present water-free model was developed for simulations in the NVT ensemble. An equilibrated 512-lipid system of 0.58 nm2 area per lipid was duplicated in two dimensions to obtain a 257 A˚ by 233 A˚ rectangular box. The resulting 2048-lipid system required density equilibration to minimize the occurrence of lipid tails popping out of the bilayer in heterogeneous patches. The fact that lipid tails do not stick out of the bilayer in 512-lipid simulations of 0.58 nm2 area per lipid suggests that the boundaries have a non-negligible influence in the smaller periodic box. Patches of disrupted lipids were manually deleted over the course of short equilibration simulations at 310 K totaling 250 ns, finally resulting in a stable 1888-lipid system. The 1888-lipid system’s area per lipid of 0.63 nm2 coincides with the upper limit of available experimental data for 1:1 DOPC/ DOPE.48 The closed structure of MsbA (3B5X.pdb) was then positioned in the equilibrated bilayer with the aid of the OPM database.50 Lipids with protein van der Waals overlaps were manually deleted during 300 ns of bilayer equilibration with the protein held fixed, to obtain a 1842 lipid 1 MsbA system (967:875 DOPC/ DOPE ratio). RMSD alignment with this system was used to generate starting configurations for simulations of the inward facing open conformer (3B5W.pdb) and hydrolysis intermediate (3B60.pdb) structures. In the construction of subsequent systems for simulation, lipids with close (5 A˚) Ca protein contacts were moved to random positions in the bilayer bulk rather than deleted to maintain a constant density. To aid comparison between the three ortholog structures, all simulations were performed with the MsbA Vibrio cholerae sequence, in which undefined sidechain atoms were placed via homology using the automated SWISS-MODEL.51 Dynamics trajectories were visualized with VMD.52 The lateral mean square displacement (MSD) was computed in the xy plane using the g_msd function in Gromacs, which corrects for periodic boundary conditions. Square displacements in choline, ethanol-

amine, or protein center of mass positions were ensemble averaged restarting every 200 ps in the trajectory. Lateral diffusion coefficients (D) were obtained from the time dependence fit over the 0–25 ns regime: MSD 5 4Dt. To examine annular lipids, the time-averaged number densities of tail and head groups were computed in gridcount53 after superimposing the TM helices of each snapshot with the final configuration. Only aligned frames from the second half of each full-length trajectory were used.

RESULTS AND DISCUSSION Model development and assessment

In addition to integrating out the water degrees of freedom, the CG scheme employed invokes three primary approximations. High temperature atomistic reference simulations were performed on a peptide–lipid water soup to attain faster binding and exchange than could be achieved with a larger peptide-bilayer simulation. The small reference system was needed to obtain converged MS-CG potentials for the set of protein–lipid CG site types. While depth-dependent lipid–lipid interactions are crucial for capturing bilayer structure, there omission in the peptide–lipid system will need to be examined in future validations of the model.21 Second, protein– protein interactions rely on a model for soluble protein folding derived from reference simulations of peptide aggregation in water. Figure 3 compares representative

Figure 3 Correspondence between interaction potentials generated from atomistic force matching of protein in the presence (solid) and absence (dash) of lipid molecules. Attractive pairwise potentials generated from a peptide– lipid soup simulation of 2Trpzip 1 3DOPC compare well with unfolding simulations of aqueous Trpzip.20 Nonbonded potentials are aligned by their attractive wells for the following protein–protein CG particle interactions: Ca–Ca (green), apolar–apolar (black), polar–polar (cyan), and positive–negative (pink). [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

PROTEINS

2181

A.B. Ward et al.

Figure 4 Tabulated nonbonded potentials employed between CG site types for protein–lipid interactions. (E) Potentials are shown for protein alpha carbons paired with the six lipid site types: choline (blue), ethanolamine (pink), phosphatidyl (red), glyceryl (green), ester (cyan), and tail (black) groups. Similarly, potentials are shown for the lipid sites paired with positive (A), negative (B), polar (C), and apolar (D) protein sidechain sites. [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

protein–protein potentials generated from the 2Trpzip 1 3DOPC system with those from a system of three Trpzip peptides in the same size water box (10% protein solution). Despite their good agreement, extrapolation of these interactions to a protein within the bilayer is nevertheless an essential assumption of the model. A third approximation involves the selection of derived potentials from the Trpzip/DOPC and Trpzip/DOPE mixtures. Force-matched potentials for apolar-head group interactions exhibited attractive strengths of order kcal mol21 (1 kcal 5 4.184 kJ). Preliminary CG simulations of the MsbA 1 bilayer system with such potentials demonstrated instability of the bilayer near the protein, due to prominent interactions between lipid head groups and apolar sidechains of the transmembrane helices. Mindful of the similarity in size (volume exclusion) and polarity (attractive well depth) of the five chemical site types assigned for protein and lipid, the solution employed henceforth was to replace selected potentials with those previously derived for lipid bilayer head–tail interactions.29 The apolar–choline, apolar–ethanolamine, and apolar–phosphate potentials were replaced by repulsive interactions developed for the ‘‘SM’’ tail groups: SM– choline, SM–ethanolamine, and SM–phosphate, respectively. An atypically strong interaction was also observed between ethanolamine and negative sidechains (5 kcal well depth), and it was replaced with the phosphate–eth-

2182

PROTEINS

anolamine pair potential. The resulting set of selected protein–lipid potentials is shown in Figure 4. CG simulations of the 2Trpzip 1 3DOPC and 2Trpzip 1 3DOPE systems were used to compare their aggregation profiles with atomistic reference simulations. Figure 5 shows the radial distribution functions for unique protein–lipid site pairings. Because the CG lipid–lipid potentials were developed from lipid bilayer simulations, high accuracy in the lipid–water soup simulations was not expected for all pairwise correlations. Protein elastic network

To prevent misfolding at the physiological temperature used in bilayer simulations (310 K), proteins treated with the model of Hills et al. require structural restraints to stabilize the Ca backbone conformation of the native state.20 A Ca elastic network was used to restrain MsbA in each of three conformations. Stabilizing harmonic interactions can be placed between protein backbone residues within a given distance of each other in the reference structure. While elastic networks can reproduce conformational transitions in normal mode analysis,3 they exhibit less flexibility in molecular dynamics where the bonds are unbreakable. Four different elastic network models were examined for their effect on backbone RMSD fluctuation in molecular dynamics simulations

Coarse Grain Lipid–Protein Molecular Interactions

Figure 5 Correspondence between atomistic (dash) and CG (solid) simulations of a peptide–lipid soup. The pair distribution functions between nonbonded CG sites are shown for protein sites paired with choline (blue), ethanolamine (pink), phosphatidyl (red), glyceryl (green), ester (cyan), and tail (black) groups. Protein–lipid correlations are shown for backbone alpha carbons (A) and apolar (B), polar (C), positive (D), and negative (E) sidechain sites. Lipid–lipid interactions were developed from a lipid bilayer simulation rather than this system. [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

(Table I). Alpha carbon pairs separated by two or more bonds and within a given cutoff distance in the crystal structure were assigned a harmonic potential with a uniform force stiffness constant, KSPRING: U 5 KSPRING(rij r0ij)2/2. In vacuo simulations of closed MsbA demonstrate that a cutoff of 10.5 A˚ and KSPRING 5 150 kJ mol–1 nm– 2 results in the appropriate balance of flexibility and stability in the TM helices and mixed a/b NBDs. Similar parameters have been successfully employed in simulations with the MARTINI model.30 Subsequent MsbA 1 bilayer simulations of the open, closed or hydrolysis intermediate structures employed backbone restraints of this form for all native contacts in the given starting crystal structure. The compact closed structure did not undergo conformational transitions because the harmonic spring energy increases with increasing contact distance.54 Application of an elastic network to the open MsbA structure, on the other hand, did allow for conformational change since no restraints were assigned between the two NBDs and there was sufficient flexibility in the hinge formed by extracellular Loops 2 and 3 (EL2/EL3). Similar to solvent-free simulations of adenylate kinase with the Hills et al. model,20 an in vacuo simulation starting with inward facing open MsbA resulted in rapid association of

the two NBDs via closing of the EL2/EL3 hinge. Within 15 ns the two monomers adopt a compact conformation with the NBD–NBD interface fully formed. The relative arrangement of the NBDs resembled the close packing of the hydrolysis intermediate structure containing the canonical ATP sandwich,31 with the Ca RMSD alignment dropping swiftly from 28 A˚ down to 10 A˚. Figure 6 shows the RMSD of the transmembrane region during the simulation. All 12 TM helices became tightly packed, resembling the TMDs of the closed MsbA structure. As with the open to closed transition of adenylate kinase,

Table I

Ca RMSD (A˚) of Elastic Network Modelsa Cutoff (
)

KSPRING (kJ mol–1 nm–2)

8 9 10.5d 11

500 1000 150d 500

TMDsb RMSD (
) 5.7 1.2 2.0 0.8

   

0.3 0.1 0.2 0.1

NBDsb 11.6 2.0 4.0 1.2

   

1.2 0.4 0.9 0.3

NBD Ac 3.5 1.0 1.6 0.6

   

0.2 0.2 0.2 0.04

All Ca 7.1 1.4 2.5 0.8

   

0.4 0.2 0.3 0.1

a

100 ns in vacuo CG simulations of closed MsbA. Mean RMSD  standard deviation after alignment with TM helices in crystal structure. c RMSD of NBD of monomer A in its own reference frame. d Parameters used for simulations of MsbA 1 bilayer system. b

PROTEINS

2183

A.B. Ward et al.

Figure 6 Inward facing open to closed conformational transition in MsbA. CG simulations were performed with Ca restraints favoring the open conformer in vacuo (dark gray) as well as in equilibrated (black) and unequilibrated (light gray) MsbA 1 bilayer systems (Fig. 7). The Ca RMSD of the transmembrane region was computed relative to the open (A) and closed (B) crystal structures.

rapid domain collapse can be attributed to the dominant influence of hydrophobic surface tension in solvent-free simulations.55 The final RMSD is indicative of the NBDs binding in a nonnative orientation, underscoring the limited use of the elastic network in molecular dynamics. Future improvements in the model will require encoding multiple conformations with breakable Lennard–Jones54 or Gaussian well56 type bonds. Simulations starting from the open conformation were next performed in the presence of the bilayer to examine its effect on conformational stability. For the MsbA 1 bilayer system equilibration was performed with the protein fixed in the open conformer, allowing lateral diffusion of lipids into the available space between the two TMDs of the MsbA dimer (Fig. 7). The diffusion of phospholipids into the large substrate binding chamber is not surprising and may influence the mechanism of substrate entry and translocation. Diffusion of lipopolysaccharide substrate along the inner leaflet of the membrane and into the internal chamber of MsbA has previously been implicated in the alternating access transport mechanism.31,32 After 200 ns equilibration, backbone constraints were turned off to allow for examination of transporter dynamics. In contrast to the in vacuo MsbA simulation, the open conformer remained stable using the same CG model 1 flexible Ca elastic network within the bilayer. The effect of lipid bilayer equilibration on conformational dynamics is illustrated in Figure 6, which also shows the open to closed conformational transition

2184

PROTEINS

observed in a MsbA 1 bilayer simulation in which phospholipids were not given time to equilibrate in the internal chamber. The stability of the open inward facing packing arrangement of TM helices can be seen in Figure 8. Hinge bending motions were not observed in the extracellular loops and the transmembrane region exhibited minimal RMSD fluctuation. Each NBD experienced minimal intradomain vibration, though rigid body motions placed the relative orientation of the NBDs 11  2 A˚ RMSD from that of the open crystal structure. The stabilizing effect of the bilayer on protein conformation was also evident in CG simulations of the outward facing hydrolysis intermediate conformation, in which helices TM1/TM2 are separated from TM3/TM6 resulting in accessibility to the outer leaflet of the membrane. In vacuo MsbA simulation with an elastic network stabilizing the hydrolysis intermediate conformer exhibited transient partial closing of the accessible chamber due to association of TM1 helices from each monomer near Ser46, three turns from extracellular Loop 1. During the 100 ns simulation, seven transitions to and from this lower energy structure were observed for a population fraction of 70%. The two conformational basins exhibited TMD RMSDs of 2.4 and 3.9 A˚ from the crystal structure. By comparison, a 200 ns

Figure 7 Equilibration of open MsbA 1 bilayer system. The protein backbone was held fixed in the open conformer to allow lipid diffusion into the transmembrane pore during bilayer equilibration. Snapshots are shown before (A) and after (B) 200 ns of equilibration with protein in blue and lipids colored as in Figure 2 with the tails depicted yellow. [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

Coarse Grain Lipid–Protein Molecular Interactions

Figure 8 Structural fluctuations observed in CG simulations of the equilibrated open transporter conformer. The transmembrane region exhibited a Ca RMSD of 3.4  0.5 A˚ from the crystal structure (pink surface). The NBDs experienced intradomain fluctuations of 1.4  0.1 A˚ RMSD and interdomain fluctuations of 5.9  1.8 A˚ (calculated in the reference frame of the transmembrane region). Superimposed backbone traces of TM1/ TM2, NBD/TM3/TM6 and TM4/TM5 of monomer A are colored olive, cyan and green, respectively, with b-strands in yellow. A cutaway view of the final lipid coordinates is shown. [Color figure can be viewed in the online issue, which is available at wileyonlinelibrary.com.]

MsbA 1 bilayer simulation of the conformer did not display this transition, exhibiting an all Ca mean RMSD of 1.8 A˚. Improved conformational stability in the presence of lipid bilayer justifies its explicit incorporation compared to implicit computational models of the membrane.57 Nevertheless, the propensity for closing of the accessible chamber in the hydrolysis intermediate state may be relevant to completion of the conformational cycle upon nucleotide release.34,36,37,58 Lipid–protein interactions

Final MsbA 1 bilayer simulations did not employ restraints to enforce the orientation or position of protein within the membrane. MsbA freely diffused in the bilayer for a total of 0.5, 1, and 0.2 ls with Ca elastic networks stabilizing the open, closed, and hydrolysis in-

termediate conformers, respectively. Figure 9(A–D) shows that the most stable orientation of each conformation is approximately normal to the bilayer plane. The balance of CG interactions effectively provided a harmonic restoring force for tilting of the protein with respect to the bilayer normal [Fig. 9(D)]. No significant change in the depth of proteins within the membrane was observed despite considerable lateral diffusion of their center of mass. Though parameterization of the present model made no assumptions regarding protein architecture, the correct insertion depth and orientation of other membrane proteins may not similarly be a local energy minimum. That the insertion of a large multidomain transporter was successful lends credibility to the overall sidechain and lipid polarities. Capturing the precise insertion energetics of amino acids in a water-free bilayer will be necessary for future improvements to the model. PROTEINS

2185

A.B. Ward et al.

Figure 9 Freely diffusing protein. The tilt angle of MsbA relative to the bilayer normal was determined from the vector spanning the NBD and TMD centers of mass for simulations of the open (A), closed (B), and hydrolysis intermediate (C) conformations. The corresponding Boltzmann-inverted probability distributions were computed using U 5 –kT ln[p(y)/sin y] (D). Linearity in the MSD for lateral protein diffusion enabled the estimation of lower and upper bounds for the diffusion coefficient (E). Two consecutive simulations were performed for the closed conformer (dark gray).

To examine the interaction between lipid and protein, and make comparisons to previous models, their diffusion was measured. The CG bilayer model of Lu and Voth exhibited a lateral lipid diffusion coefficient of D 5 400 lm2 s21, measured from the linear slope in MSD versus time for the phosphate positions.29 This is in contrast to nanosecond atomistic bilayer simulations, which are slow to reach the linear MSD regime and where all lipids do not escape their initial surroundings.48,59 In close agreement with experiment, atomistic simulation of the liquid crystal phase yields diffusion constants of  8 lm2 s21, requiring on average 1 ls for a lipid to travel 5 nm.60,61 Although absolute timescales cannot readily be compared, the order of magnitude effective speedup in the CG dynamics can be attributed to the reduction in degrees of freedom and smoother landscape.6,62 Present CG simulations employed a Langevin frictional force. Bilayer lipid diffusion coefficients for DOPC and DOPE in the 1888-lipid system lacking protein were 262  6 lm2 s21 and 250  1 lm2 s21, respectively. A range of smaller coefficients were observed in the different

2186

PROTEINS

MsbA 1 bilayer simulations: DDOPC 5 232–248 lm2 s21, DDOPE 5 220–232 lm2 s21. The reduction in lateral lipid diffusion in the presence of membrane protein has been observed in atomistic simulations,63 and is likely a general phenomenon for the concerted diffusion of solvent particles in the vicinity of a colloid. Fluorescence correlation spectroscopy was recently used to measure the lateral diffusion of membrane proteins within large vesicles.64 Fluorescent labeled integral membrane proteins similar in size to MsbA, such as MscL and MscS, exhibited diffusion constants of  4 lm2 s21 at concentrations below 103/lm2. The MSD of MsbA’s center of mass was monitored in each of four simulation runs [Fig. 9(E)]. Linearity was observed on the order of nanoseconds, enabling the determination of lower and upper bounds for the diffusion coefficient: DMsbA 5 7–10 lm2 s21. The factor of 2 effective speedup in translation of the large CG protein is relatively modest. Simulations with the MARTINI model have demonstrated a factor of 4 speedup in lateral diffusion for membrane-inserted WALP peptide.65 To examine phospholipid localization in the vicinity of protein, probability maps were computed for the positions of lipid head and tail group CG sites. A threedimensional grid was constructed in the MsbA reference frame from RMSD alignment of the TMDs in each snapshot and the local lipid density was computed relative to the bulk density of homogenous bilayer. Lipid density isosurfaces for the open, closed, and hydrolysis intermediate conformations are compared in Figure 10. For each structure, a single well-defined solvation shell is observed coating the transmembrane accessible protein surface including its internal chamber. Judging from the magnitude of the density relative to the bulk,66 the annular lipids are weakly bound with a favorable association free energy <1 kcal. The energy barrier for dissociation of annular lipids must be >1 kcal mol21 given the residence times observed. Other than for the lipids explicitly shown in Figure 10, which remained bound to the same protein site for most of the simulation, a given binding site exchanged lipids at intervals typically ranging from 10 to 100 ns, occasionally even 2 or 200 ns. Only a single phospholipid tail was usually bound, the other tail being completely flexible in the bilayer leaflet. The tail binding sites are found in grooves in the protein surface. In the case of the closed conformer, binding grooves are normal to the bilayer and reside between the exposed sidechains of adjacent TM helices [Fig. 10(B), inset]. For the two MsbA conformations in which the TMDs are not normal to the bilayer, distorted7 modes of binding transverse to the helix axes are present [Fig. 10(A,C)]. Of the lipids that remained bound throughout the simulation, some inserted between the protein backbones of adjacent TM helices. Functional roles for the specific binding of buried ‘‘nonannular’’ lipids have been discussed elsewhere,7 but their prevalence and nonspecific

Coarse Grain Lipid–Protein Molecular Interactions

Figure 10 Annular lipid density in CG simulations of open (A), closed (B), and hydrolysis intermediate (C) MsbA dimer conformers. Densities were computed in the reference frame of the protein, which freely diffused in the bilayer. Isosurfaces are drawn at 3.4 and 1.8 times the bulk density of lipid head (purple mesh) and tail (yellow mesh) groups, respectively. Shaded in light blue is the 1.9 head group density contour illustrating the uniform hydrophobic thickness of the bilayer excepting select local deformations. Hydrophobic, polar, positive, and negative protein residues are colored white, green, blue, and red, respectively. The amphipathic N-terminal elbow helix ending in Lys24 resides at different depths in the membrane in each conformation. Lipids that remain bound for the last half of the simulation are shown in ball and stick.

binding in this study may be an artifact of the reduced steric bulk of the low resolution model.67 In simulations of closed MsbA, a few nonannular lipids were even observed to make their way into the internal substrate binding chamber. Lastly, a nonannular lipid became statically bound under the TM1 elbow helix ending in Lys24 despite the fact that the N-terminus extends past the hydrophobic thickness of the bilayer. Single lipid tails were also transiently exposed from the bilayer, particu-

larly for the closed and hydrolysis intermediate structures that have more nonpolar residues sticking outside the membrane than the open conformer. Failure to delete all overlapping lipids during the embedding of protein into the bilayer setup resulted in a few lipids slowly diffusing along the aqueous protein surface in the water-free simulations. Distinct from the annular and nonannular lipids, nonradially symmetric clusters of head groups formed in the PROTEINS

2187

A.B. Ward et al.

vicinity of protein. Ellipsoid pockets of high probability can be identified spanning a distance of 5–15 A˚ from the cylindrical protein surface (Fig. 10). Some head group clusters penetrate the hydrophobic bilayer thickness, up to  5 A˚ in the case of the open conformer. The binding of phosphatidyl groups to basic sidechains at the membrane–water interface of voltage sensor domains was found to result in a 10 A˚ radial pattern of decreased hydrophobic thickness.9 The large clusters of zwitterionic choline and ethanolamine head groups did not interact directly with charged MsbA sidechains, but rather seem to coincide with regions of the protein surface that either are not normal to the bilayer or expose nonpolar sidechains to the head group region of the bilayer. Cluster location was thus dependent on MsbA conformation and TMD orientation. The lipid tails emanating from such pockets were dynamic but frequently adopted orientations parallel to the membrane–water interface to accommodate the irregular protein shape. Hydrophobic mismatch theory predicts a graded change in bilayer hydrophobic thickness to accommodate nonpolar protein surfaces.68 All CG simulations, however, exhibited a uniform bilayer thickness of 31 A˚ defined by the most probable distance between ester coordinates, closely matching MsbA’s 32 A˚ hydrophobic thickness. For complex protein architectures like MsbA including its amphipathic elbow helix, which changes orientation at the membrane interface during conformational cycling, hydrophobic matching and adaptation may proceed via local deformations rather than an alteration of bilayer thickness.8 Attachment of molecular reporter groups will shed light on such conformational interplay of polar and nonpolar accessibility at the membrane–water interface.32

CONCLUSION A reduced representation of sidechain sterics and polarity was developed and then combined with a waterfree model bilayer to investigate lipid–protein interactions in an ABC transporter. Elastic network models of diverse conformational states enabled molecular dynamics simulations of effectively long timescales and diffusive processes for both lipid and protein. While the stability of MsbA in the model membrane was demonstrated, the insertion of other protein architectures will need to be rigorously assessed, along with sidechain partitioning energetics.21,69 An extensive shell of weakly bound single lipid tails was observed to solvate the nonpolar protein surface. The average lipid exchange time between bulk and annulus was in accord with the original experimental measurement of 1027 s for spin labeled phospholipids binding to at least 22 annular sites on calcium–magnesium– ATPase.70 Disturbances in lipid bilayer structure were monitored in the vicinity of membrane protein as a func-

2188

PROTEINS

tion of conformation, TMD orientation and membrane accessible surface. The head group region exhibited local deformations at the membrane–water interface in response to nonpolar exposed surface and TM helix orientations deviating from the bilayer normal. Analogously, bilayer deformation caused by the irregular shape of a transmembrane protease has been implicated in the mechanism of substrate entry.71 Snorkeling of transmembrane arginine and lysine sidechains was seen to result in 4 A˚ decreased hydrophobic thickness out to a distance of 10 A˚ from the protein surface.72 Aquaporin has been cocrystallized with an annulus of polar lipids, as have a growing number of proteins such as KscA and adrenergic receptor with specifically bound anionic lipids or cholesterol crucial for their function.7,40,73,74 A question that remains is to what extent function is dictated by individual lipid binding events or general properties of the bilayer bulk and lipid composition. Phospholipid–protein interactions may play a competitive or cooperative role in polyspecific binding recognition and translocation events in lipid flippase transporters, and further exploration is needed to advance our understanding of multidrug efflux and resistance. The study presented here provides a platform for investigating the functional consequences of direct molecular lipid interactions with membrane proteins. ACKNOWLEDGMENTS RDH thanks the Voth Group for the MS-CG parameters and Lanyuan Lu for valuable discussion. Calculations utilized the NSF Rocks Cluster Toolkit. REFERENCES 1. Hills RD, Jr, Brooks CL, III. Insights from coarse-grained Go models for protein folding and dynamics. Int J Mol Sci 2009;10:889– 905. 2. Schug A, Onuchic JN. From protein folding to protein function and biomolecular binding by energy landscape theory. Curr Opin Pharmacol 2010;10:709–714. 3. Bahar I, Lezon TR, Bakan A, Shrivastava IH. Normal mode analysis of biomolecular structures: functional mechanisms of membrane proteins. Chem Rev 2010;110:1463–1497. 4. Apajalahti T, Niemela P, Govindan PN, Miettinen MS, Salonen E, Marrink SJ, Vattulainen I. Concerted diffusion of lipids in raft-like membranes. Faraday Discuss 2010;144:411–430. 5. Marrink SJ, de Vries AH, Tieleman DP. Lipids on the move: simulations of membrane pores, domains, stalks and curves. Biochim Biophys Acta-Biomembr 2009;1788:149–168. 6. Risselada HJ, Marrink SJ. The molecular face of lipid rafts in model membranes. Proc Natl Acad Sci USA 2008;105:17367–17372. 7. Lee AG. Biological membranes: the importance of molecular detail. Trends Biochem Sci 2011;36:493–500. 8. Sonntag Y, Musgaard M, Olesen C, Schiott B, Moller JV, Nissen P, Thogersen L. Mutual adaptation of a membrane protein and its lipid bilayer during conformational changes. Nat Commun 2011;2:304. 9. Mokrab Y, Sansom MS. Interaction of diverse voltage sensor homologs with lipid bilayers revealed by self-assembly simulations. Biophys J 2011;100:875–884.

Coarse Grain Lipid–Protein Molecular Interactions

10. Marrink SJ, Risselada HJ, Yefimov S, Tieleman DP, de Vries AH. The MARTINI force field: coarse grained model for biomolecular simulations. J Phys Chem B 2007;111:7812–7824. 11. Perlmutter JD, Sachs JN. Interleaflet interaction and asymmetry in phase separated lipid bilayers: molecular dynamics simulations. J Am Chem Soc 2011;133:6563–6577. 12. Schafer LV, de Jong DH, Holt A, Rzepiela AJ, de Vries AH, Poolman B, Killian JA, Marrink SJ. Lipid packing drives the segregation of transmembrane helices into disordered lipid domains in model membranes. Proc Natl Acad Sci USA 2011;108: 1343–1348. 13. Smirnova YG, Marrink SJ, Lipowsky R, Knecht V. Solvent-exposed tails as prestalk transition states for membrane fusion at low hydration. J Am Chem Soc 2010;132:6710–6718. 14. Parton DL, Klingelhoefer JW, Sansom MS. Aggregation of model membrane proteins, modulated by hydrophobic mismatch, membrane curvature, and protein class. Biophys J 2011;101:691–699. 15. Yesylevskyy SO, Schafer LV, Sengupta D, Marrink SJ. Polarizable water model for the coarse-grained MARTINI force field. PLoS Comput Biol 2010;6:e1000810. 16. Monticelli L, Kandasamy SK, Periole X, Larson RG, Tieleman DP, Marrink SJ. The MARTINI coarse-grained force field: extension to proteins. J Chem Theory Comput 2008;4:819–834. 17. Bereau T, Deserno M. Generic coarse-grained model for protein folding and aggregation. J Chem Phys 2009;130:235106. 18. Liwo A, He Y, Scheraga HA. Coarse-grained force field: general folding theory. Phys Chem Chem Phys 2011;13:16890–16901. 19. Bond PJ, Sansom MSP. Insertion and assembly of membrane proteins via simulation. J Am Chem Soc 2006;128:2697–2704. 20. Hills RD, Jr, Lu L, Voth GA. Multiscale coarse-graining of the protein energy landscape. PLoS Comput Biol 2010;6:e1000827. 21. Singh G, Tieleman DP. Using the Wimley-White hydrophobicity scale as a direct quantitative test of force fields: the MARTINI coarse-grained model. J Chem Theory Comput 2011;7: 2316–2324. 22. Bond PJ, Wee CL, Sansom MS. Coarse-grained molecular dynamics simulations of the energetics of helix insertion into a lipid bilayer. Biochemistry 2008;47:11321–11331. 23. Murtola T, Bunker A, Vattulainen I, Deserno M, Karttunen M. Multiscale modeling of emergent materials: biological and soft matter. Phys Chem Chem Phys 2009;11:1869–1892. 24. Izvekov S, Voth GA. Solvent-free lipid bilayer model using multiscale coarse-graining. J Phys Chem B 2009;113:4443–4455. 25. Wang ZJ, Deserno M. A systematically coarse-grained solvent-free model for quantitative phospholipid bilayer simulations. J Phys Chem B 2010;114:11207–11220. 26. Izvekov S, Voth GA. A multiscale coarse-graining method for biomolecular systems. J Phys Chem B 2005;109:2469–2473. 27. Noid WG, Chu JW, Ayton GS, Krishna V, Izvekov S, Voth GA, Das A, Andersen HC. The multiscale coarse-graining method. I. A rigorous bridge between atomistic and coarse-grained models. J Chem Phys 2008;128:244114. 28. Thorpe IF, Goldenberg DP, Voth GA. Exploration of transferability in multiscale coarse-grained peptide models. J Phys Chem B 2011;115:11911–11926. 29. Lu L, Voth GA. Systematic coarse-graining of a multicomponent lipid bilayer. J Phys Chem B 2009;113:1501–1510. 30. Periole X, Cavalli M, Marrink SJ, Ceruso MA. Combining an elastic network with a coarse-grained molecular force field: structure, dynamics, and intermolecular recognition. J Chem Theory Comput 2009;5:2531–2543. 31. Ward A, Reyes CL, Yu J, Roth CB, Chang G. Flexibility in the ABC transporter MsbA: alternating access with a twist. Proc Natl Acad Sci USA 2007;104:19005–19010. 32. Zou P, McHaourab HS. Alternating access of the putative substratebinding chamber in the ABC transporter MsbA. J Mol Biol 2009;393:574–585.

33. Campbell JD, Biggin PC, Baaden M, Sansom MS. Extending the structure of an ABC transporter to atomic resolution: modeling and simulation studies of MsbA. Biochemistry 2003;42:3666–3673. 34. Weng JW, Fan KN, Wang WN. The conformational transition pathway of ATP binding cassette transporter MsbA revealed by atomistic simulations. J Biol Chem 2010;285:3053–3063. 35. Aittoniemi J, de Wet H, Ashcroft FM, Sansom MS. Asymmetric switching in a homodimeric ABC transporter: a simulation study. PLoS Comput Biol 2010;6:e1000762. 36. Becker JP, Van Bambeke F, Tulkens PM, Prevost M. Dynamics and structural changes induced by ATP binding in SAV1866, a bacterial ABC exporter. J Phys Chem B 2010;114:15948–15957. 37. Gyimesi G, Ramachandran S, Kota P, Dokholyan NV, Sarkadi B, Hegedus T. ATP hydrolysis at one of the two sites in ABC transporters initiates transport related conformational transitions. Biochim Biophys Acta 2011;1808:2954–2964. 38. Jones PM, George AM. Molecular-dynamics simulations of the ATP/apo state of a multidrug ATP-binding cassette transporter provide a structural and mechanistic basis for the asymmetric occluded state. Biophys J 2011;100:3025–3034. 39. Oliveira AS, Baptista AM, Soares CM. Conformational changes induced by ATP-hydrolysis in an ABC transporter: a molecular dynamics study of the Sav1866 exporter. Proteins 2011;79: 1977–1990. 40. Hite RK, Li Z, Walz T. Principles of membrane protein interactions with annular lipids deduced from aquaporin-0 2D crystals. EMBO J 2010;29:1652–1658. 41. Zhou M, Morgner N, Barrera NP, Politis A, Isaacson SC, MatakVinkovic D, Murata T, Bernal RA, Stock D, Robinson CV. Mass spectrometry of intact V-type ATPases reveals bound lipids and the effects of nucleotide binding. Science 2011;334:380–385. 42. Hess B, Kutzner 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. 43. Kaminski GA, Friesner RA, Tirado-Rives J, Jorgensen WL. Evaluation and reparametrization of the OPLS-AA force field for proteins via comparison with accurate quantum chemical calculations on peptides. J Phys Chem B 2001;105:6474–6487. 44. Berger O, Edholm 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. 45. Martinez-Seara H, Rog T, Karttunen M, Reigada R, Vattulainen I. Influence of cis double-bond parametrization on lipid membrane properties: how seemingly insignificant details in force-field change even qualitative trends. J Chem Phys 2008;129:105103. 46. Lu L, Izvekov S, Das A, Andersen HC, Voth GA. Efficient, regularized, and scalable algorithms for multiscale coarse-graining. J Chem Theory Comput 2010;6:954–965. 47. Sugita Y, Okamoto Y. Replica-exchange molecular dynamics method for protein folding. Chem Phys Lett 1999;314:141–151. 48. de Vries AH, Mark AE, Marrink SJ. The binary mixing behavior of phospholipids in a bilayer: a molecular dynamics study. J Phys Chem B 2004;108:2454–2463. 49. Das A, Andersen HC. The multiscale coarse-graining method. V. Isothermal-isobaric ensemble. J Chem Phys 2010;132:164106. 50. Lomize MA, Lomize AL, Pogozheva ID, Mosberg HI. OPM: orientations of proteins in membranes database. Bioinformatics 2006;22:623–625. 51. Arnold K, Bordoli L, Kopp J, Schwede T. The SWISS-MODEL workspace: a web-based environment for protein structure homology modelling. Bioinformatics 2006;22:195–201. 52. Humphrey W, Dalke A, Schulten K. VMD: visual molecular dynamics. J Mol Graph 1996;14:33–38. 53. Beckstein O, Sansom MS. Liquid-vapor oscillations of water in hydrophobic nanopores. Proc Natl Acad Sci USA 2003;100:7063– 7068.

PROTEINS

2189

A.B. Ward et al.

54. Daily MD, Phillips GN, Jr., Cui Q. Interconversion of functional motions between mesophilic and thermophilic adenylate kinases. PLoS Comput Biol 2011;7:e1002103. 55. Chen J, Brooks CL, III. Implicit modeling of nonpolar solvation for simulating protein folding and conformational transitions. Phys Chem Chem Phys 2008;10:471–481. 56. Lammert H, Schug A, Onuchic JN. Robustness and generalization of structure-based models for protein folding and function. Proteins 2009;77:881–891. 57. Feig M. Implicit membrane models for membrane protein simulation. Methods Mol Biol 2008;443:181–196. 58. Zou P, Bortolus M, McHaourab HS. Conformational cycle of the ABC transporter MsbA in liposomes: detailed analysis using double electron-electron resonance spectroscopy. J Mol Biol 2009;393:586–597. 59. Kneller GR, Baczynski K, Pasenkiewicz-Gierula M. Consistent picture of lateral subdiffusion in lipid bilayers: molecular dynamics simulation and exact results. J Chem Phys 2011;135:141105. 60. Coppock PS, Kindt JT. Atomistic simulations of mixed-lipid bilayers in gel and fluid phases. Langmuir 2009;25:352–359. 61. Wohlert J, Edholm O. Dynamics in atomistic simulations of phospholipid membranes: nuclear magnetic resonance relaxation rates and lateral diffusion. J Chem Phys 2006;125:204703. 62. Rodgers JM, Sorensen J, de Meyer FJ, Schiott B, Smit B. Understanding the phase behavior of coarse-grained model lipid bilayers through computational calorimetry. J Phys Chem B 2012;116:1551– 1569. 63. Niemela PS, Miettinen MS, Monticelli L, Hammaren H, Bjelkmar P, Murtola T, Lindahl E, Vattulainen I. Membrane proteins diffuse as dynamic complexes with lipids. J Am Chem Soc 2010;132:7574–7575. 64. Ramadurai S, Holt A, Krasnikov V, van den Bogaart G, Killian JA, Poolman B. Lateral diffusion of membrane proteins. J Am Chem Soc 2009;131:12650–12656.

2190

PROTEINS

65. Ramadurai S, Holt A, Schafer LV, Krasnikov VV, Rijkers DT, Marrink SJ, Killian JA, Poolman B. Influence of hydrophobic mismatch and amino acid composition on the lateral diffusion of transmembrane peptides. Biophys J 2010;99:1447–1454. 66. Raman EP, Yu W, Guvench O, Mackerell AD. Reproducing crystal binding modes of ligand functional groups using site-identification by ligand competitive saturation (SILCS) simulations. J Chem Inf Model 2011;51:877–896. 67. Lammert H, Wolynes PG, Onuchic JN. The role of atomic level steric effects and attractive forces in protein folding. Proteins 2012;80:362–373. 68. Marsh D. Energetics of hydrophobic matching in lipid–protein interactions. Biophys J 2008;94:3996–4013. 69. Gumbart J, Roux B. Determination of membrane-insertion free energies by molecular dynamics simulations. Biophys J 2012;102: 795–801. 70. East JM, Melville D, Lee AG. Exchange rates and numbers of annular lipids for the calcium and magnesium ion dependent adenosinetriphosphatase. Biochemistry 1985;24:2615–2623. 71. Bondar AN, del Val C, White SH. Rhomboid protease dynamics and lipid interactions. Structure 2009;17:395–405. 72. Bond PJ, Sansom MS. Bilayer deformation by the Kv channel voltage sensor domain revealed by self-assembly simulations. Proc Natl Acad Sci USA 2007;104:2631–2636. 73. Hanson MA, Cherezov V, Griffith MT, Roth CB, Jaakola VP, Chien EY, Velasquez J, Kuhn P, Stevens RC. A specific cholesterol binding site is established by the 2.8 A structure of the human beta2-adrenergic receptor. Structure 2008;16:897–905. 74. Valiyaveetil FI, Zhou Y, MacKinnon R. Lipids in the structure, folding, and function of the KcsA K1 channel. Biochemistry 2002;41:10771–10777.