HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Full Text (PDF) Supplemental Research Data Alert me when this article is cited Alert me if a correction is posted Citation Map Services Email this article to a friend Similar articles in this journal Similar articles in PubMed Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via Google Scholar Google Scholar Articles by Lomize, A. L. Articles by Mosberg, H. I. Search for Related Content PubMed PubMed Citation Articles by Lomize, A. L. Articles by Mosberg, H. I. Social Bookmarking                   What's this?

Positioning of proteins in membranes: A computational approach

¹ College of Pharmacy
² College of Literature, Science and the Arts, University of Michigan, Ann Arbor, Michigan 48109-1065, USA

(RECEIVED January 30, 2006; FINAL REVISION February 24, 2006; ACCEPTED February 24, 2006)

	Abstract

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

 
A new computational approach has been developed to determine the spatial arrangement of proteins in membranes by minimizing their transfer energies from water to the lipid bilayer. The membrane hydrocarbon core was approximated as a planar slab of adjustable thickness with decadiene-like interior and interfacial polarity profiles derived from published EPR studies. Applicability and accuracy of the method was verified for a set of 24 transmembrane proteins whose orientations in membranes have been studied by spin-labeling, chemical modification, fluorescence, ATR FTIR, NMR, cryo-microscopy, and neutron diffraction. Subsequently, the optimal rotational and translational positions were calculated for 109 transmembrane, five integral monotopic and 27 peripheral protein complexes with known 3D structures. This method can reliably distinguish transmembrane and integral monotopic proteins from water-soluble proteins based on their transfer energies and membrane penetration depths. The accuracies of calculated hydrophobic thicknesses and tilt angles were 1 Å and 2°, respectively, judging from their deviations in different crystal forms of the same proteins. The hydrophobic thicknesses of transmembrane proteins ranged from 21.1 to 43.8 Å depending on the type of biological membrane, while their tilt angles with respect to the bilayer normal varied from zero in symmetric complexes to 26° in asymmetric structures. Calculated hydrophobic boundaries of proteins are located 5 Å lower than lipid phosphates and correspond to the zero membrane depth parameter of spin-labeled residues. Coordinates of all studied proteins with their membrane boundaries can be found in the Orientations of Proteins in Membranes (OPM) database:http://opm.phar.umich.edu/.

Keywords: membrane protein; tilt; database; modeling; channel; receptor; transporter; rhodopsin

	Introduction

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

 
Thousands of membrane-associated proteins have been deposited in the Protein Data Bank (PDB; Berman et al. 2000), and their number is rapidly growing. However, the precise spatial positions of these proteins in membranes are unknown. Membrane proteins are unique because they function in the highly anisotropic environment of a lipid bilayer, which is characterized by complex polarity gradients and a heterogeneous molecular composition in different regions and leaflets. The positioning of proteins in the lipid matrix may affect their biological activity, folding, thermodynamic stability, and binding with surrounding macromolecules and substrates (White and Wimley 1999; Booth et al. 2001; Bowie 2001; DeGrado et al. 2003; Engelman et al. 2003; Hong and Tamm 2004; Lee 2004). Hence, the orientations of many peptides and proteins in membranes have been studied by a variety of experimental techniques including chemical modification, spin-labeling, paramagnetic or fluorescence quenching, X-ray scattering, neuron diffraction, electron cryomicroscopy, NMR, and polarized infrared spectroscopy (Frillingos et al. 1998; Hristova et al. 1999; London and Ladokhin 2002; de Planque and Killian 2003; Hubbell et al. 2003; Opella and Marassi 2004; Tatulian et al. 2005). However, since the amount of such experimental data is limited, this problem should also be addressed computationally to keep up with the expanding flow of structures in the PDB.

The arrangement of a protein with respect to the membrane can be defined by its shift along the bilayer normal (d), rotational and tilt angles ( and ), and thickness of its membrane-spanning region (D = 2z₀, Fig. 1). Although the orientation of a transmembrane (TM) protein in a lipid bilayer can be assessed manually (Lee 2003), development of automated methods is necessary to provide more objective, reproducible, and accurate results. The existing computational approaches range from elaborate molecular dynamic (MD) simulations of proteins with explicit water and lipids (Ash et al. 2004; Roux et al. 2004; Gumbart et al. 2005) to simplified approaches that minimize protein transfer energy from water to a hydrophobic slab, which serves as a crude approximation of the membrane hydrocarbon core (Yeates et al. 1987; Rees et al. 1989). In the latter case, the transfer energy can be estimated by using various hydrophobicity scales of whole residues (Zucic and Juretic 2004), a normalized nonpolar accessible surface area (Tusnady et al. 2004), or atomic solvation parameters derived from partition coefficients of model organic compounds between water and nonpolar solvents (Basyn et al. 2003). One such computational method has recently been applied to create the PDB_TM database that provides an up-to-date list of all TM peptides and proteins from the PDB (Tusnady et al. 2005). However, the hydrophobic boundaries of proteins in PDB_TM were not compared with relevant experimental studies and their accuracy is uncertain.

View larger version (51K): [in this window] [in a new window]   Figure 1. Schematic representation of a TM protein in a hydrophobic slab. Parameters that define the arrangement of a TM protein in the membrane hydrocarbon core: d, shift along the bilayer normal; , tilt angle; , rotational angle; and D = 2z0, hydrophobic thickness of the protein.

	Results

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

 
Development of the method
Choice of atomic solvation parameters
The results for TM proteins were strongly dependent on the choice of atomic solvation parameters that have been applied for calculations of protein transfer energy (Table 1). These parameters can be derived from partition coefficients of uncharged solutes between water and octanol, cyclohexane, or other nonpolar solvents (Eisenberg and McLachlan 1986; Ducarme et al. 1998; Efremov et al. 1999; Lomize et al. 2004). To choose the appropriate parameter set, we calculated the orientations of TM proteins by our program using three alternative scales and compared the results with experimental studies (Table 2). The comparison shows that 1,9-decadiene and cyclohexane-based parameters produce nearly identical results that are consistent with the experimental data, while the octanol scale performs poorly. In all subsequent calculations, we used parameters for decadiene, because this solvent was shown to be the best model of the bilayer interior in studies of membrane permeability barriers (Xiang and Anderson 1994a, b; Mayer and Anderson 2002). We also found that two slightly different parameter sets should be applied for proteins in detergents and bilayers (Table 1). The results of the calculations with decadiene-based parameters ("lipid bilayer scale") were more consistent with experimental studies for mechanosensitive MscL channel, F-type Na⁺-ATP synthase, and rhodopsin in bilayers, while hexadecene-based parameters ("detergent scale") were more suitable for reproducing EPR data for rhodopsin in dodecyl maltoside (see below).

Treatment of internal cavities in TM proteins
It was also important that the transfer energy does not include any contributions of atoms that face internal polar cavities of TM proteins and do not directly interact with surrounding bulk lipid. Without this, the orientations of -barrels and many -helical transporters would be calculated incorrectly, because they have large interior channels or funnels filled by polar and charged residues. This problem was recognized and addressed previously in calculations with TMDET and Garlic (Tusnady et al. 2004; Zucic and Juretic 2004) and in Monte Carlo simulations of -barrels (Basyn et al. 2003). Our algorithm for excluding atoms in the funnels or channels was different from that in the previously published methods (see Materials and Methods).

Main features of the method
The following approximations were found to be necessary and sufficient for reproducing the experimental data: (1) a lipid bilayer is represented as a planar hydrophobic slab with adjustable thickness and a narrow interfacial area with a sigmoidal polarity profile; (2) a protein is considered as a rigid body with flexible side chains whose transfer energy is minimized with respect to four variables; (3) transfer energy is calculated at an all-atom level using atomic solvation parameters determined for the water-decadiene system; (4) explicit electrostatic interactions were neglected, while including ionization penalties for charged residues that are considered neutral in the nonpolar environment; (5) contributions of pore-facing atoms in TM proteins are automatically eliminated. The computational model obtained is exceptionally simple, because it depends only on five atomic solvation parameters (for N, O, S, and sp2 and sp3 carbons), a constant defining the size of the interfacial region, and ionization energies of charged groups. None of these parameters are adjustable, but rather, independently derived from various experimental sources. This method was implemented into the program PPM 1.0 (Positioning of Proteins in Membranes).

Precision of the method
Parameters from different crystal forms of the same protein were calculated to estimate the precision of the method (Supplemental Material, Table 1). Deviations of parameters D and were within 1 Å and 2°, respectively, in complexes with identical numbers of TM subunits. Deviations were primarily caused by different conformations of flexible side chains and nonregular loops in crystal structures. Orientations of proteins can fluctuate significantly if their energy surfaces are shallow. The corresponding maximal variations (±) of D and parameters were calculated within 1 kcal/mol around the global minimum of transfer energy, as has been suggested previously (Yeates et al. 1987). Obtained ranges are indicated in all Tables. They are larger than the precisions of the corresponding parameters estimated from comparison of different crystal forms.

Verification of the method
Results were verified through all available experimental data for 24 TM proteins of known 3D structure whose spatial positions in bilayers have been experimentally studied: rhodopsin (1gzm), bacteriorhodopsin (1py6), sensory rhodopsin II (1h2s), photosynthetic reaction centers from two species (1rzh and 1dxr), cytochrome c oxidase (1v55), Na⁺-ATPase (1yce), Ca²⁺-ATPase(1iwo, 2agv, 1wpe, 1t5s, 1su4, 1wpg), phospholamban (1zll), lactose permease LacY (1pv6), protein translocase SecY (1rh5), Na⁺/H⁺ antiporter (1zcd), K⁺-channel KcsA (1r3j), MscL mechanosensitive channel (1msl), acetylcholine receptor (2bg9), outer membrane proteins OmpA (1qjp [PDB] ), OmpX (1qj8), OmpLA (1qd6), OmpF (1hxx), ferric enterobactin receptor receptor FepA (1fep), ferric hydroxamate uptake receptor FhuA (1qfg), cobalamine transporter BtuB (1nqe), -hemolysin (7ahl), and gramicidin A (1grm). A number of methods, such as chemical modification, spin-labeling, fluorescence spectroscopy, ATP FTIR, NMR, X-ray scattering, neutron diffraction, electron cryo-microscopy, and hydrophobic matching studies were used to determine hydrophobic thicknesses or tilts of these proteins, to locate their membrane-embedded segments, and to evaluate penetration depths and environments of their residues in lipid bilayers or detergents (Supplemental Material).

Comparison with hydrophobic thicknesses of matching lipid bilayers
The calculated hydrophobic thicknesses of 12 TM proteins agree with the corresponding experimental values obtained from site-directed spin-labeling studies (BtuB transporter and FepA receptor), EM data (Na⁺ ATPase), X-ray scattering of photoreceptor membranes (rhodopsin) or hydrophobic matching experiments (other proteins), as can be seen from comparison of D_calc(dcd) and D_exp values in Table 2. The hydrophobic matching studies determine the bilayer thickness that provides optimal functional activity of a protein (Lee 2004) or maximum protein–lipid binding affinity (Lee 2003), or identify lipids whose temperature of phase transition is not affected by the presence of the protein (Dumas et al. 1999). The hydrocarbon thicknesses of lipid bilayers are obtained by subtracting 10 Å from their phosphate-to-phosphate distances determined by X-ray scattering (Lewis and Engelman 1983a; Nagle and Tristram-Nagle 2000). Thus, the calculated hydrophobic boundaries are located 5 Å from the phosphate groups toward the membrane center, i.e., at the level of the carbonyl groups of the lipid molecules.

Comparison with experimental tilt angles
NMR studies of the bacteriorhodopsin trimer show that helix A and the extracellular section of helix B are tilted with respect to the bilayer normal by 18°–22° and by less than 5°, respectively (Kamihira et al. 2005). This is consistent with the calculated tilt angles of 23° and 5°, respectively, for these helices in a trimer (1qm8). The calculated tilt of gramicidin A channel (2° ± 10°) is also in excellent agreement with solid-state NMR and infrared dichroism studies that show a nearly perpendicular arrangement of the dimer in the membrane (Nabedryk et al. 1982; Andronesi et al. 2004; Andersen et al. 2005). The calculated average tilts of -helices or -strands in TM proteins correlate well with ATR FTIR spectroscopy data (Table 3). However, the experimental values are systematically larger, which could be due to some orientational disorder under the experimental conditions. It has been noted that values of obtained by ATR FTIR spectroscopy may represent upper limits of the actual tilt angles for -helical peptides, due to such disorder (Bechinger et al. 1999; de Planque and Killian 2003).

View this table: [in this window] [in a new window]   Table 3. Average tilt angles (°) of TM -helices or -strands relative to the membrane normal calculated by PPM 1.0 (calc) and determined by ATR FTIR spectroscopy (exper)

Comparison with membrane penetration depths of individual residues
The calculated membrane-embedded portions of the regular secondary structures agree with studies of BtuB transpoter, bacteriorhodopsin, FepA receptor, and MscL and KcsA channels by spin-labeling, MscL by fluorescence, and Na+ ATPase by EM and X-ray crystallography (Table 4). The experimental and calculated penetration depths of individual spin-labeled residues are generally consistent for MscL and KcsA channels, bacteriorhodopsin, and FepA receptor (Fig. 2). However, experimental and calculated depths may deviate up to 5 Å, as observed for residue 69 in MscL channel (D_calc = 1.1 Å and D_exper = 6.0 Å). Such deviations probably appear because the depths were taken for the C-atom of the spin labeled cysteine instead of the nitroxyl group, which actually interacts with the paramagnetic quenchers. The penetration depth of the nitroxyl radical in Cys69 of MscL can be in the range of –1.0 Å to 6.5 Å, depending on four angles of the spin-labeled cysteine. The latter value is more consistent with the experiment.

View larger version (20K): [in this window] [in a new window]   Figure 2. Comparison of calculated and experimental membrane penetration depths of spin-labeled Cys residues in TM proteins. Five TM proteins were studied: MscL channel (1msl, open square), FepA receptor (1fep, black triangle), KcsA channel (1r3j, gray triangle), bacteriorhodopsin (1py6, black diamond), and BtuB porin (1nqe, black circle). The experimental distances are taken from the original publications (Altenbach et al. 1990, 1994; Greenhalgh et al. 1991; Klug et al. 1997; Perozo et al. 1998, 2001; Fanucci et al. 2002) but counted relative to the depth where parameter is equal to zero. Points with zero depth correspond to residues that have been identified as first or last in the membrane-embedded segments of -helices of -strands. Calculated depths are defined as distances from C-atoms of the corresponding residues to the closest boundary plane.

View this table: [in this window] [in a new window]   Table 4. Comparison of calculated and experimentally determined membrane-embedded portions of -helices or -strands

Comparison with environments of residues
The water- or lipid-facing environments of individual residues in TM proteins can be mapped by different chemical probes. Lactose permease LacY (1pv6) has been more extensively studied by chemical modification than any other TM protein. A total of 393 single-Cys mutants of LacY were modified by bifunctional reagents to identify residues that could be involved in intermolecular cross-linking (Guan et al. 2002; Ermolova et al. 2003). Only residues located in regions of sufficiently high polarity could produce the reactive thiolate anion required for the formation of intermolecular disulfide bonds. Indeed, our calculations show that most of the residues susceptible to cross-linking are situated in water-exposed periplasmic and cytoplasmic loops or in a relatively narrow (5 Å) layer of the hydrophobic core where the dielectric permittivity may be intermediate between that in lipid and water (Fig. 3A). These layers are parallel to the calculated interfacial planes, thus confirming that these planes were identified correctly. However, some of the reactive residues (Trp78, Phe398) are situated 7–10 Å from the surface, which is most probably because they occupy polar sites at the C-termini of TM helices close to Lys74, Lys188, and Lys289 where local dielectric constant may be higher.

View larger version (47K): [in this window] [in a new window]   Figure 3. Comparison of calculated membrane core boundaries with experimental data for lactose permease (1pv6) in lipid bilayer (A,B), for rhodopsin (1gzm) in native membrane (C) and in detergent (D), and for OmpX (1qj8) in detergent (E). The boundaries are indicated by blue dots for the inner membrane side and red dots for the outer membrane side. (A) Cys-substituted residues in lactose permease, which were modified by cross-linking thiosulfonate agents (C-atoms are colored red). (B) Residues of lactose permease that can be modified by NEM are colored red, residues inaccessible to NEM are colored blue, and lipid-accessible residues identified in spin-labeling studies are colored green. (C) Residues of rhodopsin that were modified by hydrophilic (red) and hydrophobic (blue) chemical probes in native photoreceptor membranes. (D) Spin-labeled residues of rhodopsin in detergent assigned to nonpolar and polar environments by EPR are colored blue and red, respectively (residues with undefined environment are colored gray). (E) NH and CH3 groups of OmpX that form NOEs with hydrophobic tails of DHPC are colored blue and green, respectively. NH groups that interact with head groups of DHPC are colored red. (F) A schematic representation of a protein in a native membrane (left) or in detergent (right).

Vertebrate rhodopsins have also been studied in great detail. Environments of many rhodopsin residues were characterized in intact photoreceptor membranes (Davison and Findlay 1986a, b) or in n-dodecyl--D-maltoside (DM) (Hubbell et al. 2003). The results are consistent with the data for the native membranes, primarily with the results of chemical modification of ovine rhodopsin by hydrophobic and hydrophilic probes, which interact with residues exposed to nonpolar or polar environments, respectively (Barclay and Findlay 1984; Davison and Findlay 1986a, b). The hydrophobic probe was shown to modify all 14 Cys, Trp, Tyr, and His residues within the calculated hydrophobic boundaries (blue in Fig. 3C) and a few residues situated just outside the boundary in the lipid head group area (His65, Lys66, Tyr74, Lys231, Cys316) (Davison and Findlay 1986a). All polar residues that were modified by nonpermeable hydrophilic probe applied from the intracellular side (Barclay and Findlay 1984) are located outside the calculated hydrophobic slab (red in Fig. 3C).

The accessibilities to paramagnetic quenchers of two residues in Na⁺/H⁺ antiporter (1zcd) and eight residues in sensory rhodopsin II (1h2s) (Wegener et al. 2000; Hilger et al. 2005) are consistent with locations of the membrane boundaries calculated for the corresponding proteins.

Comparison with studies of TM proteins in detergents
The hydrophobic dimensions of TM proteins have also been studied in detergents using neutron diffraction with contrast variation in crystals (two photosynthetic reaction centers, trimeric porin OmpF and monomeric phosphoripase OmpLA), solution NMR (OmpX) and spin-labeling (rhodopsin). Detergent molecules form monolayers around nonpolar surfaces of TM proteins and most of them are oriented perpendicular to the protein surface, unlike lipids in bilayers (Fig. 3F; le Maire et al. 2000). In crystals of TM proteins, these monolayers look like regular rings with dimensions of 15–20 Å in the direction perpendicular to TM domains and 15–30 Å parallel to them (Roth et al. 1989, 1991). The latter values are in agreement with hydrophobic thicknesses of the corresponding proteins calculated by PPM 1.0 (Table 5). Moreover, the calculated membrane boundary planes of the proteins closely correspond to the borders of the detergent monolayer. For example, these planes pass through the aromatic rings of Trp78, Trp98, Phe109, and Phe122 residues in monomeric phospholipase A (OmpLA) in agreement with neutron diffraction (Snijder et al. 2003).

Results of the calculations with "detergent" and "bilayer" scales (Table 1) were nearly identical for all proteins studied in detergents, except rhodopsin. To reproduce the experimental conditions, we used the crystal structure of rhodopsin with extended helix V (1gzm), removed its C-terminal palmitates, and applied the "detergent" solvation parameters. This led to an increased hydrophobic thickness (from 32.4 Å to 36.9 Å) and tilt angle (from 8° to 16°) (Fig. 3D). The expanded membrane boundaries are in much better agreement with spin-labeling data in DM (Hubbell et al. 2003) than with chemical modification data in native membranes (Davison and Findlay 1986a, b). The EPR studies identified a number of interfacial residues that were buried from water when substituted by a spin-labeled cysteine (V63, P71, V137, H152, K231, T251, and N310). The C-atoms of these residues are indeed situated within the expanded boundaries calculated with "detergent" parameters, but well outside the boundaries obtained with "membrane" parameters (Supplemental Material, Table 3). Most importantly, several residues in the last turn of helix V (227–231) were shown to be coated with the detergent (Hubbell et al. 2003), but are accessible to water in the native membrane.

Application of the method to TM proteins from the PDB
After successful testing of the method for 24 well-studied TM proteins, it was applied to all other TM proteins deposited in the PDB. The calculations were conducted for 109 TM protein complexes (80 -helical, 28 -barrels, and gramicidin A dimer), 32 representative integral monotopic and peripheral proteins selected from the literature, and a control set of 20 water-soluble proteins with the highest hydrophobicity score in PDB_TM. Any protein that did not traverse the membrane after the optimization was interpreted as peripheral or monotopic, and the maximal membrane penetration depth of its atoms is calculated instead of its hydrophobic thickness.

Figure 4 shows that TM, integral monotopic and peripheral proteins occupy separate areas of the plot of hydrophobic thickness (D = 2z₀) versus transfer energy (G_transf), and therefore, they can be easily distinguished based on these two parameters. All peripheral and monotopic proteins have penetration depths of <15 Å. Sixteen hydrophobic water-soluble proteins from PDB_TM have transfer energies in the range from 0 to –0.5 kcal/mol. However, there are four exceptions that can be interpreted as probable membrane-associated proteins. Among them are V_H antibody domain resistant to aggregation (1ohq, G_transf = –7.2 kcal/mol), extracellular domain of bone morphogenetic protein receptor (1es7, G_transf = –5.3 kcal/mol), gephyrin domain that links glycine receptor and tubulin (1t3e, G_transf = –3.4 kcal/mol), and caseine kinase (1rqf, G_transf = –1.6 kcal/mol). Indeed, interactions with membranes were shown to be functionally important for gephyrin and bone morphogenetic receptor (Sebald et al. 2004; Sola et al. 2004). Thus, our program can be applied for automatic identification and discrimination of TM, integral monotopic and water-soluble proteins, although the threshold between peripheral and nonmembrane proteins is sometimes blurred.

View larger version (18K): [in this window] [in a new window]   Figure 4. Hydrophobic thickness or membrane penetration depth (Dcalc) vs. transfer energy (Gtransf) plot for transmembrane (-helical and -barrel), integral monotopic, and peripheral proteins.

D and parameters of TM proteins were prone to fluctuations within 1 kcal/mol around the global minimum of transfer energy. These fluctuations were usually smaller than 2 Å and 4°, respectively. However, the fluctuations were larger for proteins with a smaller TM perimeter (Fig. 5).

View larger version (18K): [in this window] [in a new window]   Figure 5. Dependence of fluctuations of hydrophobic thickness (A) and tilt angle (B) on the size of TM proteins. Transfer energy divided by hydrophobic thickness is roughly proportional to the length of the outer perimeter in TM proteins.

In contrast, the hydrophobic thicknesses of F- and V-type ATPases and lipid flippases were unusually large, 36 Å (1a91, 1c17, 1yce, 2bl2, 1pf4, and 1z2r). This is expected to produce a significant hydrophobic mismatch between these proteins and their host bilayers. The mismatch may facilitate large-scale movements of these proteins in the lipid bilayers, because it reduces protein–lipid binding affinities (Lee 2003).

Some potential problems can also be detected based on the calculated tilt angles. These angles were close to zero for all TM complexes with noncrystallographic symmetry except cytochrome b6f from M. laminosus (1vf5) and fumarate reductase dimer from E. coli (1kf6), which both have of 3°. The non-zero overall tilt of cytochrome b6f complex appears due to significant differences between the symmetry-related subunits (RMSD of C-atoms 1 A). The dimer of fumarate reductase is loosely packed and therefore was suggested to be non-native (Iverson et al. 1999).

It is noteworthy that some TM proteins can form non-native dimers or trimers in crystals, for example, rhodopsin (1u19 and 1gzm), bacteriorhodopsin (1py6), lactose permease (1pv6), OmpA (1qjp), OmpX (1qj8), and fatty acid transporter FadL (1t16 and 1t1l). Hydrophobic thicknesses of such non-native oligomers are usually reduced (see Supplemental Material). Therefore, it is important to know the complete and correct quaternary structure of a multimeric complex to calculate its position in the membrane.

	Discussion

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

 
The positioning of proteins in membranes is an important problem that has been previously addressed by various methods. Our computational approach to this problem is exceptionally simple, because it neglects all fine details of the protein–lipid interactions in the membrane interfacial areas and accounts only for the hydrophobic burial of the proteins in the membrane core. We found that it is important to implement an appropriate atomic solvation parameter set, include ionization penalties for charged residues, and exclude contributions of pore residues that are not involved in interactions with surrounding bulk lipids. Since the results of the calculations are generally consistent with the experimental studies of 24 TM proteins, the main underlying assumptions of the method seem to be reasonable, as discussed below.

Description of the membrane
The bilayer core was approximated by a hydrophobic slab with narrow interfacial areas, based on EPR studies of spin-labeled lipid analogs (Marsh 2002). This crude approximation is appropriate, because the positioning of TM proteins in membranes is dominated by hydrophobic forces (Booth et al. 2001; Engelman et al. 2003), although electrostatic and other interactions at the lipid interface can also play an important role, especially for weakly bound peripheral peptides and proteins (Murray et al. 1997; White and Wimley 1999; Cho and Stahelin 2005). The calculated boundary planes of the proteins are located 5 Å deeper than lipid phosphates judging from the comparison of hydrophobic thicknesses of the proteins and lipid bilayers (Table 2). The effective concentrations of polar and nonpolar paramagnetic probes are approximately equal at these boundaries ( = 0) judging from the spin-labeling studies of TM proteins (Fig. 2). The borders between nonpolar and polar regions at the surfaces of protein complexes are usually well approximated by planes. Only a few nonpolar residues remained outside the calculated planar boundaries and few charged groups penetrated inside the calculated membrane-spanning regions. Some tendency for local hydrophobic thinning was observed in the structures of chloride channels and mitochondrial cytochrome c oxidase.

The concept of hydrophobic matching
The hydrophobic thickness was applied as a free variable, because the lipid bilayers are known to adjust their local thickness to match the nonpolar areas of large TM proteins (Dumas et al. 1999). Such adjustment costs relatively little energy, 0.025 kcal/mol per fatty acyl chain carbon atom (Williamson et al. 2002; Lee 2003), compared to 0.7 kcal/mol required for transfer of a CH₂ group from nonpolar solvent to water. The corresponding contractions or expansions of bilayers upon addition of TM proteins have been observed in biological membranes (Mitra et al. 2004). Our results demonstrate that TM proteins from inner bacterial membranes have more diverse hydrophobic thicknesses (from 22.5 Å to 37 Å) than TM proteins from other membrane types (Table 6). This may be due to a higher elasticity of the inner bacterial membranes facilitated by their heterogeneous lipid composition and lack of cholesterol (Denich et al. 2003). Indeed, the depletion of Escherichia coli membranes of proteins results in a substantial decrease of the thickness (from 27.5 Å to 23.5 Å; Table 6), while the parameters of the apical plasma membranes remained constant (Mitra et al. 2004). Thus, certain types of biological membranes may be relatively rigid due to the presence of cholesterol or dolichol/sphingomyelin (apical plasma membranes), lipids with phytanyl chains (archeabacteria), or lipopolysaccarides (outer membranes of Gram-negative bacteria). In all such cases, one could expect more uniform hydrophobic thicknesses in all proteins from the same type of membranes, as is actually observed in our results.

Implicit solvation model and hydrophobicity scale
An implicit solvation approach was appropriate because the interior of lipid bilayers is fluid. Transfer energy of TM proteins from water to the membrane interior was calculated using atomic solvation parameters developed for the water-decadiene system. The decadiene can serve as a good approximation of the bilayer interior, as follows from the experimental studies of membrane permeability barriers (Xiang and Anderson 1994a, b; Mayer and Anderson 2002), our previous work (Lomize et al. 2004), and results described here. In contrast, the experimental hydrophobic thicknesses of TM proteins were reproduced poorly with the octanol scale (Table 2). This is not surprising, since the octanol solution contains a significant amount of water (2.3M; Leo et al. 1971) and has high dielectric constant (10–20), unlike the hydrocarbon core of lipid bilayers (Benz et al. 1975; Radzicka and Wolfenden 1988). It was noted that water–octanol transfer energies of neutral solutes correlate poorly with membrane permeability barriers (Walter and Gutknecht 1986).

Comparison with other computational methods
Our method, PPM 1.0, represents a significant improvement upon other simple methods, such as optimization of normalized nopolar accessible surface areas (TMDET; Tusnady et al. 2004), Monte Carlo simulations with atomic solvation parameters (IMPALA; Basyn et al. 2003), or manual assessment (Lee 2003), because it provides results that are in better agreement with experimental data (Table 7). Surprisingly, the manually assessed hydrophobic thicknesses are close to those calculated by PPM 1.0. There are just a few cases of moderate discrepancies between the manually assessed and PPM's values, such as bacteriorhodopsin, bc1 complex, MscL channel, and Ca²⁺-ATPase. However, the parameters calculated by IMPALA and especially by TMDET (PDB_TM) differ significantly from ours, sometimes by more than 10 Å. It is noteworthy that the thicknesses of all -barrels from Gram-negative bacteria are very uniform when calculated with PPM 1.0 (22–25 Å), as was previously suggested (Tamm et al. 2004). However, the -barrel thicknesses calculated with IMPALA vary from 14 Å to 29 Å, while those in the PDB_TM database vary from 17.5 Å to 34 Å. Comparison of our method with the results of MD simulations with explicit lipids is less straightforward, because MD does not produce the hydrophobic thickness as an intrinsic property of a protein. The dynamically averaged tilt angles of TM proteins in MD studies are more or less consistent with our results. For example, the tilt of OmpA -barrel was 8° ± 6° and 5°–10° when generated by PPM 1.0 and MD (Bond et al. 2002), respectively. However, the agreement was less than perfect for the -barrel of OmpT and -helical dimer of glycophorin A. Both proteins are oriented nearly perpendicularly to the membrane plane in our results ( = 2° ± 5° and = 2° ± 10°, respectively), but are more tilted in MD studies ( = 20° and = 8°–10°, respectively) (Petrache et al. 2000; Baaden and Sansom 2004).

	Conclusion

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

	Materials and methods

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

 
Calculation of transfer energy
A protein was considered as a rigid body that freely floats in the planar hydrocarbon core of a lipid biliayer. The transfer energy of a protein, G_transfer, was calculated as a function of variables d, z₀, , and in a coordinate system whose Z-axis coincides with the membrane normal (Fig. 1):

where ASA_i is accessible surface area of atom i, _i^W–M is the solvation parameter of atom i (transfer energy of the atom from water to membrane interior expressed in kcal/mol per A²), and f(z_i) is interfacial water concentration profile. ASA were determined using the subroutine SOLVA from NACCESS (obtained from S.J. Hubbard and J.M. Thornton, University College London), with radii of Chothia (1975) and without hydrogens.

Solvation parameters have been derived specifically for lipid bilayers (Table 1) and normalized by the effective concentration of water, which changes gradually along the bilayer normal in a relatively narrow region between the lipid head group region and the hydrocarbon core. We used a sigmoidal water concentration profile as determined in EPR studies on spin-labeled phospholipids (Marsh 2002):

The characteristic distance of this profile was chosen as 0.9 Å (see Results).

Transfer energy was calculated only for atoms that are exposed to water or lipid at the outer perimeter of TM protein complexes. Contributions from exposed pore-facing atoms were excluded. The "lipid-facing" atoms were automatically defined as those that are not staggered by other protein atoms when "looked at" from the inner longitudinal axis of the TM -bundle or -barrel. This axis was calculated as the sum of individual TM secondary structure vectors whose signs were chosen to provide an identical general orientation of all vectors relative to the membrane. The protein was first rotated to superimpose the found longitudinal axis with the Z-axis. Then, each atom, i, was considered as staggered by atom j if the following conditions were satisfied. (1) |z_i – z_j| < 2 Å; (2) r_j > r_i + 2 Å, where r_j and r_i are distances from atoms i and j to Z-axis, and (3) the distance from atom j to a line, which is perpendicular to the Z-axis and passes through atom i, is <3.6 Å (sum of two van der Waals radii). For large TM complexes, one central axis "to look from" was insufficient; it was translated into different positions around the central point (x = 0, y = 0), and each atom that was not staggered by other atoms while looking from at least one position of the axis was assigned as lipid-facing.

Ionization/protonation energies of charged residues were described by the Henderson–Hasselbalch equation (Fersht 1999).

These energies are small at pH = 7: 6.9, 4.7, 4.9, 4.0, and 0.6 kcal/mol for Arg, Lys, Asp, Glu, and His, respectively, as follows from average pK_a values of these residues in proteins, i.e., 12.0, 10.4, 3.4, 4.1, and 6.6, respectively (Fersht 1999; Edgcomb and Murphy 2002; Forsyth et al. 2002). The ionization penalties were included as a part of the solvation parameters of N or O atoms in the charged groups, except those that are buried in the protein interior, because interactions and ionization states of the buried charged groups are not expected to change during immersion into the lipid phase. An ionizable group was considered as buried if it had ASA of polar atoms less than 1 Å, or if it formed at least two hydrogen bonds in the protein structure.

Global energy optimization
The transfer energy surfaces are typically shallow and have numerous local minima. Therefore, we used a reliable deterministic search strategy that included two steps: (1) grid scan to determine a set of low-energy combinations of variables z₀, d, , and (those with relative energies <4 kcal/mol); and (2) local energy minimization by the Davidon-Fletcher-Powell method starting from each low-energy point. Partial derivatives of the energy with respect to z₀, d, , and variables were calculated analytically. The appropriate grid scan steps were chosen empirically as 0.5 Å for z₀ and d, and 5° and 2° for and , respectively. Performance of the global search procedure was verified by starting from different spatial orientations of the complexes to make sure that they converge to essentially the same solution.

After global energy minimization, we determined minimal and maximal values of each variable that could be achieved in an interval of relative energies from 0 to 1 kcal/mol (x_min and x_max). This was done using values of energy calculated and stored during a grid scan and subsequent local energy minimizations. The x_min and x_max parameters provided an amplitude of possible fluctuations of the variable around the position of global minimum (x₀). The amplitudes were calculated as (x_max – x_min) and then used as ± values. Thus, they represent maximal rather than mean deviations.

We found that the results of calculations can be significantly affected by the orientation of interfacial side chains that are charged, accessible to the solvent, and situated close to the calculated boundaries. Many such side chains are flexible and poorly resolved (have high B-values), or even omitted in PDB entries. These side chains can easily change orientation upon association with the membrane. At present, there are no proven computational methods for automated side chain packing at the membrane interface, which must be done simultaneously with optimizing the orientation of the entire protein in the bilayer. Therefore, rearrangement of side chain rotamers was done semiautomatically (using PPM 1.0 and molecular modeling modules of QUANTA) as follows. First, the protein orientation was calculated for the PDB structure with original side chain rotamers. Secondly, we identified all solvent-exposed charged side chains (usually Lys and Arg) that were situated close to the calculated boundaries and could be rotated away from the hydrophobic core without creating interatomic overlaps or reducing the number of hydrogen bonds ( angles were kept close to ±60° or 180°). Remarkably, the majority of charged side chains were already turned away from the hydrophobic slab in the crystal structures of TM proteins, so only a few of them needed to be adjusted (see Supplemental Material). Third, conformations of the found side chains were readjusted, all charged residues with missing atoms near the boundaries were reconstructed, and membrane boundaries were recalculated. Fourth, the resulting structure was accepted if it had lower transfer energy than the original structure.

Protein data set
Structures of TM proteins were taken from the current release of PDB (Supplemental Material, Table 5). We temporarily omitted theoretical models, entries with backbone coordinates only, some EM-based structures, NMR models derived from orientational constraints, and nonfunctional structures (such as peptide fragments, monomeric units of protegrin, mellitin, alamethicin, and zervamicin, double helices of gramicidin A). A representative structure (usually the one with the highest resolution) was selected for every protein, i.e., we did not consider all mutants, complexes with different ligands, or independent crystal determinations. It was important to conduct all calculations for biologically relevant quaternary structures (biological units) rather than for individual polypeptide chains or non-native oligomers. The quaternary structures were taken from the Protein Quaternary Structure (PQS) database (Henrick and Thornton 1998), except several known non-native complexes such as the antiparallel dimer of rhodopsin. All cofactors, lipids, detergents, ions, and water molecules were excluded from the calculations, except those in photosynthetic complexes where they form a significant part of the protein core (photosynthetic reaction centers, photosystems I and II, and light-harvesting and cytochrome b6f complexes), and in the complex of lipid flippase with lipopolysaccharide (1z2r). In all other complexes, the incorporation of cofactors or crystallized lipids only led to minor changes of the calculated parameters. Some disordered N- or C-terminal segments (5–7 in 1pw4, 62–69 in 1afo, and 1030–1036 in 1oye) and three misplaced N-terminal helices of KvAP channel (1orq) were omitted in final calculations.

All regular secondary structures were first determined by the DSSP algorithm incorporated in QUANTA. Then, all gaps in the middle of TM -helices and -strands were eliminated. All -aneurisms, short fragments of 3₁₀ helices, and -bulges were interpreted as parts of continuous secondary structures. Each regular secondary structure was represented by a vector passing through its starting and ending points. These points were calculated as average positions of several C atoms at the beginning or end of the corresponding -helix or -strand (the number of C atoms for the averaging was chosen as 7, 4, or 1 depending on length of the secondary structure).

Electronic supplemental material
Supplemental materials include PDB codes of all proteins used in this study, hydrophobic thicknesses and tilt angles calculated for different crystal forms of TM proteins, a list of 24 TM proteins used for testing the methods, calculated and ESR-determined water-inaccessible segments of rhodopsin, and a list of flexible side chains whose conformers are adjusted to fit the membrane boundaries.

	Footnotes

 
Supplemental material: see www.proteinscience.org

Reprint requests to: Andrei L. Lomize, College of Pharmacy, University of Michigan, Ann Arbor, MI 48109-1065, USA; e-mail: almz{at}umich.edu; fax: (734) 763-5595.

Article and publication are at http://www.proteinscience.org/cgi/doi/10.1110/ps.062126106.

Abbreviations:TM, transmembrane; OPM, Orientations of Proteins in Membranes (database); PPM, Positioning of Proteins in Membranes, (software); PDB, Protein Data Bank; ATR FTIR spectroscopy, attenuated total reflection Fourier transform infrared spectroscopy; MD, molecular dynamics; EM, electron cryo-microscopy; NEM, N-ethylmaleimide; DM, n-dodecyl--D-maltoside; DHPC, 1,2-dihexanoyl-sn-glycero-3-phosphatidylcholine; RMSD, root-mean-square deviation; G_transf, transfer free energy; D, hydrophobic thickness; , tilt angle; , atomic solvation parameter.

	Acknowledgments

 
We thank Drs. Simon Hubbard, Vladimir Maiorov, and Simon Sherman for provided software, Dr. Kim Henrick for discussion of protein quaternary structure, and Drs. Eugene Krissinel, and Gabor Tusnady for explanations in regard to the SSM server and the PDB_TM database. This work was supported by NIH grant DA003910 (HIM) and an Upjohn Research Award from the College of Pharmacy, University of Michigan (ALL).

	References

TOP Abstract Introduction Results Discussion Conclusion Materials and methods References

 
Abraham M.H., Chadha H.S., Whiting G.S., Mitchell R.C. 1994. Hydrogen bonding. 32. An analysis of water-octanol and water-alkane partitioning and the logP parameter of Seiler J. Pharm. Sci. 83: 1085–1100.[Medline]

Akitake B., Anishkin A., Sukharev S. 2005. The "dashpot" mechanism of stretch-dependent gating in MscS J. Gen. Physiol. 125: 143–154.[Abstract/Free Full Text]

Altenbach C., Marti T., Khorana H.G., Hubbell W.L. 1990. Transmembrane protein structure: spin labeling of bacteriorhodopsin mutants Science 248: 1088–1092.[Abstract/Free Full Text]

Altenbach C., Greenhalgh D.A., Khorana H.G., Hubbell W.L. 1994. A collision gradient method to determine the immersion dept of nitroxides in lipid bilayers: application to spin-labeled mutants of bacteiorhdopsin Proc. Natl. Acad. Sci. 91: 1667–1671.[Abstract/Free Full Text]

Andersen O.S., Koeppe R.E., Roux B. 2005. Gramicidin channels IEEE Trans Nanobioscience 4: 10–20.[CrossRef][Medline]

Andronesi O.C., Pfeifer J.R., Al-Momani L., Ozdirekcan S., Rijkers D.T.S., Angerstein B., Luca S., Koert U., Killian J.A., Baldus M. 2004. Probing membrane protein orientation and structure using fast magic-angle-spinning solid-state NMR J. Biomol. NMR 30: 253–265.[CrossRef][Medline]

Arkin I.T., Rothman M., Ludlam C.F.C., Aimoto S., Engelman D.M., Rothschild K.J., Smith S.O. 1995. Structural model of the phospholamban ion channel complex in phospholipids membranes J. Mol. Biol. 248: 824–834.[CrossRef][Medline]

Ash W.L., Zlomislic M.R., Oloo E.O., Tieleman D.P. 2004. Computer simulations of membrane proteins Biochim. Biophys. Acta 1666: 158–189.[Medline]

Baaden M. and Sansom M.S.P. 2004. OmpT: Molecular dynamics simulations of an outer membrane enzyme Biophys. J. 87: 2942–2953.[Medline]

Barclay P. and Findlay J.B.C. 1984. Labelling of the cytoplasmic domains of ovine rhodopsin with hydrophilic chemical probes Biochem. J. 220: 75–84.[Medline]

Bass R.B., Strop P., Barclay M., Rees D.C. 2002. Crystal structure of Escherichia coli MscS, a voltage-modulated and mechanosensitive channel Science 298: 1582–1587.[Abstract/Free Full Text]

Basyn F., Spies B., Bouffioux O., Thomas A., Brasseur R. 2003. Insertion of X-ray structures of proteins in membranes J. Mol. Graph. Model. 22: 11–21.[CrossRef][Medline]

Bechinger B., Ruysschaert J.M., Goormaghtigh E. 1999. Membrane helix orientation from linear dichroism of infrared attenuated total reflection spectra Biophys. J. 76: 552–563.[Medline]

Benz R., Frohlich O., Lauger P., Montal M. 1975. Electrical capacity of black lipid films and of lipid bilayers made from monolayers Biochim. Biophys. Acta 394: 323–334.[Medline]

Berman H.M., Westbrook J., Feng Z., Gilliland G., Bhat T.N., Weissig H., Shindyalov I.N., Bourne P.E. 2000. The Protein Data Bank Nucleic Acids Res. 28: 235–242.[Abstract/Free Full Text]

Bishop R.E. 2005. The lipid A palmitoyltransferase PagP: molecular mechanisms and role in bacterial pathogenesis Mol. Microbiol. 57: 900–912.[CrossRef][Medline]

Blaurock A.E. and Wilkins M.H.F. 1972. Structure of retinal photoreceptor membranes Nature 236: 313–314.[CrossRef][Medline]

Bond P.J., Faraldo-Gomez J.D., Sansom M.S.P. 2002. OmpA: A pore or not a pore? Simulation and modeling studies Biophys. J. 83: 763–775.[Medline]

Booth P.J., Templer R.H., Meijberg W., Allen S.J., Curran A.R., Lorch M. 2001. In vitro studies of membrane protein folding Crit. Rev. Biochem. Mol. Biol. 36: 501–603.[CrossRef][Medline]

Bowie J.U. 1997. Helix packing in membrane proteins J. Mol. Biol. 272: 780–789.