Loading [MathJax]/jax/output/SVG/jax.js
Research article Topical Sections

Molecular dynamics study of homo-oligomeric ion channels: Structures of the surrounding lipids and dynamics of water movement

  • Molecular dynamics simulations were used to study the structural perturbations of lipids surrounding transmembrane ion channel forming helices/helical bundles and the movement of water within the pores of the ion-channels/bundles. Specifically, helical monomers to hexameric helical bundles embedded in palmitoyl-oleoyl-phosphatidyl-choline (POPC) lipid bilayer were studied. Two amphipathic α-helices with the sequence Ac-(LSLLLSL)3-NH2 (LS2), and Ac-(LSSLLSL)3-NH2 (LS3), which are known to form ion channels, were used. To investigate the surrounding lipid environment, we examined the hydrophobic mismatch, acyl chain order parameter profiles, lipid head-to-tail vector projection on the membrane surface, and the lipid headgroup vector projection. We find that the lipid structure is perturbed within approximately two lipid solvation shells from the protein bundle for each system (~15.0 Å). Beyond two lipid “solvation” shells bulk lipid bilayer properties were observed in all systems. To understand water flow, we enumerated each time a water molecule enters or exited the channel, which allowed us to calculate the number of water crossing events and their rates, and the residence time of water in the channel. We correlate the rate of water crossing with the structural properties of these ion channels and find that the movements of water are predominantly governed by the packing and pore diameter, rather than the topology of each peptide or the pore (hydrophobic or hydrophilic). We show that the crossing events of water fit quantitatively to a stochastic process and that water molecules are traveling diffusively through the pores. These lipid and water findings can be used for understanding the environment within and around ion channels. Furthermore, these findings can benefit various research areas such as rational design of novel therapeutics, in which the drug interacts with membranes and transmembrane proteins to enhance the efficacy or reduce off-target effects.

    Citation: Thuy Hien Nguyen, Catherine C. Moore, Preston B. Moore, Zhiwei Liu. Molecular dynamics study of homo-oligomeric ion channels: Structures of the surrounding lipids and dynamics of water movement[J]. AIMS Biophysics, 2018, 5(1): 50-76. doi: 10.3934/biophy.2018.1.50

    Related Papers:

    [1] Daniela Meleleo, Cesare Sblano . Influence of cholesterol on human calcitonin channel formation. Possible role of sterol as molecular chaperone. AIMS Biophysics, 2019, 6(1): 23-38. doi: 10.3934/biophy.2019.1.23
    [2] Lars H. Wegner . Cotransport of water and solutes in plant membranes: The molecular basis, and physiological functions. AIMS Biophysics, 2017, 4(2): 192-209. doi: 10.3934/biophy.2017.2.192
    [3] Mathieu F. M. Cellier . Evolutionary analysis of Slc11 mechanism of proton-coupled metal-ion transmembrane import. AIMS Biophysics, 2016, 3(2): 286-318. doi: 10.3934/biophy.2016.2.286
    [4] István P. Sugár . A generalization of the Shell Theorem.
    Electric potential of charged spheres and charged vesicles surrounded by electrolyte. AIMS Biophysics, 2020, 7(2): 76-89. doi: 10.3934/biophy.2020007
    [5] Wenyuan Xie, Jason Wei Jun Low, Arunmozhiarasi Armugam, Kandiah Jeyaseelan, Yen Wah Tong . Regulation of Aquaporin Z osmotic permeability in ABA tri-block copolymer. AIMS Biophysics, 2015, 2(3): 381-397. doi: 10.3934/biophy.2015.3.381
    [6] István P. Sugár, Parkson Lee-Gau Chong . Monte Carlo simulations of the distributions of intra- and extra-vesicular ions and membrane associated charges in hybrid liposomes composed of negatively charged tetraether and zwitterionic diester phospholipids. AIMS Biophysics, 2017, 4(2): 316-336. doi: 10.3934/biophy.2017.2.316
    [7] Jacques Fantini, Francisco J. Barrantes . How membrane lipids control the 3D structure and function of receptors. AIMS Biophysics, 2018, 5(1): 22-35. doi: 10.3934/biophy.2018.1.22
    [8] Sudarat Tharad, Chontida Tangsongcharoen, Panadda Boonserm, José L. Toca-Herrera, Kanokporn Srisucharitpanit . Local conformations affect the histidine tag-Ni2+ binding affinity of BinA and BinB proteins. AIMS Biophysics, 2020, 7(3): 133-143. doi: 10.3934/biophy.2020011
    [9] Ken Takahashi, Yusuke Matsuda, Keiji Naruse . Mechanosensitive ion channels. AIMS Biophysics, 2016, 3(1): 63-74. doi: 10.3934/biophy.2016.1.63
    [10] Chia-Wen Wang, Meng-Han Lin, Wolfgang B. Fischer . Cholesterol affected dynamics of lipids in tailor-made vesicles by ArcVes software during multi micro second coarse grained molecular dynamics simulations. AIMS Biophysics, 2023, 10(4): 482-502. doi: 10.3934/biophy.2023027
  • Molecular dynamics simulations were used to study the structural perturbations of lipids surrounding transmembrane ion channel forming helices/helical bundles and the movement of water within the pores of the ion-channels/bundles. Specifically, helical monomers to hexameric helical bundles embedded in palmitoyl-oleoyl-phosphatidyl-choline (POPC) lipid bilayer were studied. Two amphipathic α-helices with the sequence Ac-(LSLLLSL)3-NH2 (LS2), and Ac-(LSSLLSL)3-NH2 (LS3), which are known to form ion channels, were used. To investigate the surrounding lipid environment, we examined the hydrophobic mismatch, acyl chain order parameter profiles, lipid head-to-tail vector projection on the membrane surface, and the lipid headgroup vector projection. We find that the lipid structure is perturbed within approximately two lipid solvation shells from the protein bundle for each system (~15.0 Å). Beyond two lipid “solvation” shells bulk lipid bilayer properties were observed in all systems. To understand water flow, we enumerated each time a water molecule enters or exited the channel, which allowed us to calculate the number of water crossing events and their rates, and the residence time of water in the channel. We correlate the rate of water crossing with the structural properties of these ion channels and find that the movements of water are predominantly governed by the packing and pore diameter, rather than the topology of each peptide or the pore (hydrophobic or hydrophilic). We show that the crossing events of water fit quantitatively to a stochastic process and that water molecules are traveling diffusively through the pores. These lipid and water findings can be used for understanding the environment within and around ion channels. Furthermore, these findings can benefit various research areas such as rational design of novel therapeutics, in which the drug interacts with membranes and transmembrane proteins to enhance the efficacy or reduce off-target effects.


    1. Introduction

    Biological membranes contain various types of lipid, proteins, and other molecules whose functions are essential for life. To understand biological membranes, it is critical to understand the interactions of their main components, lipids, proteins, and other molecules. It is well known that the structure and function of the membrane proteins, which impart many of the functions to biological membranes, are highly sensitive to the intermolecular forces that act within the lipid bilayer, such as polar and nonpolar interactions [1]. It is also well known that different bilayer environments play a critical role in altering protein structure and stability. A particular class of membrane structures are ion channels which are essential in maintaining normal cell functions such as sustaining osmotic equilibrium, facilitating bioenergetics, and providing the means for transmitting environmental signals [2,3]. Ion channels are targets for many therapeutics, and many classes of ion channels have been implicated in either genetic or acute diseases [4,5]. By offering insight into the local molecular environment around and within ion channels, possible therapeutics such as new antivirals or treatments to channel-related diseases could be developed.

    Ion channels are often multimeric proteins containing many transmembrane (TM) segments. To understand general ion-channel properties, a minimalist approach was adapted to comprehend the function of ion channels which uses synthetic peptides that consist of repetitive sequences with features and topology found to be important in more complex and natural channel-forming proteins [6]. These canonical synthetic peptides are amenable to study, and the conclusions of such studies should be applicable to the entire class of homo-oligomeric ion-channels including naturally occurring channels. These peptides have been successful in studying ion-channel structure, assembly and gating [6,7]. The two synthetic peptides we use are Ac-(LSLLLSL)3-NH2 (LS2), and Ac-(LSSLLSL)3-NH2 (LS3), where leucine (L) is nonpolar with high hydrophobicity and helix-forming propensity, and serine (S) is polar and represents the hydrophilic side of the helix. In addition to the amphipathic nature of the peptide found in ion-channels, the length of the peptide (21 residues ∼30 Å) is designed to be comparable to the thickness of the hydrocarbon portions of the lipid bilayer minimizing the hydrophobic mismatch [8,9]. While these peptides have a simple sequence, they are known to form ion channels in membranes and are voltage gated similar to biologically active ion channels [10,11]. The sequence pattern of L/S residues forming LS2 and LS3 peptides are topologically similar, each with a hydrophilic face which would form the lumen of the pore and hydrophobic face interacting with lipids or other hydrophobic helices. This topology can be seen in many naturally occurring peptide sequences such as those found in virus [12], e.g. M2 influenza virus [13,14,15,16,17] and VPU [18], fungus e.g. Alamethicin [19] and mammalian e.g. channel-lining segments of nicotinic AchR [20].

    Previous work has shown that the LS3 peptides forms monovalent channels that is cation-selective and prefers to be a hexameric channel, placing LS3 in the group of channels with midsized pores which include ion-selective nicotinic receptors and anion-selection GABA and glycine receptors [21]. Unlike LS3, LS2 peptides prefer tetrameric bundles and are proton selective [10,22,23], behaving similar to the influenza A viral M2 protein [16,17,24]. Understanding the movement of water traveling within these simple molecular assemblies of the LS2 and LS3 channels can garner information on pore size preferences that can link these simple peptides with naturally occurring peptides that have similar channel behaviors.

    The mobile lipids next to the more rigid membrane proteins can adjust to a changing environment without compromising the integrity of the membrane. Therefore the interplay between membrane proteins and lipids are crucial. Further, lipids are essential in protein assembly or oligomerization of multi-subunit complexes, and they help maintain channel stability, control protein insertion and folding processes [1,25,26]. Numerous studies have been conducted on the biochemical and biophysical importance of protein-lipid interactions for assembly, stability, and function of membrane proteins, however, atomistic detail of protein-lipid interactions are required to understand the effects of lipid structure on the aggregation of transmembrane proteins [1]. Despite considerable efforts, it is experimentally difficult to map details of protein-lipid interactions since the lipid bilayer environment is a complex two-dimensional liquid crystalline system [27]. Investigators have used Molecular Dynamics (MD) simulations to shed light on biomembrane phenomenon such as helix aggregation, local motions and interactions with lipid bilayers [27]. Here, we aim to build upon our previous study of ion-channel protein-protein interactions [9] to offer insight into the interactions of TM proteins with lipid molecules and water flow through the ion-channels.


    2. Materials and methods

    MD simulations were carried out on the LS2 and LS3 peptides that form ion channels. The details of the methods are described in our previous paper on peptide-peptide interactions [9], for completeness we give the details of the simulations here. Models ranging from one to six helix bundles were generated using an idealized structure with the TM bundle parallel to the bilayer normal and the serine residues of each α-helices lining the pore, forming a N-fold symmetry where N is the number of helices in the bundle. Each bundle were embedded in a pre-equilibrated lipid bilayer consisting of 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) (Figure 1) and placed between two lamellae of waters. A total of 200 lipid molecules (100 per leaflet) and 11,228 water molecules were used for each system, for a total of 60,484 atoms for the environment, with 360 atoms per helix for LS2 and 334 atoms per helix for LS3.

    Simulations were carried out using the NAMD package [28] with the AMBER03 force field [29] and TIP3P waters [30]. Periodic boundary conditions were applied in all three spatial directions. Particle Mesh Ewald (PME) method was used for long-ranged electrostatic forces with a grid of 90 × 90 × 90 [31]. A time step of 1.0 fs was used. An initial equilibration was done using the NVT ensemble for 5 ns under a constant temperature of 310 K where the peptide backbone alpha-carbons were fixed to remove all “bad” interactions. After these initial simulations, an additional 10 ns was performed without any constraints under the NVT ensemble. The calculations were then continued for at least 75 ns using an NPγT ensemble where a surface tension (γ) of 45 mN/m was applied [32,33,34], which reproduced a surface area of 65 ± 2 Å2/headgroup [35]. This is in accord with experimental value of 68.3 ± 1.5 Å2/headgroup [36]. The thickness of the bilayer was found to be 40 ± 1 Å, relatively in good agreement with the experimental value of 37 Å [36]. A constant pressure of 1 atm was maintained using Langevin piston with a period of 200 fs and a decay time of 50 fs. The temperature was maintained at 310 K using Langevin damping with a damping coefficient of 5 ps−1. The simulation unit cell fluctuates separately in each dimension and had an average size of roughly 84 Å × 76 Å × 94 Å, which is large enough such that periodic images of the peptides do not interact.

    Figure 1. Top: TM ion channel forming helices, where (A) is LS2 and (B) is LS3. Each helix is depicted in a semitransparent representation where cyan represents the backbone, red represents polar residues, and green represents nonpolar residues. Each panel depicts the same protein TM sequence rotated 180° with the hydrophilic face on the right and the hydrophobic face on the left. Bottom: 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) with labeled atom names, where the palmitoyl tail is from atoms C32-C316, and the oleoyl tail is from atoms C22-C218.


    3. Results and discussion

    The topology of LS3 and LS2 peptides differs where LS3 consists of more polar serine residues than LS2 (Figure 1A, B). The different structures resulting from packing and tilting of the oligomers results in commensurate lipid structural and dynamics changes near the oligomers. We quantify the collective bilayer response hydrophobic mismatch and meniscus formation, and give more detail of the molecular structural changes with lipid order parameters and vector projections.

    As well as the response of the lipid, these oligomers form ion-channels, which allow for water to transverse the bilayer. Due to the difference in the number of polar residues, we initially hypothesize that the structure of the pore of the LS3 ion channels would allow for more water to flow through since more polar interactions with water can be formed. However, our results suggest that the polar topography of the helix does not dictate the flow of water, but rather the structural arrangement of the bundle control the flow of water. We have quantified the structure packing of the bundles in a previous publication [9]. Here we highlight and expand on structural features related to water flow. We find that only bundles with more than three helices exhibit water flow and thus form ion channels. To examine the movement of water, we enumerate each water-crossing event by following the water movement through the hydrophobic region of the bilayer.


    3.1. Oligomer structure and helix vector projection

    The structural and orientation of the helices within the oligomers have been discussed in a previous paper [9], here we highlight structure of the bundles with N > 4 and the pore lumen through which water transverses. The average vector projections (representing the hydrophobic to hydrophilic direction) for each helix are reported in Figure 2. For each helix, this is defined as the XY projection of the vector from the center-of-mass (COM) of the backbone atoms of the center residue (11th) to the COM of the third or fifth serine side chain atoms for LS2 or LS3, respectively. This produces a vector pointing toward the hydrophilic face of the peptide, which highlights the structural arrangement of each helix and the alignment of the serine residues with respect to the lumen of the pore (Figure 2).

    Figure 2. LS2 (blue) and LS3 (red) average protein structure (A) and vector projections (B): (1) LS2 Tetramer; (2) LS2 Pentamer; (3) LS2 Hexamer; (4) LS3 Tetramer; (5) LS3 Pentamer; (6) LS3 Hexamer.

    Base on Figure 2, it is evident that a pore is formed in all bundles with N > 4 (N = 4-6 peptides) where the hexameric bundle has the largest pore, and the tetramer bundle has the smallest pore. The LS2 tetramer and hexamer bundle forms a stable pore, while the pentamer bundle show signs of instability where helix 3 (H3) is interacting strongly with helix 2 (H2). In comparison, all of the LS3 bundles (N = 4-6) have a stable pore lining arrangement with the hexamer conformation having the largest pore. We find that all adjacent helices within each bundle have at least 1 serine-serine contact with LS3 having more contacts than LS2 as expected.


    3.2. Perturbations of lipid structure

    We analyzed the lipid-protein interactions of embedded bundles comprised of 1-6 peptides. Specifically, we investigate the hydrophobic mismatch between the protein and the lipid bilayer, and its effects on the lipid structures. We analyzed the lipid structure using the head-to-tail orientations, acyl-chain order parameters, and the headgroup orientations. We expand on each of these analyses and discuss the results in turn below. Anticipating our results, we find distinctive lipid perturbation patterns surrounding bundles containing 1-3 peptides (lower order bundles) compared to bundles with 4-6 peptides (higher order bundles). In addition, we determined that the membrane proteins only perturb the lipid structure for 2 lipid solvation shells.


    3.2.1. Hydrophobic mismatch and meniscus formation

    Hydrophobic mismatch occurs when the hydrophobic portion of the proteins do not match the acyl core of the membrane lipids, thus exposing the hydrophobic surfaces to an hydrophilic environment which is either the lipid head groups or the aqueous solution. In order to avoid such exposure and to minimize the free energy of the system, proteins and/or lipids typically adopt different conformations to reduce or eliminate the mismatch. Proteins may tilt to change their effective hydrophobic length, self-associate/aggregate, or adapt their structure. Lipids can react by stretching, compressing, or distorting their acyl chains to increase or decrease lipid thickness to match the proteins. Numerous studies at the molecular level have been conducted to elucidate the effects caused by the protein-lipid hydrophobic mismatch. Results indicated that the extension of lipid perturbation depends on several criteria such as protein size, the degree of the hydrophobic mismatch, the curvature of the protein surface in contact with the lipid hydrocarbon chains, and the temperature setting during the simulation [37,38,39,40]. Most models assume that the lipid's fatty acyl chains in the vicinity of a membrane protein adjust their conformation to match the hydrophobic thickness of the protein with minimal local perturbation.

    To investigate the hydrophobic mismatch at the molecular level, we plotted the height of the lipid head group to the center of lipid as a function of the radial distance from the COM of the TM protein (Figure 3). The height of the lipid head group is defined as the distance from the phosphorus atom to the COM of the bilayer projected along the lipid bilayer normal (Z-axis), averaged over all lipid molecules in the two leaflets and over the configurations sampled over time. Helix tilt angles are measured as the angle between the principal eigenvector of the moment of inertia (MOI) tensor of each bundle (Table 1) and the normal of the bilayer-water interface [41].

    Table 1. Total simulation time and the average tilt angle of the bundle with respect to the bilayer's norm, calculated by using the moment of inertia tensor and averaged over the last 50 ns of trajectory.
    Systems LS2 Simulation Time (ns) LS3 Simulation Time (ns) LS2 Tilt Angle (degrees) LS3 Tilt Angles (degrees)
    Monomer 75 75 22.3 ± 2.0° 5.89 ± 2.9°
    Dimer 80 90 28.4 ± 2.1° 5.06 ± 2.0°
    Trimer 135 80 21.7 ± 2.0° 10.1 ± 2.7°
    Tetramer 80 75 3.5 ± 2.9° 2.3 ± 1.3°
    Pentamer 80 120 29.7 ± 1.8° 18.9 ± 2.1°
    Hexamer 80 80 24.0 ± 5.8° 19.0 ± 4.2°
     | Show Table
    DownLoad: CSV
    Figure 3. Lipid head-group height (Z-project of the distance from the phosphorus atom and the COM of the bilayer) averaged over time and two leaflets as a function of its radial distance from the monomer/bundle. The cold colors (aqua, blue, purple) represent LS2 systems and hot colors (yellow, red, orange) represent LS3 systems. (A) monomer to trimer systems; (B) tetramer to hexamer systems.

    Simulation of a pure POPC bilayer system yielded an average lipid head group height of 20.3 ± 2.1 Å. For the protein/lipid systems, Figure 3 shows that the lipid head-group heights range from a minimum of ∼17.5 Å to 19.0 Å in systems with three or less peptides (lower order bundles), and from a minimum of ∼16.0 to 19.7 Å in systems with four or more peptides (higher order bundles). The general trend found was that as the distance from the bundle/monomer increases, the membrane becomes unperturbed as expected, but in the vicinity closest to the bundles/monomers (less than 10 Å for lower order, and less than 15 Å for higher order), the bilayer thins slightly. At areas closest to the bundles/monomers, the lipid chains expand out along the bilayer interface to provide a thinner bilayer with a small decrease in bilayer thickness. This is due to the hydrophobic thickness of the bilayer being somewhat greater than that of the protein helix length therefore it decreases to match with the protein and minimize the free energy of the system. However, since the meniscus of the bilayer with embedded peptides is within the standard deviation of 2.1 Å of the height of a pure bilayer (random fluctuations of the pure bilayer are more than the mismatch), we conclude that the hydrophobic mismatch effect is not significant in the simulated systems.

    In terms of protein conformation, we have observed significant tilting of peptides in most systems (Table 1). The tilt angles of the bundle are found to be greater in the LS2 systems than that in the LS3 systems (Table 1). An extreme example can be observed in the monomer system where the LS2 peptide has a tilt angle of 22° while the LS3 peptide remained fairly vertical with a tilt angle of 6°. However, there is no clear correlation between the protein tilt angles (Table 1) and the lipid height meniscus (Figure 3). In addition, the length of the LS2 and LS3 peptides (∼30 Å) is slightly shorter than the height of the hydrophobic core of lipid bilayer. Tilting of the peptides would not lead to matching hydrophobic lengths but the opposite and more of a mismatch.

    We reason that the tilting of the peptides is not caused by meniscus formation in response to hydrophobic mismatch where the perturbation on protein orientation and structure are not caused by membrane curvature, but rather direct protein-lipid interactions. To elaborate, the hydrophobic interaction of the membrane core with the leucine residues is competing with the hydrophilic interactions between the serine residues and lipid head groups. Thus, if the peptide contained only nonpolar leucine residues (LLLLLLL)3, then the peptide would be located at the center of the bilayer in a parallel orientation (tilt angle of 90°) leaving it in a parallel orientation that is not transmembrane. On the other hand, a peptide containing only polar serine residues (SSSSSSS)3 will be expelled from the bilayer's core and reside in the aqueous phase. Hence, combinations of polar and nonpolar residues are required to maintain a transmembrane conformation. Experimentally it has been shown that the LS3 monomer is not transmembrane and prefers an orientation that is parallel to the membrane interface [42]. However, within a simulation time of 75 ns, the LS3 monomer is in a meta-stable state where it remains transmembrane, and we see the monomer rocking back and forth in this local minimum. Since LS3 has more polar residues than LS2, LS3 interacts more strongly with the lipid headgroups to keep its vertical orientation, leading to a smaller tilt angle. It has been estimated that a single polar amino acid such as serine, can stabilize the transmembrane peptides by 1-2 kcal/mol [43].


    3.2.2. Hydrocarbon lipid tail order parameter

    Lipid acyl chain order parameter is a direct measurement of the hydrocarbons chain orientation and mobility. The deuterium order parameter can be obtained through 2H NMR experiments from the quadrupole splitting [44,45,46]. We calculated the lipid deuterium order parameters (Sbilayer) of POPC acyl tails (Figure 1) using:

    Sbilayer=123cos2θ1 (1)

    where θ is the angle between the molecular axis (carbon-hydrogen bond vector of methylene segment) with respect to the bilayer normal. We calculate acyl chain order parameters on both the palmitoyl and oleoyl tails (Figure 1). Lipid molecules within and beyond two solvation shells (10.0 Å for lower order bundles and 15.0 Å for higher order bundles) with respect to the protein bundle were separated and the order parameters were calculated separately. Order parameters of lipids more than two solvation shells away from the bundles were found to be similar to that of a pure bilayer. We thus conclude that lipid perturbation by the transmembrane protein system does not extend beyond two solvation shells. On the other hand, Figure 4 shows the lipid deuterium order parameters measured for lipids that are within two solvation shells of the bundles are affected. Each system is affected to a different extent by the presence of the protein bundles.

    Figure 4. Lipid acyl chain deuterium order parameter of lower order bundles (A & C) and higher order bundles (B & D) for both the palmitoyl (A & B) and oleoyl (C & D) tails within about two solvation shells (10 Å for 1-3 bundles and 15 Å for 4-5 bundles) are shown. The cold colors (aqua, blue, purple) represent LS2 systems, hot colors (yellow, red, orange) represent LS3 systems and the black line represent the order parameter of a pure POPC bilayer.

    Typically, a higher order parameter corresponds to a more ordered system. The saturated lipid tail of POPC displayed a general trend of decreasing lipid order towards the core of the bilayer (Figure 4A, B). The unsaturated lipid chain is not monotonic, where there are two minimum ranges indicating that the unsaturated acyl chains are more disordered. The first range is between C22-C25 atoms of the acyl chain in the vicinity of the glycerol group. Second, the acyl tail becomes more disorder near the point of unsaturation between C29-C210 (Figure 4C, D). Experimental deuterium NMR studies of lipid-only system yielded a similar trend for both saturated and unsaturated acyl chains [46,47]. The order parameters of the fatty lipid tails are significantly affected by the presence of the lower order bundles (N = 1-3). The lipids surrounding the LS3 peptides are more ordered than that of the LS2 peptides, which likely results from the LS2 peptides tilting allowing more space resulting in increased mobility and decrease order of the lipid tails. Clearly observed in the palmitoyl tail (Figure 4A), the hydrocarbons order parameters in a pure bilayer lie in between LS2 and LS3 systems. In the oleoyl tail, the same trend is observed from C211-C218 atoms, towards the bilayer core. However, due to the nearby glycerol group and the presence of the double bond, the oleoyl tail order parameters from C22-C210 are affected by the protein bundles somewhat differently. All systems are more disordered as compared to a pure bilayer system. As shown in Figure 4, in higher order bundles (N = 4-6), the lipid order parameters were found to be fairly close to the pure POPC bilayer and with each other for both the palmitoyl and oleoyl tails. The fact that the presence of these higher order bundles has little effect on the order parameters of lipid tails is probably due to less accessibility of the polar residues to the fatty acyl chains in these higher order bundles.


    3.2.3. Lipid tail vector field and vector length

    To further quantify the lipid tail ordering we use a lipid tail vector projection (LTVP) and average it over the simulation trajectory. A LTVP is defined as the vector from the central glycerol carbon (C2 in Figure 1) to the respective terminal carbon (C218 or C316 in Figure 1) projected onto the XY plane. In an equilibrated pure lipid system or an unperturbed lipid such as an area away from the protein bundles, the average LTVP over time would be expected to vanish indicating no preferred orientation of tail directions. Therefore, large deviation from zero indicates a preferred direction of the lipid tails.

    A graphical representation of the radial and angular distribution of the LTVPs with respect to the COM of the protein bundle/monomer is shown in Figure 5. The XY plane of the membrane were divided based on distance and angle from the COM of the bundle with a bin size of 5.0 Å for distance and 40.0° for angle. The colors of the arrows represent the density of lipids in the bin area, where purple has the lowest density, and yellow with the highest density. The origin is the COM of the bundles/monomers, the z-axis is perpendicular to the bilayer plane, and the x-axis is from the COM of the bundle to the COM of helix 1 or in the monomer case from the COM of the helix to the COM of all leucine residues. Figure 6 quantifies and Figure 5 displays the LTVP lengths as a function of the distance from the COM of the bundle/monomer averaged over the angular space and two leaflets (upper and lower).

    The acyl chains in the vicinity of the protein bundles generally have longer LTVPs indicating lipid perturbation by the protein bundles. As a result, these lipids show preferred orientation and/or tilting of their tails from a vertical orientation along the bilayer normal to a more extended orientation along the membrane surface. Typically, the lipid tail extends away from the bundle. This is especially true on the polar residue side of the peptide as evidenced in the LS3 monomer and dimer cases (Figure 5D, E). Occasionally, some lipid tails reach towards the hydrophobic face of the bundle. This is mostly due to the tilting of the protein bundles leaves voids around the bundles therefore lipids extend their tails to try to fill the void. Figure 6 shows that the vector lengths die down rather quickly for lipids that are more than 10.0 Å (lower order bundle) or 15.0 Å (higher order bundle) away from the center of the protein bundles. The fact indicates that lipids are only perturbed locally within 10.0-15.0 Å of the homo-oligomeric bundle systems which is consistent with the order parameter results.

    Figure 5. LTVPs around: (A) LS2 monomer; (B) LS2 dimer; (C) LS2 Tetramer; (D) LS3 monomer; (E) LS3 dimer; (F) LS3 Tetramer conformation. The arrow represents the vectors from the glycerol to terminal carbon groups projected on the XY plane, while the listed colors represent the density of the lipids in the order of least dense to most dense: Purple (lowest) > Blue > Green > Red > Orange > Yellow (highest).

    Second, we observe that acyl chains of the lipid tails tilt more (longer vectors) in the vicinity of the protein in response to the tilting of the lower ordered protein bundles. This is evident as we compare the protein tilt angles (Table 1) versus the lipid tail orientation (Figures 5 and 6A) between the LS2 and LS3 lower order bundle systems. The LTVPs of the LS2 monomer and dimer systems are significantly longer than those of the LS3 monomer and dimer systems, respectively. For example for the first solvation shell (∼5 Å to the COM of protein), while the LS2 dimer has an averaged LTVP length of 9.0 Å, the LS3 dimer only has an average vector length of 4.2 Å. Therefore the amphipathic peptides, when they act alone or in lower order bundle, perturb the lipid structure to a higher extent when they are more tilted. This is probably due to the fact that the higher degree of tilting would result in more lipids getting in contact with the hydrophilic faces of these lower order bundles. As the size of the bundles grow and the hydrophilic faces of the peptides get more shielded, the effect of tilting on lipid orientation decreases. For example, the differences in LTVP lengths between the LS2 and LS3 trimer systems are much smaller than that of the monomer and dimer systems.

    Figure 6. LTVP length averaged over time and both leaflets with respect to the lipid's radial distance from the center of the monomer/bundle. The cold colors (aqua, blue, purple) represent LS2 systems and hot colors (yellow, red, orange) represent LS3 systems. (A) monomer to trimer systems; (B) tetramer to hexamer systems.

    Third, in higher order bundles (N = 4-6 peptides), the radial and angular distributions of the LTVPs are not significantly influenced by the tilting of the peptides since the polar residues are shielded from the lipid core. In the tetramer case, the vector field of the LS2 tetramer (Figure 5C) shows that the LTVPs are shorter in length, less directional and are more evenly distributed around the peptide bundle than that of the LS3 tetramer. Our analysis of protein structure and dynamics [9] indicates that LS3 tetramer and pentamer bundles have tighter packing arrangements than the corresponding LS2 bundles due to more hydrophilic contacts within the pore. As a result, the dynamics of the peptides in the LS3 bundles are slower than that of the LS2 systems. This will likely influence the dynamics of the lipids nearby so that lipid tails close to the LS3 tetramer and pentamer are hindered. Hence, the packing arrangements of higher order bundles contribute to the perturbation of the lipid head-to-tail displacement.


    3.2.4. Lipid head group vector field

    Similar to the LTVPs that quantify the tail orientations, we have also analyzed the orientation of the head groups. Specifically, the vectors defined as from the COM of the choline group (approximately N in Figure 1) to the COM of the phosphate group (approximately P in Figure 1) is calculated and projected onto the membrane plane (XY). The graphical representation of the angular and radial distribution of this vector projection (LHVP as in Lipid Head Vector Projection) is generated using the same procedure as discussed above for the LTVP. Specifically, the two dimensional (XY) space around the LS2 and LS3 bundles/monomers were divided into bins with a bin size of 5.0 Å for distance and 40.0° for angles. The average LHVP for each bin over time and space are calculated. Figure 7 shows examples of the LHVP vector field around the protein bundles for lower (monomer and dimer) and higher (tetramer) order bundles for both LS2 and LS3 systems. Further, the radial distributions of the LHVP lengths as a function of the radial distance from the COM of the bundles/monomers are plotted in Figure 8.

    As shown in Figure 8, the perturbation of lipid by the protein bundles in terms of lipid head vector orientation is also localized. However, the perturbation extends slightly further out from the center of the bundles. Specifically, the LHVP lengths level off around 15 Å for lower order bundles and 20 Å for higher order bundles. In both cases, the perturbation in terms of head vector orientation extends out 5 Å more than that of the tail vector orientation. In addition, unlike the LTVPs that are mostly affected by specific protein-lipid interactions, the LHVPs are likely to be affected by the electrical field generated by the hydroxyl groups of the serine residues.

    Figure 7. LHVPs around: (A) LS2 monomer; (B) LS2 dimer; (C) LS2 tetramer; (D) LS3 monomer; (E) LS3 dimer; (F) LS3 tetramer conformation. The peptides are colored green for the leucine residues and red for the serine residues. The arrows represent the vectors from the head group choline to phosphate projected onto the XY plane, where the listed colors represent the density of the lipids from least dense to most dense: Purple (lowest) > Blue > Green > Red > Orange > Yellow (highest).

    Figure 8. LHVP length averaged over time and both leaflets with respect to the distance from the center of the monomer/bundle. The cold colors (aqua, blue, purple) represent LS2 systems and hot colors (yellow, red, orange) represent LS3 systems. (A) monomer to trimer systems; (B) tetramer to hexamer systems.

    For example, in the LS3 monomer system, there is a group of LHVPs on the opposite side of the hydrophilic face and pointing away from the peptide. These choline (positive) to phosphate (negative) vector projections align well with the dipole vectors (from the more positive H to the more negative O atoms) of the hydroxyl groups of the LS3 peptide, which is almost vertically sitting in the membrane. In the LS2 monomer system, the hydroxyl groups are more randomly distributed due to the tilting of the peptide therefore the LHVPs are more randomly distributed. The behaviors of the LHVPs in the dimer systems, especially the LS3 dimer, also support the idea that electrical field is the main source for lipid perturbation in terms of head group orientation. As shown in Figure 7B, E, the LHVPs in the dimer systems are generally short and less directional. This is in agreement with the fact that the dipole moment in the dimer systems are likely small due to the cancelation of charges between the hydroxyl groups from the two peptides when they directly facing each other and forming hydrogen bonds. The effect is more evident in the LS3 case since the two peptides are sitting in the membrane almost vertically with their hydrophilic faces aligning to give the maximum serine-serine contacts [9].

    In the LS3 trimer and all higher order bundle systems, the general pattern of the LHVPs are that they are pointing away from the center of the bundle. The distribution of the LHVPs seems to correlate well to the symmetry of the bundle and this again supports the conclusion that the orientation of the head group vectors is mostly affected by the electrical field. For example, the LS2 tetramer has the highest symmetry in terms of equal distance from each helix to the center of the bundle [9], similar tilt angle for individual helix and similar interhelical distance for each pair, etc. Therefore, the LHVP field of the LS2 tetramer system (Figure 7C) shows a very interesting symmetric pattern revolve around the center of the tetramer bundle, which can be explained by an electrical field caused by the symmetric arrangement of the helical bundle.


    3.3. Water dynamics


    3.3.1. Pore radius profile

    The pore radius profiles for all systems were determined using the HOLE program [48,49]. HOLE measures the maximum pore diameter of ion channels by using a Monte Carlo simulated annealing procedure. Traveling through the z-axis of the center of the pore, HOLE determines the largest sphere that can be accommodated through the channel without overlapping with the van der Waals surface of any atoms of the peptides. To compare the change in pore size, the pore radius profile was calculated during the production run at 25 ns and at the end of the run at 100 ns for all systems (Figure 9).

    Hexamer: In good agreement with the vector projections (Figure 2B), the hexameric bundle has the largest pore radius. A steady stream of water can be observed within the channel in contrast to the other conformations. The larger pore radius found within the hexameric bundle suggests that ions can pass through the pore. In particular, the interaction distance for two TIP3P water molecule is ∼3.15 Å [50], which would allow ions to pass through the pore while remaining partially hydrated throughout.

    Pentamer: The average pore radius during the last 70 ns of simulation time for LS2 was 1.94 ± 0.55 Å and 1.53 ± 0.61 Å for LS3. Both LS2 and LS3 bundles have a narrow region (pore minimum) average of 1.02 ± 0.4 Å for LS2 and 0.427 ± 0.2 Å for LS3, suggesting a “closed” channel. The water positions on the pore radius profile (Figure 9B) at 25 ns versus 100 ns give a picture of the dramatic decline of water as the pentamer relaxed from its initial to final configuration. Compared to the LS2 pentameric bundle, which has 40 water molecules sitting within the bundle, LS3 has significantly less water molecules (20) sitting within the LS3 bundle. While it is evident that the LS3 bundle formed a “closed” channel, it is not so clear for the LS2 bundle. As discussed later, few water molecules was observed crossing through the LS2 channel. This is because the channel is dynamic and the motions of the channel allow for water to move into cavities formed by the channels and eventually cross the channel.

    Figure 9. Pore radius profiles and positions of water molecules along the z-axis (bilayers' norm) at 25 ns and 100 ns. (A) Tetramer; (B) Pentamer; (C) Hexamer. Red for LS2 systems and green for LS3 systems.

    Tetramer: It was observed that the arrangements of the peptides in tetrameric bundles play a pivotal on allowing water to cross. As shown from the pore radius profile (Figure 9A), LS2 tetramer remained in an open conformation at 100 ns and has an average pore radius of 2.3 ± 0.45 Å. Since a contact distance of ∼3.15 Å is needed for waters to interact with one another, this suggest that the LS2 tetrameric bundle contain a small water column, which would hinder the movement of ions to cross through the pore. This is in congruence with studies done by Sansom et al. [23], where the average pore radius profile was 2.2 Å. Due to the restriction of water movement through the small pore, an existing theory suggests that hydrogen's are pass along from one or more water “wires” formed through the interior of the pore in a Grottüs manner [51] which explains why LS2 is a proton selective channel. The LS3 tetramer remained in a closed conformation from its equilibrated structure throughout the simulation with a pore radius less than 1.0 Å (0.78 ± 0.34 Å). In addition, no water was observed within the LS3 tetrameric bundle whereas multiple water molecules were observed within and crossing through the LS2 tetrameric bundle.

    The volume of each bundle was calculated based on the average pore radius and a pore height of ∼30 Å [9]. The number of waters located through the channel was estimated and referred as the ideal number of waters that can occupy the volume of each bundle (e.g. density of water time the volume of the pore). The number of ideal waters was compared to the actual number of water sitting through the lining of the pore at the final structure at 100 ns (Table 2). The numbers of water seen within the pores are in excellent agreement (within ∼10%) to the ideal number of water that can occupy the pore considering the simplicity of the ideal model. Differences in pore profile, protein dynamics, as well as water structure and dynamics seem to cancel out to give a number very close to the ideal model.

    Table 2. The average pore radius profile of each bundle during the last 70 ns with standard deviation, ideal number of water that can occupy the bundle, and the actual number of water molecules inside bundles.
    # of Peptide Systems Avg. Pore Radius (Å) Ideal # of waters Observed # of waters
    4 LS2 2.30 ± 0.45 50 45
    LS3 0.78 ± 0.34 6 4
    5 LS2 1.94 ± 0.55 36 40
    LS3 1.53 ± 0.61 22 20
    6 LS2 4.72 ± 0.56 210 225
    LS3 3.88 ± 0.39 142 143
     | Show Table
    DownLoad: CSV

    3.3.2. Water crossing events

    It is crucial to understand water movement through an ion channel since water provides a pathway for ions to travel, which in turn affects ion selectivity and permeation rates. The flow of water going through the bundles was analyzed by tracking water diffusion from one side of the membrane to the other from −15 to 15 Å on the z-axis. An example of such diffusion is shown in Figure 10 for the LS2 tetrameric bundle in comparison to the LS2 pentameric bundle. Summarized in Table 3 are the numbers of crossing events and the direction of water flow. The flux of the channel is roughly zero where the numbers of waters traveling in a direction from top to bottom (downward direction) and from bottom to top (upward direction) is somewhat equal to one another. This zero flux is expected as there is no driving force in any preferred direction regardless the fact that the channel is structurally asymmetric. Specifically, in order to get asymmetric flow (non-zero water flux), one needs a gradient of either chemical, electrical, or other type. A structural asymmetry is not enough. We do not have any gradients in our system, the periodicity of the simulation box ensures that water above and below the membrane are connected, and we impose no electrical or other gradient on the system. Another way one could get non-zero water flux is dynamic coupling of the motions of the water to the channel. For example, waves could force water through in a prefered direction, but we find no evidence that there is such a dynamic coupling.

    Table 3. Total number of water crossing events, during the last 70 ns of simulation.
    # of Peptide Systems Total # of crossing downward crossing upward crossing
    tetramer LS2 252 127 125
    LS3 0 0 0
    pentamer LS2 16 8 8
    LS3 0 0 0
    hexamer LS2 3206 1605 1601
    LS3 2110 1000 1110
     | Show Table
    DownLoad: CSV

    For LS2 bundles, most water crossed through the hexamer arrangement, followed by the tetramer, and last the pentamer conformation (Table 3). In Figure 2B, the vector projection of the arrangement of the polar residues indicated that LS2 form a pore in the hexameric and tetrameric conformation, while the pentameric conformation contained a helix that is not aligned with the rest of the bundle that partially block the bundle. Thus, more water was observed traveling through the tetramer and hexamer bundle, while few waters (1-3 waters) crossed through the pentamer bundle (Figure 10). While the flow or water correlates with the pore radius, the pore radius isn't monotonic with numbers of oligomers in the bundle.

    Figure 10. Water molecules that crossed are shown traveling from 20 Å to −20 Å on the z-axis within the (A) LS2 tetrameric bundle and (B) LS2 pentameric bundle.

    In the case of LS3, most water traveled through the hexamer bundle, followed by pentamer and no water crossed through the tetramer bundle (Table 3). As expected, the pore is largest in the hexamer conformation, and smallest in the tetramer conformation (Figure 9). Unexpectedly from a topological point of view, LS2 allowed for more water to cross than LS3 in their respectable bundle arrangement despite the fact that LS3 peptides have more hydrophilic residues than LS2. However, as stated above, the water crossing events tracks with the pore radius (Figure 9). In the case of the tetrameric bundle, LS2 have an open “rectangle” conformation with a larger pore while LS3 have a close “diamond” structural arrangement [9]. No water was observed crossing through the LS3 tetrameric bundle indicating a closed channel while water was observed crossing though the LS2 tetrameric bundle.

    The LS2 pentamer bundle has large fluctuations in the pore radius (Figure 9B) and where all adjacent helices polar serine residues are not within a contact distance of one another [9]. Although the average narrowest region of the pore is 1.02 Å, the average pore radius of the bundle is 1.94 ± 0.55 Å where a few waters (1-3 water atoms) crossed through the pore at every nanosecond. On the contrary, the LS3 pentameric bundle has at least one serine-serine residues that are within contact with one another yet it resulted in a closed conformation that does not allow for any water to travel through the pore.

    Surprisingly, LS2 hexameric bundle has the largest pore (Figure 9C) as we initially thought that LS3 hexamer would have the largest pore because it had the most hydrophilic residues. Although LS3 contains more polar residues than LS2, it is apparent that the structural arrangements between helices play a larger role in water movement. Since adjacent helix within the LS2 bundle are not interacting as strongly as the LS3 bundle the LS2 bundle forms a much wider pore allowing for more water to cross. Hence, regardless of the number of polar residues in each helix, the movement of water is driven by structural arrangement of the bundle.


    3.3.3. Transition time of crossing and non-crossing events

    In addition to analyzing the numbers of water that crossed through the channel, the probability distribution of the water crossing events and the amount of time it took for water to travel through the channel were determined during the last 70 ns of simulation time (Figure 11). We find that water molecules traveled through the hexameric bundle the fastest, followed by the tetrameric and pentameric bundle in congruence with the number of water that crossed. Aside from analyzing the transition time of water that crossed through the “open” bundles, the transition time of “Non-crossing” events was also determined for all bundles. “Non-crossing” events refer to waters that entered the pore and traveled at least midway through the pore but do not cross. These water molecules exited and entered the pore on the same side of the bundle (Figure 12).

    Figure 11. (A) The probability distribution of the transition time that it takes for a water molecule to cross through the channel (solid line) fitted with a chi-squared distribution (dashed line). (B) Probability distribution of the number of water crossing events traversing through the “open” bundles (solid line) fitted with a Poisson distribution (dashed line).

    The probability of a water molecule crossing through the channel is hypothesized to be independent of its history, i.e. each water-crossing event is independent of both other water crossing events and past water crossing events. Such a process is characterized as a “Poisson process” and the total number of water crossing events within the simulation obeys the “Poisson distribution” (Figure 11B). Further information and derivations regarding the Poisson distribution can be found in probability theory textbooks such as in Larsen and Marx [52]. The probability distribution of water crossing events at every nanosecond is a good fit to a Poisson distribution, which is consistant with our hypothesis and indicates that the crossing events of water is indeed a stochastic process. The probability distribution of the transition time it took for water to cross and the transition time of “non-crossing” events were fitted to a chi-squared distribution (Figures 11A and 12). A chi-squared distribution is the sum of independent random variables (k) that measures for the goodness of fit of an observed distribution to a theoretical one [53,54]. The transition time of both the crossing and non-crossing events are in good agreement with the chi-squared distribution, which indicates that the movement of water is diffusive. Further, no openings and closing events were observed for each bundle whether by water crossing events or protein movement (Figures 11 and 12). In addition, the time it takes for water to cross through the entire pore and the time it takes for water to travel mid-way through the pore but do not cross are similar, which also is consistent with they hypothesis that water is diffusively traveling through the bundles with no cooperatively between waters or between water and the protein motions (Figures 11A and 12).

    Figure 12. Probability distribution of the amount of time it took for a “non-crossing” events to occur (solid line) fitted with a chi-squared distribution (dashed line) for all systems. The cold colors (violet, blue, cyan) represent LS2 systems and the warm colors (red, orange, yellow) represent LS3 systems.

    Since the LS2 hexameric bundle has a significantly larger average pore than the LS3 hexameric bundle, more water can travel through the LS2 bundle more readily and with a smaller transition time for “non-crossing” events (Table 4). Since LS3 forms a slightly smaller bundle due to forming stronger interactions with adjacent helices, water sits longer in the pore (higher transition time of “non-crossing” events), which does not allow water to flow as readily as the LS2 bundle.

    Table 4. Average pore radius profile and pore minimum radius during the last 70 ns with standard deviations, the average number of water that crossed through the pore with standard deviations, and both the rate of crossing and “non-crossing” events.
    # of Peptide Systems Avg. pore radius (Å) Avg. pore minimum radius (Å) Avg. # of crossing events per ns Crossing time (ns) (k) Non-crossing time (ns) (k)
    tetramer LS2 2.30 ± 0.45 1.54 ± 0.25 4 ± 2 1.51 (14) 1.36 (8)
    LS3 0.78 ± 0.34 0.175 ± 0.11 0 0 1.36 (9)
    pentamer LS2 1.94 ± 0.55 1.02 ± 0.39 2 ± 2 n/a* 1.13 (9)
    LS3 1.53 ± 0.61 0.427 ± 0.21 0 0 1.36 (8)
    hexamer LS2 4.72 ± 0.56 4.03 ± 0.23 43 ± 8 0.63 (8) 0.68 (5)
    LS3 3.88 ± 0.39 3.25 ± 0.26 26 ± 6 0.88 (9) 0.90 (6)
    *Requires more sampling; k is the number of degrees of freedom of the chi-squared distribution.
     | Show Table
    DownLoad: CSV

    For the pentameric bundle, LS3 formed a closed pore that allowed no water to traverse through the bundle. LS2, however, formed an asymmetric pore that allows few waters to cross, but not enough crossing events were sampled to accurately measure the transition time. From the probability distribution of “non-crossing” events, water traveled at a slower rate through the pentameric bundles compared to the hexameric and tetrameric bundles (Table 4). Also, the distribution of “non-crossing” events is not as pure “Gaussian” compared to all other systems. This is due to two narrow regions (<1.5 Å) found at the ends of pore where the pore is more confined at the N terminal than the C terminal (Figure 9B). Hence, it took less time for water to travel midway through the pore on the C terminal side than the N terminal, which in turn alter the shape of the distribution. In addition, it was observed that some “non-crossing” water molecules sit within the pore for long periods of time (>3 ns) due to the asymmetry of pore size. This can be visualized in Figure 10B, where there are few waters crossed through the LS2 pentameric bundle, but many waters was observed sitting within the pore for a long time. This give rise to the “tailing” at the end of the distribution of the “non-crossing” events transition time (Figure 12).

    For the tetrameric bundle, LS3 formed a closed pore where no water molecules traveled through the channel. Water traveled at a much slower rate through the LS2 tetramer bundle than the LS2 hexameric bundle due to the decrease in pore size (Table 4) as expected. The chi-squared distribution fitted to the transition time distribution of crossing events in the LS2 tetrameric bundle does not agree as strongly as the fit in the hexameric systems (Figure 11B). Since LS2 tetramer bundle formed a smaller pore, waters traveling through the bundle partially block the port and slow down the rate for other water molecules to travel, which resulted in a less diffusive water motion. The “non-crossing” transition time of water were found to be the same for waters traveling midway through the tetrameric bundles of LS2 and LS3.

    Figure 13. Number of water crossing events as a function of average pore area, and fitted with a linear function (dashed line).

    A linear fit was conducted on the average pore area and number of water crossing events, where the x-intercept (average pore area) was calculated to be 9.0 (Figure 13). This means that a minimum average pore radius of 1.70 Å is needed for a steady stream of water to cross through the bundle. The linear fit also provides a general equation where an estimate of water crossing events can be determined base on the average pore area.


    4. Conclusions

    Systems of homo-oligomeric transmembrane helical bundles ranging from N = 1-6 amphipathic peptides embedded in POPC lipid bilayer were explored using MD simulations. Our main objective was is to characterize environment around and within the bundles. Around the oligomers the degree of lipid perturbations induced by the transmembrane protein bundle was investigated. Within the oligomers we investigates the water dynamics and ion-channel characteristics of the bundles.

    Around the bundles, our analysis shows that lipids are only locally perturbed up to two solvation shells around the bundles/monomers. The properties of lipids three and more solvation shells away from the bundles/monomers are similar to those of the pure lipid bilayer systems. It was also found that the conformation of the LS2 and LS3 bundles, such as packing arrangements and symmetry etc., play a significant role in the perturbation of the lipids. We observe no significant meniscus formation of the lipid bilayer, i.e. no significant hydrophobic mismatch, in all systems. Therefore, lipid perturbation by the LS2 and LS3 monomer/bundles is mostly due to specific protein-lipid interactions. It is no surprise that they differ as the nature of the protein (LS2 vs. LS3) and the conformation of the bundles (size, symmetry, etc.) change. In lower order bundles that do not form a pore, since LS2 peptides tilted more and their polar residues are more exposed to the lipid hydrophobic core, a larger extent of lipid perturbations is observed. In higher order bundles (N = 4-6 and LS3 trimer), polar residues of both LS2 and LS3 peptides are shielded from the bilayer hydrophobic core and therefore lipid perturbation are not significantly affected by the tilt of the bundle. Finally, the lipid perturbation of the head groups are found to correlate mostly with the structure of the protein bundle and its resulting electric field by the hydroxyl group of the serine residues.

    Within the bundles, the movements of water traveling through the LS2 and LS3 bundles were examined. We calculated the average number contacts between serine-serine residues of adjacent helices for all bundles and analyzed the pore radius profile of each bundle from its initial and final structure where the positions of water through the bundle were identified. In all systems, it was found that LS2 had more water crossing events than LS3 in their respective bundle. The size of the pore correlates with the amount of water that can traverse through the pore. We find that the increase in polar residues will lead to stronger interactions amongst helices, leading to smaller pore formation that will decrease the amount of water crossing events.

    For LS2, the hexameric arrangement allowed for a great amount of water to travel through the bundle. So much so, that any ions can readily travel through the pore. In the pentameric arrangement, the channel allowed for the least number of water crossing events with the longest transition time of “non-crossing” events. LS2 is known to be a proton channel that requires a steady stream of water to flow through the pore, but the pore should have some selectivity. It is evident that the LS2 prefers to be in a tetrameric arrangement, which is in agreement with works done by Randa et al. [23] base on the transition time of water crossing through the channel, and the number of water crossing events. It is interesting to point out that from the vector projections it can be observed that the polar residues of the LS2 peptides are not facing towards the center of the bundle, but form a dimer of dimer arrangement [9]. Despite this, a steady stream of water was observed crossing through the bundle (Figure 10). For LS3, in a tetrameric and pentameric bundle, little or no water was observed traveling through the pore, hence the only ion channel conformation is hexameric. This is also in agreement with works done by Randa et al. [23] that claim that LS3 ion channel is in a hexameric packing arrangement.

    From our work, we demonstrated that the structural arrangement of the LS2 and LS3 peptides strongly correlate with the bilayer structure and movement of water through the pore. The lipids are perturbed within two solvation shells of the oligomers. The movement of water is a diffusive process through the pore. The transition time of both crossing and “non-crossing” events fit well to a chi-squared distribution that supports our argument. Further, we calculated that an average pore size of 1.7 Å is needed in order to have a steady stream of water flow within any channel. By comparing pore sizes of these simple channels to natural and complicated channel, an estimate of the number of water and the amount of time it takes for water to cross through the channel can be made, this in turn can shed light on the permeation rate of ions crossing through the channel.


    Acknowledgements

    The authors acknowledge that this research was supported in part by a gift from the H.O. West Foundation, a grant from the University of the Sciences in Philadelphia, a grant from the National Institute of Health (R15GM075990), Pennsylvania Department of Health (PDH) (SAP4100057688), a grant from the National Science Foundation (CHE1229564), and a graduate fellowship from NASA's Harriett G. Jenkins Pre-doctoral Fellowship Project. We also gratefully acknowledge Temple Owlsnest computer facility and Dr. Axel Kohlmeyer for computational resources.


    Conflict of interest

    All authors declare no conflict of interest in this paper.


    [1] Lee AG (2004) How lipids affect the activities of integral membrane proteins. BBA-Biomembranes 1666: 62–87. doi: 10.1016/j.bbamem.2004.05.012
    [2] Pohorille A, Schweighofer K, Wilson MA (2006) The origin and early evolution of membrane channels. Astrobiology 5: 1–17. doi: 10.1017/S1473550406002886
    [3] Hille B (2001) Ion Channels of Excitable Membranes, 3 Eds., Sinauer.
    [4] Kaczorowski GJ, Mcmanus OB, Priest BT, et al. (2008) Ion channels as drug targets: The next GPCRs. J Gen Physiol 131: 399–405. doi: 10.1085/jgp.200709946
    [5] Ackerman MJ, Clapham DE (1997) Ion channels-basic science and clinical disease. N Engl J Med 336: 1575–1586. doi: 10.1056/NEJM199705293362207
    [6] Lear JD, Wasserman ZR, Degrado WF (1988) Synthetic amphiphilic peptide models for protein ion channels. Science 240: 1177–1181. doi: 10.1126/science.2453923
    [7] Kienker PK, Degrado WF, Lear JD (1994) A helical-dipole model describes the single-channel current rectification of an uncharged peptide ion channel. Proc Natl Acad Sci USA 91: 4859–4863. doi: 10.1073/pnas.91.11.4859
    [8] Petrache HI, Zuckerman DM, Sachs JN, et al. (2002) Hydrophobic matching mechanism investigated by molecular dynamics simulations. Langmuir 18: 1340–1351. doi: 10.1021/la011338p
    [9] Nguyen THT, Liu Z, Moore PB (2013) Molecular dynamics simulations of homo-oligomeric bundles embedded within a lipid bilayer. Biophys J 105: 1569–1580. doi: 10.1016/j.bpj.2013.07.053
    [10] Howard KP, Lear JD, Degrado WF (2002) Sequence determinants of the energetics of folding of a transmembrane four-helix-bundle protein. Proc Natl Acad Sci USA 99: 8568–8572. doi: 10.1073/pnas.132266099
    [11] Arseneault M, Dumont M, Otis F, et al. (2012) Characterization of channel-forming peptide nanostructures. Biophys Chem 162: 6–13. doi: 10.1016/j.bpc.2011.12.001
    [12] Fischer WB (2005) Viral Membrane Proteins: Structure, Function, and Drug Design, In: Protein Rev, Kluwer Academic/Plenum Publishers.
    [13] Wang J, Kim S, Kovacs F, et al. (2001) Structure of the transmembrane region of the M2 protein H+ channel. Protein Sci 10: 2241–2250.
    [14] Stouffer AL, Acharya R, Salom D, et al. (2008) Structural basis for the function and inhibition of an influenza virus proton channel. Nature 451: 596–599. doi: 10.1038/nature06528
    [15] Schnell JR, Chou JJ (2008) Structure and mechanism of the M2 proton channel of influenza A virus. Nature 451: 591–595. doi: 10.1038/nature06531
    [16] Kovacs FA, Cross TA (1997) Transmembrane four-helix bundle of influenza A M2 protein channel: Structural implications from helix tilt and orientation. Biophys J 73: 2511–2517. doi: 10.1016/S0006-3495(97)78279-1
    [17] Acharya R, Carnevale V, Fiorin G, et al. (2010) Structure and mechanism of proton transport through the transmembrane tetrameric M2 protein bundle of the influenza A virus. Proc Natl Acad Sci USA 107: 15075–15080. doi: 10.1073/pnas.1007071107
    [18] Moore PB, Zhong Q, Husslein T, et al. (1998) Simulation of the HIV-1 Vpu transmembrane domain as a pentameric bundle. FEBS Lett 431: 143–148. doi: 10.1016/S0014-5793(98)00714-5
    [19] Woolley GA, Wallace BA (1992) Model ion channels: Gramicidin and alamethicin. J Membrane Biol 129: 109–136.
    [20] Opella SJ, Marassi FM, Gesell JJ, et al. (1999) Structures of the M2 channel-lining segments from nicotinic acetylcholine and NMDA receptors by NMR spectroscopy. Nat Struct Mol Biol 6: 374–379. doi: 10.1038/7610
    [21] Akerfeldt KS, Kienker PK, Lear JD (1996) Structure and conduction mechanisms of minimalist ion channels. Compr Supramol Chem 10: 659–686.
    [22] Gratkowski H, Lear JD, Degrado WF (2001) Polar side chains drive the association of model transmembrane peptides. Proc Natl Acad Sci USA 98: 880–885. doi: 10.1073/pnas.98.3.880
    [23] Randa HS, Forrest LR, Voth GA, et al. (1999) Molecular dynamics of synthetic leucine-serine ion channels in a phospholipid membrane. Biophys J 77: 2400–2410. doi: 10.1016/S0006-3495(99)77077-3
    [24] Oiki S, Danho W, Madison V, et al. (1988) M2 δ, a candidate for the structure lining the ionic channel of the nicotinic cholinergic receptor. Proc Natl Acad Sci USA 85: 8703–8707. doi: 10.1073/pnas.85.22.8703
    [25] Carruthers A, Melchior DL (1986) How bilayer lipids affect membrane protein activity. Trends Biochem Sci 11: 331–335. doi: 10.1016/0968-0004(86)90292-6
    [26] Palsdottir H, Hunte C (2004) Lipids in membrane protein structures. BBA-Biomembranes 1666: 2–18. doi: 10.1016/j.bbamem.2004.06.012
    [27] Lindahl E, Sansom MSP (2008) Membrane proteins: Molecular dynamics simulations. Curr Opin Struct Biol 18: 425–431. doi: 10.1016/j.sbi.2008.02.003
    [28] Phillips JC, Braun R, Wang W, et al. (2005) Scalable molecular dynamics with NAMD. J Comput Chem 26: 1781–1802. doi: 10.1002/jcc.20289
    [29] Duan Y, Wu C, Chowdhury S, et al. (2003) A point-charge force field for molecular mechanics simulations of proteins based on condensed-phase quantum mechanical calculations. J Comput Chem 24: 1999–2012. doi: 10.1002/jcc.10349
    [30] Hu Z, Jiang J (2010) Assessment of biomolecular force fields for molecular dynamics simulations in a protein crystal. J Comput Chem 31: 371–380.
    [31] Darden T, York D, Pedersen L (1993) Particle mesh Ewald: An N.log(N) method for Ewald sums in large systems. J Chem Phys 98: 10089–10092.
    [32] Feller SE, Yin D, Pastor RW, et al. (1997) Molecular dynamics simulation of unsaturated lipid bilayers at low hydration: Parameterization and comparison with diffraction studies. Biophys J 73: 2269–2279. doi: 10.1016/S0006-3495(97)78259-6
    [33] Jojart B, Martinek TA (2007) Performance of the general amber force field in modeling aqueous POPC membrane bilayers. J Comput Chem 28: 2051–2058. doi: 10.1002/jcc.20748
    [34] Taylor J, Whiteford NE, Bradley G, et al. (2009) Validation of all-atom phosphatidylcholine lipid force fields in the tensionless NPT ensemble. BBA-Biomembranes 1788: 638–649. doi: 10.1016/j.bbamem.2008.10.013
    [35] Rosso L, Gould IR (2007) Structure and dynamics of phospholipid bilayers using recently developed general all-atom force fields. J Comput Chem 29: 24–37.
    [36] Kučerka N, Tristram-Nagle S, Nagle JF (2006) Structure of fully hydrated fluid phase lipid bilayers with monounsaturated chains. J Membrane Biol 208: 193–202. doi: 10.1007/s00232-005-7006-8
    [37] Nielsen SO, Ensing B, Ortiz V, et al. (2005) Lipid bilayer perturbations around a transmembrane nanotube: A coarse grain molecular dynamics study. Biophys J 88: 3822–3828. doi: 10.1529/biophysj.104.057703
    [38] De Planque MR, Killian JA (2003) Protein-lipid interactions studied with designed transmembrane peptides: Role of hydrophobic matching and interfacial anchoring (Review). Mol Membr Biol 20: 271–284. doi: 10.1080/09687680310001605352
    [39] Nyholm TKM, Oezdirekcan S, Killian JA (2007) How protein transmembrane segments sense the lipid environment. Biochemistry 46: 1457–1465. doi: 10.1021/bi061941c
    [40] Sonne J, Jensen MO, Hansen FY, et al. (2007) Reparameterization of all-atom dipalmitoylphosphatidylcholine lipid parameters enables simulation of fluid bilayers at zero tension. Biophys J 92: 4157–4167. doi: 10.1529/biophysj.106.087130
    [41] Venturoli M, Smit B, Sperotto MM (2005) Simulation studies of protein-induced bilayer deformations, and lipid-induced protein tilting, on a mesoscopic model for lipid bilayers with embedded proteins. Biophys J 88: 1778–1798. doi: 10.1529/biophysj.104.050849
    [42] Chung LA, Lear JD, Degrado WF (1992) Fluorescence studies of the secondary structure and orientation of a model ion channel peptide in phospholipid vesicles. Biochemistry 31: 6608–6616. doi: 10.1021/bi00143a035
    [43] Choma C, Gratkowski H, Lear JD, et al. (2000) Asparagine-mediated self-association of a model transmembrane helix. Nat Struct Mol Biol 7: 161–166. doi: 10.1038/72440
    [44] Douliez JP, Leonard A, Dufourc EJ (1995) Restatement of order parameters in biomembranes: Calculation of C-C bond order parameters from C-D quadrupolar splittings. Biophys J 68: 1727–1739. doi: 10.1016/S0006-3495(95)80350-4
    [45] Douliez JP, Ferrarini A, Dufourc EJ (1998) On the relationship between C-C and C-D order parameters and its use for studying the conformation of lipid acyl chains in biomembranes. J Chem Phys 109: 2513–2518. doi: 10.1063/1.476823
    [46] Seelig J, Niederberger W (1974) Two pictures of a lipid bilayer. A comparison between deuterium label and spin-label experiments. Biochemistry 13: 1585–1588.
    [47] Seelig J, Waespesarcevic N (1978) Molecular order in cis and trans unsaturated phospholipid bilayers. Biochemistry 17: 3310–3315. doi: 10.1021/bi00609a021
    [48] Smart OS, Breed J, Smith GR, et al. (1997) A novel method for structure-based prediction of ion channel conductance properties. Biophys J 72: 1109–1126. doi: 10.1016/S0006-3495(97)78760-5
    [49] Smart OS, Neduvelil JG, Wang X, et al. (1996) HOLE: A program for the analysis of the pore dimensions of ion channel structural models. J Mol Graphics 14: 354–360. doi: 10.1016/S0263-7855(97)00009-X
    [50] Jorgensen WL, Chandrasekhar J, Madura JD, et al. (1983) Comparison of simple potential functions for simulating liquid water. J Chem Phys 79: 926–935. doi: 10.1063/1.445869
    [51] Zhong Q, Jiang Q, Moore PB, et al. (1998) Molecular dynamics simulation of a synthetic ion channel. Biophys J 74: 3–10. doi: 10.1016/S0006-3495(98)77761-6
    [52] Larsen RJ, Marx ML (2017) An Introduction to Mathematical Statistics and Its Applications, Pearson, 742.
    [53] Heijmans RDH, Pollock DSG, Satorra A (2000) Innovations in Multivariate Statistical Analysis, Springer US, 298.
    [54] Canal L (2005) A normal approximation for the chi-square distribution. Comput Stat Data An 48: 803–808. doi: 10.1016/j.csda.2004.04.001
  • This article has been cited by:

    1. Domenico Lombardo, Pietro Calandra, Maria Teresa Caccamo, Salvatore Magazù, Luigi Pasqua, Mikhail A. Kiselev, Interdisciplinary approaches to the study of biological membranes, 2020, 7, 2377-9098, 267, 10.3934/biophy.2020020
    2. Michael A. Wilson, Andrew Pohorille, Structure and Computational Electrophysiology of Ac-LS3, a Synthetic Ion Channel, 2022, 126, 1520-6106, 8985, 10.1021/acs.jpcb.2c05965
  • Reader Comments
  • © 2018 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(5702) PDF downloads(1004) Cited by(2)

Article outline

Figures and Tables

Figures(13)  /  Tables(4)

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog