Scholarly article on topic 'The Influence of Fatty Acids on the GpA Dimer Interface by Coarse-Grained Molecular Dynamics Simulation'

The Influence of Fatty Acids on the GpA Dimer Interface by Coarse-Grained Molecular Dynamics Simulation Academic research paper on "Biological sciences"

CC BY
0
0
Share paper
Academic journal
IJMS
OECD Field of science
Keywords
{""}

Academic research paper on topic "The Influence of Fatty Acids on the GpA Dimer Interface by Coarse-Grained Molecular Dynamics Simulation"

Int. J. Mol. Sci. 2014, 15, 14247-14268; doi:10.3390/ijms150814247

OPEN ACCESS

International Journal of

Molecular Sciences

ISSN 1422-0067

www.mdpi.com/journal/ijms

Article

The Influence of Fatty Acids on the GpA Dimer Interface by Coarse-Grained Molecular Dynamics Simulation

Nadine Flinner, Oliver Mirus and Enrico Schleiff *

Cluster of Excellence Macromolecular Complexes, Center of Membrane Proteomics,

Department of Biosciences, Molecular Cell Biology of Plants, GU Frankfurt am Main,

60439 Frankfurt, Germany; E-Mails: nadine-flinner@gmx.de (N.F.);

o.mirus@bio.uni-frankfurt.de (O.M.)

* Author to whom correspondence should be addressed; E-Mail: schleiff@bio.uni-frankfurt.de;

Tel.: +49-69-798-29287; Fax: +49-69-798-29286.

Received: 22 May 2014; in revised form: 14 July 2014 /Accepted: 6 August 2014 /

Published: 15 August 2014

Abstract: The hydrophobic thickness of membranes, which is manly defined by fatty acids, influences the packing of transmembrane domains of proteins and thus can modulate the activity of these proteins. We analyzed the dynamics of the dimerization of Glycophorin A (GpA) by molecular dynamics simulations to describe the fatty acid dependence of the transmembrane region assembly. GpA represents a well-established model for dimerization of single transmembrane helices containing a GxxxG motif in vitro and in silico. We performed simulations of the dynamics of the NMR-derived dimer as well as self-assembly simulations of monomers in membranes composed of different fatty acid chains and monitored the formed interfaces and their transitions. The observed dimeric interfaces, which also include the one known from NMR, are highly dynamic and converted into each other. The frequency of interface formation and the preferred transitions between interfaces similar to the interface observed by NMR analysis strongly depend on the fatty acid used to build the membrane. Molecular dynamic simulations after adaptation of the helix topology parameters to better represent NMR derived structures of single transmembrane helices yielded an enhanced occurrence of the interface determined by NMR in molecular dynamics simulations. Taken together we give insights into the influence of fatty acids and helix conformation on the dynamics of the transmembrane domain of GpA.

Keywords: Glycophorin A dimerization; dimer interface formation; bitopic transmembrane a-helix; fatty acid dependency; MARTINI force field

1. Introduction

All living cells are surrounded by membranes, which are mainly composed of glycerophospholipids, sphingolipids and sterols [1,2], which consist of a hydrophilic head group and a hydrophobic tail. The lipids are arranged in two layers within the membrane, with the tails packed against each other. Glycerophospholipids are the most prominent components of biological membranes and phosphatidylcholine (PC) is the most frequently occurring lipid in biological systems [3].

By molecular dynamics (MD) simulations it has been shown that PC lipids with saturated and unsaturated fatty acid chains separate spontaneously into an liquid ordered phase, containing the saturated lipids and cholesterol, and an liquid disordered phase, containing the lipids with unsaturated fatty acids [4,5]. A small mismatch between the saturated and unsaturated fatty acids reduces the driving force for segregation, while the increase of cholesterol leads to an enhanced driving force for the phase separation [6]. Phase separation and compartmentalization are also present in biological membranes, for example lipid rafts are observed in the eukaryotic plasma membrane, which are enriched in sphingolipids and cholesterol [7].

Native lipid bilayers are generally crowded with transmembrane proteins [7-9], which influence and which are influenced by the biophysical properties of the membrane. By MD simulations it has been shown that bitopic a-helices are sorted into the liquid disordered phase [6]. This results from the enhanced number of possible van der Waals contacts [10], independently of the hydrophobic thickness of the different phases [11]. However, some proteins, e.g., peripheral membrane proteins H-Ras or Hedgehog are sorted into the liquid ordered phase, which depends on their lipid anchor [12].

Irrespective of whether the transmembrane region is sorted into the liquid disordered or ordered phase, it is generally observed that the functional activity of membrane proteins depends, among other things, on the hydrophobic thickness and the composition of the membrane [13]. The highest activity of proteins is often observed in case that no hydrophobic mismatch between the protein and the membrane occurs [14]. To compensate for the hydrophobic mismatch, a-helical transmembrane regions can change their tilt angle, the conformation of side chains or their effective hydrophobic length [15]. Alternatively, the membrane can adapt its thickness or can undergo a midplane bending to scope with the transmembrane domain dimension [10,13]. To minimize the energetic penalty for this mismatch, membrane proteins often aggregate after occurrence of a hydrophobic mismatch [16]. Consequently, changing the fatty acid chain length would have significant effects on the packing of the a-helices [10].

Glycophorin A (GpA) is one of the model peptides used to investigate the biophysical properties of membrane domain behavior in lipid environment. GpA is a human protein of 150 amino acids with a single transmembrane helix and is localized in the plasma membrane of red blood cells. This membrane consists of phospholipids, cholesterol and glycolipids. Among the phospholipids PC, phosphatidylethanolamine, and phosphatidylserine are most abundant and the fatty acids thereof have a

length of 16 to 22 carbon atoms [17,18]. GpA carries an O- and N-glycan, which again carries sialic acid. This negative charge prevents aggregation of the red blood cells and enables the flow through the blood vessels [19]. To exhibit its function, GpA forms a homodimer as well as a heterodimer with Glycophorin B, which then forms a multimeric protein complex with other proteins like band 3 [18,20]. The structure of the homodimerizing membrane inserted segment of GpA in micelles [21,22] as well as in bicelles [22] was determined by NMR spectroscopy. The interface between the two helices consists of a GxxxG motif, a common motif for helix-helix association [23]. GpA forms a right-handed dimer with a crossing angle of approximately -40° [21]. The interface does not contain hydrogen bonds between the two monomers, while the OH group of threonine localized in the dimerization interface forms an intramolecular hydrogen bond with the backbone [21].

Although GpA dimerization has been studied in silico by MD simulations (e.g., [24-26]), the dependence of the dimeric interface on fatty acid chain length has not yet been investigated in detail. Thus, we examined the influence of different fatty acid chain lengths and saturation degrees on the stability and dynamics of the GpA homodimer by coarse-grained molecular dynamics (CG-MD) simulations with the MARTINI force field [27-29]. We performed self-assembly simulations with the transmembrane region of GpA and analyzed whether initial interface formed after assembly depends on the fatty acid composition of the membrane. Subsequently, we analyzed the fatty acid dependence of the transitions between different dimer interfaces as well. Finally, simulations on the stability of the dimer observed by NMR have been performed to explore the importance of the correct helix topology on the stability.

2. Results

2.1. Self-Assembly of the GpA-Transmembrane Region Reconstructs the NMR Interface with Low Frequency

We used the transmembrane part of GpA and performed CG-MD simulations in five different membranes to investigate the influence of the fatty acids on interface formation of the dimeric a-helical protein GpA. Each membrane consisted of lipids with different fatty acid chains, while the PC-head group was kept, because it is the most frequent head group in the plasma membrane [3], the natural environment of GpA [18,20]. We used the saturated lipids DPPC (dipalmitoyl-PC) and DSPC (distearoyl-PC) with different length, where each fatty acid is modeled by four and five CG beads representing saturated fatty acids of 15-18 or 18-21 carbons respectively. Further, a lipid with one saturated and one monounsaturated fatty acid chain (POPC: palmitoyl-oleoyl-PC) is used, with the saturated fatty acid modeled by four and the unsaturated one by five CG beads. In addition, membranes composed of the unsaturated DOPC (dioleoyl-PC with both fatty acids modeled by five CG beads) and the di-unsaturated DUPC (dilinoleyl-PC with both fatty acids modeled by four CG beads) have been included in our analysis. The hydrophobic thickness of the according membranes ranges from 2.49 ± 0.04 nm for DUPC to 3.86 ± 0.04 nm for DSPC (DUPC < DPPC (3.13 ± 0.09 nm) < POPC (3.22 ± 0.03 nm) < DOPC (3.42 ± 0.02 nm) < DSPC).

In each of the created lipid environments 300 self-assembly simulations for GpA were performed using the Martini 2.2 force field [29]. The individual simulation was terminated after formation of a stable dimer (methods) of the transmembrane regions. At first the crossing angles of the dimers was

analyzed yielding mainly negative crossing angles characteristic for right-handed dimers in the thin membranes, while in thicker membranes crossing angles around 0° have been observed as well. In contrast, positive crossing angles characteristic for left-handed dimers were hardly observed in any simulation. Consistently, the crossing angle observed by NMR analysis is -40° [21].

In general, most of the self-assembled interfaces are located at one side of the helix. Surprisingly the known dimerization interface, containing the amino acid T87, is located on the opposite side of these frequently obtained structures (Figure 1A). Next, we clustered the structures based on their interface contact plots (methods). We used all structures from self-assembly simulations mentioned in this manuscript to determine the "cluster space" and to assign identical cluster numbers to interfaces obtained with different approaches. In total 83 clusters of different interface architectures are observed.

The interfaces most frequently observed considering all short self-assemblies simulations are represented by cluster 52, which at the same time is the most frequently observed cluster for POPC, DOPC and DSPC (Figure 1B). The dimerization interface of the in cluster 52 is located at a different site than the interface of the NMR structure (Figure 1C) and has an average tilt angle near zero (-0.95° ± 6.58°; Figure 1D). The most frequent interface observed in the DPPC bilayer is represented by cluster 63, which is at the same time the second most frequently observed cluster for all short self-assemblies simulations (Figure 1B,C). The crossing angle of the dimer of the structures in this cluster is 1.17° ± 6.92° (Figure 1D). The most frequently occurring interface in the DUPC bilayer is cluster 81 (in total the third most frequently observed cluster; Figure 1B,C), which contains an asymmetric right-handed dimer with a crossing angle of -32.05° ± 4.18° (Figure 1D).

The structures forming the interface observed by NMR analysis [21] are represented by cluster 14, which has been observed for all membranes analyzed. The dimeric structures therein are right-handed with an average crossing angle of -30.81° ± 3.46°, which is slightly smaller when compared to the NMR structure (Figure 1D). Interestingly, although the interface determined by NMR analysis is generally detected with low frequency, it is most frequent in the bilayer with the smallest hydrophobic thickness (DUPC; Figure 1B).

To determine whether another interface would also explain the distances obtained by NMR experiments, the coarse grain structures from the short self-assembly simulations were backmapped into a full atomistic resolution. Subsequently we checked which of the NOE distances are fulfilled by the structures of each cluster (Figure 1E). As expected, the structures of cluster 14 match best to the NOE distances with on average 19.7 constraints fulfilled by one structure, while other clusters match some of the NOE's even better than the NMR structure. For example, there are seven NOE constraints of GpA measured in bicelles or micelles [22] which are fulfilled most often by the structures from cluster 12 (Figure 1E). This cluster represents an interface with a crossing angle of -28.56° ± 5.50° (Figure 1D) and which contains contacts for the GxxxGxxxT motif similar to that found by NMR, but which is slightly asymmetric (Figure 1C). These seven NOE constraints are also fulfilled by structures of clusters 16 (crossing angle: -26.06° ± 6.45°) and 33 (crossing angle: 29.38° ± 9.68°) (Figure 1E), where at least the glycine residues of the GxxxG motive have contacts with each other (Figure 1C). However, the crossing angle of cluster 33 is not comparable to the one observed by NMR analysis. Thus, clusters 12 and 16 are annotated as "NMR-like" and represent together with the NMR cluster 14 about 10.6% of all self-assembled structures, and thus still not a large portion of all observed interfaces.

Figure 1. Short self-assembly simulations (25 mer) (A) Density map of all self-assembled dimers in top view, superimposed on BB beads of helix A. The "*" marks the position of T87 of the first helix; (B) The frequency of the different interfaces obtained for each lipid (DUPC...blue, DPPC...green, POPC...red, DOPC...cyan, DSPC...grey) after UPGMA clustering is shown in total numbers. Only clusters containing at least three structures of different frames (0.1%) are shown, all structures from other clusters are shown in the category "oth." (others); (C) Shown are residue-residue contact maps for all structures in the corresponding cluster. The frequency of contacts is shown as grey scale with black indicating contact present in all structures of the cluster and white indicating contact not present. Framed in red are contacts observed by NMR analysis [21]. For representative CG-structures please see Figure S1; (D) The box plot shows the distribution of the crossing angles of the structures in each of the clusters shown in B. A box contains 50% of the data values and marks the 75th percentile and the 25th percentile, whereas the red line marks the median. The whiskers (horizontal lines) mark the values, which lie within the 1.5 interquartile range, which is defined as the box length; further values are marked as outliers ("+"); (E) Agreement of the backmapped structures with interhelical NOE constraints from NMR experiments. On the X-axis the single NOE constraints from four experiments (Smith51: micelle or POPC; 1afo [21]: micelle; 2kpf [22]: bicelle; 2kpe [22]: micelle) and on the T-axis the percentage of structures in the corresponding cluster which fulfill this distance are plotted. Shown are only structures with a frequency at least as high as the NMR cluster.

2.2. Interface Frequencies Are Depend on Different Fatty Acid Chain Environments

Next, the dependency of the preferred interface on the fatty acid chains was analyzed in more detail. At first, the probability for the occurrence of a novel interface after a defined number of simulations performed was determined. This analysis revealed a probability of 5% to find a new interface after 150 simulations, which in turn accounts for a 95% confidence for complete sampling. Furthermore, a 98% confidence is reached for the detection of all possible structures after 300 simulations. Interestingly, interfaces that occur only after 150 simulations are very unstable and short living, and most likely represent intermediates between two stable interfaces. This in turn indicates that we have sampled nearly the whole interface space.

Further we investigated, whether we have sampled enough to determine differences in the cluster frequencies observed in the different lipids. We calculated the autocorrelation for the function of cluster frequencies using Pearson's correlation coefficient (r) for n (n = 5 to 150) randomly selected interfaces of our 300 simulations without replacement and fitted a double exponential function to our data (Figure 2A). In general, correlation and autocorrelation functions decay or rise with several numbers of exponential terms, related to the number of underlying processes [30,31]. In our case, the first term refers to the number of simulations to detect all important clusters and the second term refers to the number of simulations to cover the distribution. The limits of the fitted functions are in the range of ~0.95 (Table 1). With n = 150 simulations r is already in the range of 0.93 (DPPC) to 0.96 (DOPC), reflecting that we have a high overlap for two distinct sets and that the sampling is nearly fully completed (Figure 2A). Thus, consistent with the probabilities to detect a new interface, we conclude that the sampling nearly covers the underlying probability distribution of the different clusters for all lipids.

Table 1. Optimized values for a, b, c, d and the limit of the fitted (auto-) correlation function.

Finally, we calculated the correlation between different lipids, again using Pearson's correlation coefficient for n (n = 5 to 300) interfaces (Figure 2B). The correlation between the three longest lipids (POPC, DOPC and DSPC) is very similar to the autocorrelation. This indicates that there is a large

a b c d Lim

DUPC-DUPC DUPC-DPPC DUPC-POPC DUPC-DOPC DUPC-DSCP DPPC-DPPC DPPC-POPC DPPC-DOPC DPPC-DSPC POPC-POPC POPC-DOPC POCP-DSPC DOPC-DOPC DOPC-DSPC DSPC-DSPC

0.51 0.10 0.42 0.02 0.93

0.37 0.07 0.18 0.01 0.56

0.24 0.09 0.12 0.01 0.36

0.18 0.11 0.07 0.02 0.25

0.18 0.10 0.09 0.02 0.27

0.55 0.07 0.38 0.01 0.94

0.58 0.07 0.24 0.01 0.82

0.58 0.08 0.22 0.01 0.80

0.57 0.07 0.22 0.01 0.79

0.59 0.09 0.36 0.02 0.95

0.73 0.08 0.20 0.01 0.93

0.67 0.08 0.25 0.01 0.92

0.65 0.13 0.31 0.02 0.96

0.73 0.09 0.19 0.01 0.93

0.59 0.11 0.35 0.02 0.95

overlap between the distributions of the different interface clusters found for these lipids. The comparison of DUPC (shortest fatty acid chain analyzed) with POPC, DOPC and DSPC shows correlation values between 0.27 and 0.36, whereas the correlation of DPPC (longer than DUPC but shorter than POPC, DOPC and DSPC) in comparison to POPC, DOPC and DSPC is much higher with a value of about 0.80. The value for the correlation between DUPC and DPPC is intermediate with 0.56.

2.3. Transitions between Dimer Conformations Are Visited by Long Self-Assembly Simulations

The results presented above suggest that the preferred interface of GpA depends on the nature of the fatty acid used for membrane construction. However, the interface found by NMR experiments is observed only with low frequency in self-assembly simulations, irrespective which fatty acid is used. This prompted the question whether the interfaces observed are only transition states to a putative final interface. To this end we conducted long self-assembly simulations (10 ^s) with focus on bilayers composed of DUPC (10 simulations), DPPC (10 simulations) and POPC (20 simulations). DOPC and DSPC containing membranes were not analyzed because the interface distribution observed in these membranes is comparable to the one in POPC (Figure 2). Additional, the long simulations enabled us to analyze the transitions between interfaces.

Figure 2. (A) Correlation and (B) Autocorrelation for the short self-assembly simulations were calculated using Pearson's correlation coefficient (7-axis) between the cluster distributions obtained from n simulations (X-axis). Plotted are the average values obtained from 1000 repetitions (sampling without replacement).

In general, all clusters observed with a frequency of at least 3% in the short self-assembly simulation of GpA in a membrane composed of one specific lipid were also observed in the long self-assembly simulation in the according membrane (compare Figure 1B with Figure 3A). This indicates that we have generally sampled the same interface space with the short and long simulations.

For DUPC and POPC even all clusters with a frequency of at least 1% were visited, indicating that the sampling is more complete in membranes composed of these two lipids. However, we observed a shift of some cluster frequencies in comparison to the short self-assembly, which is caused by a high retention time e.g., in cluster 79 or 81 (Figure 3B).

Figure 3. Long self-assembly simulations (25 mer) (A) Frequency of the different interfaces obtained in the long self-assembly simulations (10 ^s) with the 25 mer for each lipid (DUPC...blue, DPPC...green, POPC...red) after UPGMA clustering. Shown are the same clusters as in Figure 1B; (B) Retention times of different clusters in different lipids; (C) Transition map for the self-assembly simulations in DUPD and POPC. Each cluster with a frequency at least equal to the cluster of the NMR structure is represented by a circle, please note that the NMR cluster (12) is not found in POPC. Only the three most frequent transitions for each cluster are shown, where transitions occurring only once were excluded. Transitions occurring frequently only in DUPC were colored in blue, only in POPC in red and in DUPC and POPC in purple.

Similar to the results for the short self-assembly simulations, we observe differences in the cluster frequencies of different fatty acid chains. For example, the interface of cluster 52 is more frequent in POPC and, remarkably, the cluster, which contains the contacts of the NMR structure (cluster 14), occurs only in DUPC and DPPC (Figure 3A). Interestingly, all frequently occurring interfaces convert into each other—at least via some intermediate steps—meaning that the underlying graph is strongly connected.

In general, the transition maps for DPPC and POPC are similar (Figure S2) and all transitions which were frequently observed in DPPC are also present in POPC, but due to the limited sampling and higher retention times we have less confirmed transitions for DPPC. The transitions of DUPC and POPC are in general similar as well, but a significant difference with respect to the outgoing transitions from NMR-like clusters 12 and 16 exists (Figure 3C). In DUPC the most frequent transition of cluster 12 and 16 points to the NMR interface, whereas in POPC this transition is not observed. Similarly, the second or third most transition points to the other NMR-like cluster, which is also not observed in POPC. Hence, the NMR-like interfaces are observed while simulating the dimer in both lipids, while the transition to the NMR interface is only observed in DUPC. Thus, we do not only see an influence of the fatty acid chain onto the frequency of different interfaces, we also see an influence of the preferred transitions between different clusters.

2.4. Stability of the NMR Interface Depends on Helix Length and Helix Conformation

Next, we aimed to explore the reason for the underrepresentation of the NMR interface in our simulations. At first we tested the stability of the NMR dimer in the force field used. We embedded the structure determined by NMR analysis [21] truncated to the 25 mer into membranes of DUPC, DPPC and POPC and simulated each system three times using different seeds. We realized that the NMR interface is not stable for 10 p,s in any of the simulations. Interestingly, the interface changes rapidly in DUPC, however the NMR interface reforms in all three repetitions (Figure 4A). Further, it appears that the stability of the dimer interface increases in thicker bilayers and is more stable in bilayers composed of saturated fatty acid chains (DPPC) in comparison to bilayers containing unsaturated fatty acid chains (DUPC, POPC, Figure 4A). We repeated the simulation with the construct used for NMR analysis (40 amino acids) to test if the interface is stabilized by flanking residues. In addition, indeed we observed a stabilization of the NMR interface (Figure 4A vs. Figure 4B). However, as found for the 25 mer, the stability of the interface of the 40mer increase in saturated or thicker membranes (Figure 4B).

In parallel, we used the MARTINI force field version Martini2.1 [28] and GpA2010 [26], because a successful reconstruction of the NMR interface was reported while using GpA2010. For GpA2010 in five of nine simulations the interface is maintained over 10 p,s, while the same was observed for only two of nine simulations using the Martini 2.2 version (Figure 4B). In contrast, the Martini 2.1 version is not able to stabilize the interface at all, and in most cases the NMR interface is lost during the first microsecond (Figure 4B). However, the comparable behavior in maintaining the NMR interface while using the GpA2010 and Martini 2.2 force field is not reflected by the force field parameters. In general, Martini 2.1 and Martini 2.2 use by and large the same angles, dihedrals and bead types to model the properties of the side chain and the backbone [28,29]. In contrast, GpA2010 uses different bead types to describe non-bonded interactions and also the definition of the conformation of a-helix

differs drastic as GpA2010 places the backbone beads to the CA position and not the center of mass of the amino acid backbone. Through this change different bond length and bond angles are needed. Furthermore, the secondary structure is maintained by dihedral angles in Martini 2.1/Martini 2.2, while in GpA2010 an additional bond is placed between backbone beads i and i + 4. However, one major difference between Martini 2.1 and Martini 2.2 is the bond length between the backbone beads (Martini 2.1: 0.35 nm vs. Martini 2.2: 0.31 nm [28,29]). This suggests that a suitable definition of the helix is at least as important as the definition of the interaction energies between different beads.

Figure 4. The similarity (similarityNMR) of the interfaces of the structures from trajectories obtained with different lipids (DUPC...blue, DPPC.. .green, POPC...red) and (A) the 25 mer (B) the 40 mer and different force field versions (GpA2010...first panel, Martini2.1...second panel, Martini2.2...third panel) to the one determined by NMR spectroscopy in micelles (1afo, [21]) is plotted against the simulation time (time). For each lipid + force field combination three simulations started with different seeds are shown. A value of one indicates that all residue-residue contacts present in the NMR structure are reproduced, while zero indicates no reproduced contact.

2.5. Correction of the Helix Conformation Leads to Increased NMR Interface Frequency

Consequently, we compared the helix length obtained from NMR experiments of different bitopic a-helices able to form dimers with the values obtained in our simulation. Consistent with earlier reports [29] we realized that the helix length within the Martini 2.1 force field is significantly longer than found by experimental approaches (Figure 5A). In turn, the experimental helix length distribution

overlaps with the distribution in the other two force fields, although the length is slightly shorter (GpA2010) or longer (M2.2) in the force fields. This supports our thesis that the helix conformation is—next to the non-bonded interactions—the second determinant for a correct interface formation.

In order to test whether the helix topology has an impact on the dimerization behavior we used an elastic network-like approach similar to the normal elastic network used with the MARTINI force field [32] to model the helix topology. For this procedure we add additional bonds between the backbone atoms i and i + 3, i and i + 4 and i and i + 5 with an equilibrium length corresponding to the average observed in the known structures (0.49, 0.61 and 0.82 nm) and screen all combinations for force constants e {None, 500, 1,000, 5,000, 10,000} (Table 2). We realize that adding only a bond between backbone atom i and i + 3 with a force constant of 5,000 performs best in reproducing the helix length of known structures (Figure 5B). Using this additional bond we repeat the simulations for the stability of the 25 mer NMR interface and noticed that behavior of the NMR interface of the 25 mer of GpA is not altered drastically, as a reinvestigation of the NMR interface as well as the increased stability in thicker membranes can still be observed (Figure 5B). Furthermore, the retention time of the NMR interface is slightly enhanced compared to the original Martini 2.2 force field (Figure 4A vs. Figure 5B).

Figure 5. Influence of the helix conformation onto the simulation. (A) probability distribution for the helix length obtained from solved structures and from simulations using the COM position (solid lines) or the Ca position (dashed lines). The helix length is measured between backbone atom i and i + 12. Black/grey...NMR reference, red...M2.2, orange...M2.1, green...M2.2 plus bbi-bbi + 3 bond, blue...GpA2010; (B) Stability of the NMR interface of the transmembrane region using the modified M2.2 topology with additional bbi-bbi + 3 bonds. The figure is represented as described for Figure 4.

Table 2. Screen for a proper force constant for backbone atoms i and i + 3/i + 4/i + 5 for GpA in DUPC. An additional screen in DPPC results in the same minimum. Shown is the difference between the probability distribution obtained from the NMR structures and the corresponding simulation.

None 500 1000 5000 10000 bb3

None 133.07 115.7 100.2 60.34 57.02

500 137.10 123.6 114.4 72.56 62.92

1000 121.99 114.0 104.0 70.09 62.57 None

5000 83.59 84.61 82.69 73.35 71.52

10,000 75.11 77.43 76.97 76.83 79.96

None 71.72 61.61 51.89 40.38 44.57

500 92.88 83.73 70.28 51.18 50.20

1000 81.57 78.54 69.55 51.53 52.24 500

5000 58.17 60.20 60.39 61.08 63.49

10,000 64.70 66.86 67.18 73.00 76.51

None 54.63 56.27 45.74 39.17 47.65

500 87.97 79.35 69.90 50.15 50.48

1000 75.77 72.60 65.51 51.49 52.83 1000

5000 56.98 60.49 59.77 60.89 63.18

10,000 64.17 67.19 66.17 72.05 76.92

None 27.38 32.15 34.26 46.53 57.14

500 48.72 50.96 49.47 50.73 56.61

1000 49.67 49.56 49.15 51.87 59.93 5000

5000 53.08 56.04 55.63 63.25 71.04

10,000 67.34 68.75 68.40 75.04 80.79

None 32.34 33.56 37.52 55.03 64.71

500 38.60 41.83 45.41 53.83 65.94

1000 40.73 44.60 46.12 57.61 66.66 10,000

5000 59.00 60.18 59.61 69.77 75.88

10,000 71.35 71.95 73.49 79.16 84.68

Table 3. Frequencies (freq) and placing for interface clusters (nr) from all short self-assembly simulations in DUPC and DPPC obtained with different force field versions (m22, GpA2010 and m22 + bb3 - m22 together with additional bond between backbone atom i and i + 3).

DUPC DPPC

Place m22 m22 + bb3 GpA2010 m22 m22 + bb3 GpA2010

nr freq nr freq nr freq nr freq nr freq nr freq

1 81 17.0% 81 13.3% 2 31.3% 63 15.0% 52 19.7% 2 15.7%

2 76 11.7% NMR 10.3% 0 19.3% 52 11.0% 63 12.7% 0 11.7%

3 16 11.3% 76 9.3% 12 12.3% 79 8.7% 42 8.0% 7 10.0%

4 79 6.3% 33 6.7% 16 4.7% 42 7.3% 76/79 5.0% 6 6.7%

5 NMR 6.0% 16/79 6.0% 6/9 3.0% 81/12 6.0% 12 4.3% 12 6.0%

These results have been promising and thus, we repeated the short self-assembly simulations to test whether the frequency of the NMR interface is increased. Indeed the frequency of NMR interface has nearly doubled in DUPC (Table 3) and resembles the second most frequent cluster. Interestingly, DUPC is the lipid which hydrophobic length comes close to the detergent fatty acid which was used for structure determination. Furthermore, using the adopted force field we still observed a fatty acid dependent shift in cluster frequencies (Table 3). This documents that the fatty acid dependency of the interface formation and stability is reproducible and not strongly dependent on the helix topology.

3. Discussion

The low frequency of the detection of the interface observed by NMR analysis also prompted a methodological discussion. Several technical aspects might contribute to this observation. (i) The MARTINI force fields use beads with an effective size of o = 0.47 nm, which stands in contrast to the minimal distance between beads in different bitopic a-helices of dimers (smaller than 0.47 nm; Table S1). Furthermore, an inter-helical hydrogen bond between the threonine side chain of the first helix and the backbone of the second helix exists in membranes [33] that cannot be modeled with the current MARTINI force field. Here, the announced parameter for polarizable particles for the backbone [34] might provide a solution for future analysis; (ii) The helix conformation might change in thicker membranes, and the in here constrained helix topology adapted to values observed in detergent micelles would not be correct. However, arguing against a major impact of such phenomenon, our simulations show an increased retention time of the NMR interface in thick bilayers (Figures 4 and 5) parallel FRET experiments, by which the highest dimer fraction is found in membranes with 20:1 and 22:1 fatty acids, while the observed dimer fraction is lower in thinner membranes (14:1, 16:1 and 18:1) [35]. However, the correct parameterization of bonded parameter strongly influences the dimer interface formation. There is a magnitude of different elastic network methods to restrain the overall structure of proteins like SAHBNET [36], IDEN [37] or ElNeDyn [32]. Testing all was beyond the scope of this manuscript, because we were interested in the question if different fatty acid chains influence the dimerization interface and not primarily in the explicit structures formed. Nevertheless, reproducing the helix conformation known from NMR experiments by application of an additional bond, which represents a simple elastic network with a small distance cutoff (Rc = 0.49 nm) and a big spring force constant (Kspring = 5000), significantly enhanced the frequency of the NMR interface, while the effects of the fatty acids are still observable (Figure 5). Thus, the major finding presented in here concerning the fatty acid dependence of dimer formation is not dependent on the correct parameterization of bonded parameter.

The fatty acid dependence of the dimerization of GpA was explored using membranes composed of different fatty acid chains with PC head group. The use of the PC head group is rationalized as (i) conclusions on fatty acid impact requires a constant head group behavior; (ii) NMR studies on GpA have been conducted in DPC micelles or DMPC/DHPC bicelles [22]; and (iii) PC is one of the most abundant head group in the plasma membrane of erythrocytes [17,18]. The analysis of the different fatty acids was prompted by the observation that the dimer of "serine zipper" peptides is more stable in thick bilayers (POPC) when compared to thin bilayers (DLPC) or micelles [38]. In line, the GpA dimer observed by NMR was better reproduced by atomistic MD simulations in a very thin bilayer (DMPC)

than in a micelle [39]. The CG-MD analysis of the dynamic behavior of the GpA membrane anchor in bilayers composed of lipids with nearly physiological fatty acid chain length reproduced the interface observed by NMR in all membranes (Figure 1), but with highest frequency in membranes with a small hydrophobic thickness (DUPC, Figures 1 and 3). In turn, stabilization of the dimer with characteristic negative crossing angles and NMR contacts depends on proper hydrophobic thickness of the membrane. Consequently, in thin membranes (like DUPC) no stabilization occurs and in very thick membranes (like POPC) the helices are not able to form the characteristic contacts and the negative crossing angle of the NMR interface as the tilt of the monomers is too small.

However, many different interfaces have been identified by CG-MD self-assembly analysis, which parallels the proposed large conformational space for the GpA dimer based on surface-based predictions [40]. Thus, it is tempting to speculate that the observed alternative dimeric conformations reflect components of the large structural ensemble in native membranes. Supporting this hypothesis thermodynamic measurements of single amino acid mutants of GpA indicate that the GxxxG motif is not necessary for dimerization of GpA [41]. Further, different interfaces can be adopted by the same sequence as exemplified for different bitopic a-helical proteins: For example for VEGFR-2 (vascular endothelial growth factor receptor 2) was described that an inactive and an active interface can be adopted, with helices rotated along the helical axis by 180° [42]. Similarly, EphAl (Erythropoietin-producing hepatocellular receptor) forms a dominant right-handed dimerization interface, but a left-handed dimer interface was observed in NMR experiments as well [43]. The same two dimerization motifs were found for EphA2 and the preference for the interface depends on the hydrophobic thickness of the membrane [44]. Finally, APP (amyloid precursor protein) forms different interfaces of the transmembrane dimer as well, as identified by NMR measurements of the wild type transmembrane sequence in POPC:POPS:CHOL and POPC:POPS membranes [45].

Taken together, it will be very interesting to validate our findings by invitro measurements, for which different strategies might be proposed like NMR measurements of dimerizing GpA mutants with disturbed GxxxG motif (as in [42,45]) or of wild type GpA in thicker membranes (as in [45]). The best candidate for such an alternative interface is represented by cluster 81, which is the most stable cluster in long self-assemblies and which is also found after adaptation of the parameter for the definition of the helix conformation. The most important amino acids of this interface are I77, F78, M81 and A82 and so the investigation of the dimerization behavior of mutants disturbing both (the GxxxG motif and the interface of cluster 81) is an interesting experiment.

Remarkably, the occurrence of the different interfaces, their retention time and their transitions (Figures 1-4) are strongly fatty acid chain dependent. For example the retention time of the most frequent interface (cluster 81) from long self-assembly simulations is highest in DPPC (Figure 3). This finding is in agreement with the observation that the Potential of Mean Force (PMF) for GpA association is highest in intermediately thick bilayers like DPPC and lower in thicker (DOPC) and thinner (DMPC) bilayers [46]. Thus, the retention times are a good measure for the stability of the corresponding interface. This in turn reveals that the NMR interface is more stable in thick bilayers concluded from highest observed retention times (Figure 4). Worth mentioning, the results obtained by long scale simulations (Figure 4) and by atomistic simulations for 1.5 ns [47] are comparable with respect to e.g., helix tilt or crossing angle, which documents that properties of transmembrane helices mainly depend on the bilayer thickness.

Flanking amino acids provide a second contribution for the dimer stability and formation [48,49]. Indeed, this notion is consistent with the enhanced retention time for interfaces of the 40 mer in comparison to the 25 mer (Figure 4). Although it is known that the LIxxGxxxGxxxT motif of the GpA helix is sufficient for association in E. coli membranes [50,51], the flanking amino acids most likely modulate the dynamics of the assembly in vivo as well.

Thus, we conclude that the interface formation of GpA depends on the surrounding fatty acids, the amino acids flanking the transmembrane region and the helix conformation. In turn, it is important to reproduce a native environment, either in biological or in in silico experiments, to get information on the native behavior of proteins. This holds true especially for mechanisms which relay on the hydrophobic thickness of membranes like the ratio and transition between left- and right-handed interface e.g., of EphA2 [44] or EpHAl [43]. Translating our findings for the different membranes into a biological context suggests that modulation of dimerization interfaces of transmembrane proteins can be achieved by changes of the lipid composition, which is a common stresses reaction [52-54]. Such dependence of protein action and membrane dynamics has for example been described for Hik33 [55-57] or DesK [58]. Furthermore, there are several reports that the transport or enzyme activity of proteins depends on the bilayer thickness or curvature (e.g., reviewed in [59]). Remarkably, such lipid dependence has been reported for the hexose transporter band 4.5 [60], which is alike GpA localized in erythrocyte membranes. Thus, it is not unlikely that the function of GpA is modulated by membrane environment changes as well.

4. Methods

4.1. Simulation Details

Three different versions of the MARTINI force field were used to describe the protein structure: (i) the parameters used in Sengupta and Marrink 2010 (referred to as GpA2010; [26]); (ii) Martini2.1 [28]; (iii) Martini2.2 [29] parameters as created by the martinize.py script. The PC lipids with different fatty acids are standard MARTINI lipids [27]. All simulations were performed using GROMACS v4.5.5 [61] with standard MARTINI parameters: the neighbor lists are updated every 10th time step and a cutoff of 1.2 nm for the non-bonded interactions is applied (smooth switching from 0.0 to 1.2 nm for Coulomb potential and 0.9 to 1.2 nm for Lennard-Jones potential). A relative dielectric constant of 15 is used. For temperature (T = 320 K) and pressure coupling (semi-isotrope, p = 1 bar) the Berendsen thermostat and barostat are used [62]. All simulation times reported are raw simulation times and not multiplied by a special factor to consider the speed-up from coarse graining as used in other reports [63].

4.2. System Composition

The NMR structure (pdb:1afo—40 mer; [21]) was converted into a coarse grain model. The structure was truncated for simulations of the 25 amino acid long transmembrane helix (Glu72 to Arg96). For the self-assembly simulations helix A of the NMR structure was used. The stability of the dimer interface was tested by inserting the coarse grain structure based on NMR data [21], into a pre-equilibrated bilayer (protein:lipid:water ratio of 2:200:3000 (25 mer)/4500(40 mer)), followed by energy minimization and an MD run where the GpA dimer is position-restrained with a force constant

of 1000 kJmol-1nm-2. For the self-assembly the two single helices were inserted with a maximum distance from each other into a pre-equilibrated bilayer (protein:lipid:water ratio of 2:200:3000). The N- and C-terminus of the protein are uncharged and the simulation box has no net charge.

4.3. Details for the Short Self-Assembly Simulations

For short self-assembly simulations each inserted helix is rotated randomly around the helix axis to ensure that dimer formation is not influenced by the starting position. Each system is simulated for 0.3 |is and if no stable dimer is observed in the last frame continued for 0.1 |is. A dimer is considered as stable when at least 19 contacts occur between different amino acids of the two helices for a time span of 100 ns. 19 contacts is a suitable cutoff, because this number of contacts was derived from the simulation of the NMR dimer with the Martini2.2 force field (Figure S3). For the cluster analysis we used only the last frame of the short self-assembly trajectories.

4.4. Description of the Similarity to the NMR Structure and Clustering of Results

A contact map for each frame was calculated using the amino acids of the transmembrane part (25 mer). In the according figures, the X-axis indicates the residues of helix A and the 7-axis the residues of helix B. If residue i of helix A and residue j of helix B are in contact with each other with a distance shorter than 6 A the value for the position (ij) is set to one. If the distance is >6 A, zero is assigned to (ij). Distances are calculated with Yasara [64]. In the plot, the frequency of counts <6 A are given in grey scale.

For the comparison to the NMR structure, the contact map of each frame is compared to the map of the NMR dimer and a simple score representing the similarity is defined:

. . c(TT)

slmllarltyNMR = c(NMR) (1)

with c(TT) representing the number of contacts coinciding between NMR and MD structure, and with c(NMR) representing the total number of contacts between amino acids of the different helices of the NMR structure. In the transmembrane part of the NMR structure of GpA we observe 21 contacts (c(NMR) = 21).

For the clustering of all structures from the different self-assembly trajectories structures obtained every 50 ns were used. The contact maps modified to represent Boolean vectors are the basis for clustering. All frames that contain structures with less than 19 contacts were excluded from the clustering. Average linkage clustering with a cutoff of 0.5 for the cluster formation was used as implemented in SciPy [65]. A distance matrix for all contact maps is constructed using the dice dissimilarity [66] as implemented in SciPy [65].

4.5. Autocorrelation, Correlation and Fitting

For the (auto-) correlation plots we used Pearson's correlation coefficient (r):

Zn=1(x- x)(y- y)

- ~x)2J z:,(y, -y)

where is the number of randomly drawn simulations without replacement, and are frequencies of cluster i and are the averages. For testing the autocorrelation for n = 150, the first 150 interfaces were used for the distribution X and the remaining 150 interfaces for the distribution Y. In total the input orders were permuted 1000 times.

Afterwards a double exponential function is fitted onto the average correlation values using the least-squares fitting as implemented in SciPy [65]. As target function we use f = (y - eFunc (x))/ std(y) where y is the observed value, std(y)is the standard deviation of y and eFunc(x) is the calculated value with eFunc(x) = a(1 - exp(-bx)) + c(1 - exp(-dx)). The parameters a, b, c and d are estimated from the data and afterwards optimized using least-squares fitting (Table 1).

4.6. Construction of the Transition Net

The transition net is built from the long self-assembly simulations. Only clusters with a frequency, which at least corresponds to the NMR cluster, are considered for the network (clusters 12, 14, 16, 42, 43, 44, 52, 63, 76 and 81). Each cluster forms a node, which is connected via an edge if a direct or indirect transition between the two clusters is observed in the trajectories. A transition is regarded as indirect if two clusters are connected via other clusters, which are not explicitly represented in the net. Shown are the three most frequent transitions for each cluster. Transitions occurring only once are not considered, because they are not reproduced in our data.

4.7. Backmapping and Comparison to NOE's

The back mapping of the coarse grain structures to fine grain structures was performed using the back mapping tool BACKWARD [67]. In short, we used the charmm36 force field and did some energy minimizations and MD simulations with constant pressure using different time steps. The resulting structures in full atomistic resolution were compared to the maximal intermolecular distances obtained from the NMR experiments.

5. Conclusions

The self-assembly of a transmembrane segment depends on the helix length, the conformation of the helix and on the nature of the surrounding fatty acids. The latter influence the hydrophobic thickness and the fluidity of the membrane. By this property the nature of the fatty acid determines the interface formed shortly after assembly, the transitions between interfaces after dimer formation and the interface stability. The structure of the dimer interface depends to a certain extent on the force field parameter, because application of an elastic network was found to enhance the detection of the interface established by NMR experiments. Taken together, our results document on the one hand the importance of further force field development of protein parameters for the analysis of membrane domain assembly and on the other hand the requirement for correct modeling of membranes to obtain detailed information on the behavior of proteins in their native environment.

Acknowledgments

This work would not have been possible without computation time at the Center for Scientific Computing (CSC) of the Goethe University Frankfurt. This work was supported by grants from the Deutsche Forschungsgemeinschaft: CEF "Macromolecular Complexes" in the frame of the FOCUS-Project to Enrico Schleiff; the CRC of SFB 807 P17 to Enrico Schleiff/Nadine Flinner and the Center of Membrane Proteomics to Nadine Flinner is acknowledged.

Author Contributions

Enrico Schleiff, Nadine Flinner and Oliver Mirus conceptualized the project and Enrico Schleiff headed the project. Nadine Flinner performed experiments. Nadine Flinner and Enrico Schleiff were involved in writing the manuscript.

Conflicts of Interest

The authors declare no conflict of interest. References

1. Ohvo-Rekilä, H.; Ramstedt, B.; Leppimäki, P.; Slotte, J.P. Cholesterol interactions with phospholipids in membranes. Prog. Lipid Res. 2002, 41, 66-97.

2. Ikonen, E. Cellular cholesterol trafficking and compartmentalization. Nat. Rev. Mol. Cell Biol. 2008, 9, 125-138.

3. Van Meer, G.; de Kroon, A.I.P.M. Lipid map of the mammalian cell. J. Cell Sci. 2010, 124, 5-8.

4. Risselada, H.J.; Marrink, S.J. The molecular face of lipid rafts in model membranes. Proc. Natl. Acad. Sci. USA 2008, 105, 17367-17372.

5. Davis, R.S.; Sunil Kumar, P.B.; Sperotto, M.M.; Laradji, M. Predictions of phase separation in three-component lipid membranes by the MARTINI force field. J. Phys. Chem. B 2013, 117, 4072-4080.

6. Domanski, J.; Marrink, S.J.; Schäfer, L.V. Transmembrane helices can induce domain formation in crowded model membranes. Biochim. Biophys. Acta 2012, 1818, 984-994.

7. Lee, A. Membrane structure. Curr. Biol. 2001, 11, R811-R814.

8. Yeagle, P.L. Lipid regulation of cell membrane structure and function. FASEB J. 1989, 3, 1833-1842.

9. Engelman, D.M. Membranes are more mosaic than fluid. Nature 2005, 438, 578-580.

10. Lee, A.G. Lipid-protein interactions in biological membranes: A structural perspective. Biochim. Biophys. Acta 2003, 1612, 1-40.

11. Schäfer, L.V.; de Jong, D.H.; Holt, A.; Rzepiela, A.J.; de Vries, A.H.; Poolman, B.; Killian, J.A.; Marrink, S.J. Lipid packing drives the segregation of transmembrane helices into disordered lipid domains in model membranes. Proc. Natl. Acad. Sci. USA 2011, 108, 1343-1348.

12. De Jong, D.H.; Lopez, C.A.; Marrink, S.J. Molecular view on protein sorting into liquid-ordered membrane domains mediated by gangliosides and lipid anchors. Faraday Discuss. 2013, 161, 347-363; discussion 419-459.

13. Phillips, R.; Ursell, T.; Wiggins, P.; Sens, P. Emerging roles for lipids in shaping membrane-protein function. Nature 2009, 459, 379-385.

14. Holt, A.; Killian, J.A. Orientation and dynamics of transmembrane peptides: The power of simple models. Eur. Biophys. J. 2010, 39, 609-621.

15. Strandberg, E.; Esteban-Martín, S.; Ulrich, A.S.; Salgado, J. Hydrophobic mismatch of mobile transmembrane helices: Merging theory and experiments. Biochim. Biophys. Acta 2012, 1818, 1242-1249.

16. Parton, D.L.; Klingelhoefer, J.W.; Sansom, M.S.P. Aggregation of model membrane proteins, modulated by hydrophobic mismatch, membrane curvature, and protein class. Biophys. J. 2011, 101, 691-699.

17. Leidl, K.; Liebisch, G.; Richter, D.; Schmitz, G. Mass spectrometric analysis of lipid species of human circulating blood cells. Biochim. Biophys. Acta 2008, 1781, 655-664.

18. De Oliveira, S.; Saldanha, C. An overview about erythrocyte membrane. Clin. Hemorheol. Microcirc.

2010, 44, 63-74.

19. Reid, M.E. MNS blood group system: A review. Immunohematology 2009, 25, 95-101.

20. Mankelow, T.; Satchwell, T.; Burton, N. Refined views of multi-protein complexes in the erythrocyte membrane. Blood Cells Mol. Dis. 2012, 49, 1-10.

21. MacKenzie, K.R.; Prestegard, J.H.; Engelman, D.M. A transmembrane helix dimer: Structure and implications. Science 1997, 276, 131-133.

22. Mineev, K.S.; Bocharov, E.V.; Volynsky, P.E.; Goncharuk, M.V.; Tkach, E.N.; Ermolyuk, Y.S.; Schulga, A.A.; Chupin, V.V.; Maslennikov, I.V.; Efremov, R.G.; et al. Dimeric structure of the transmembrane domain of glycophorin a in lipidic and detergent environments. Acta Naturae

2011, 3, 90-98.

23. Russ, W.P.; Engelman, D.M. The GxxxG motif: A framework for transmembrane helix-helix association. J. Mol. Biol. 2000, 296, 911-919.

24. Bond, P.J.; Sansom, M.S.P. Insertion and assembly of membrane proteins via simulation. J. Am. Chem. Soc. 2006, 128, 2697-2704.

25. Psachoulia, E.; Fowler, P.W.; Bond, P.J.; Sansom, M.S.P. Helix-helix interactions in membrane proteins: Coarse-grained simulations of glycophorin a helix dimerization. Biochemistry 2008, 47, 10503-10512.

26. Sengupta, D.; Marrink, S.J. Lipid-mediated interactions tune the association of glycophorin A helix and its disruptive mutants in membranes. Phys. Chem. Chem. Phys. 2010, 12, 12987-12996.

27. Marrink, S.J.; de Vries, A.H.; Mark, A.E. Coarse grained model for semiquantitative lipid simulations. J. Phys. Chem. B 2004, 108, 750-760.

28. Monticelli, L.; Kandasamy, S.K.; Periole, X.; Larson, R.G.; Tieleman, D.P.; Marrink, S.-J. The MARTINI coarse-grained force field: Extension to proteins. J. Chem. Theory Comput. 2008, 4, 819-834.

29. De Jong, D.H.; Singh, G.; Bennett, W.F.D.; Arnarez, C.; Wassenaar, T.A.; Schäfer, L.V.; Periole, X.; Tieleman, D.P.; Marrink, S.J. Improved parameters for the martini coarse-grained protein force field. J. Chem. Theory Comput. 2013, 9, 687-697.

30. Long, F.; Du, L.-C.; Mei, D.-C. Asymmetric effects on the associated relaxation time and the correlation function of a bistable system with correlated noises. Phys. Scr. 2009, 79, 45007. doi:10.1088/0031-8949/79/04/045007.

31. Gaspard, P. The correlation time of mesoscopic chemical clocks. J. Chem. Phys. 2002, 117, 8905. doi:10.1063/1.1513461.

32. Periole, X.; Cavalli, M.; Marrink, S.-J.; Ceruso, M.A. Combining an elastic network with a coarse-grained molecular force field: Structure, dynamics, and intermolecular recognition. J. Chem. Theory Comput. 2009, 5, 2531-2543.

33. Smith, S.O.; Eilers, M.; Song, D.; Crocker, E.; Ying, W.; Groesbeek, M.; Metz, G.; Ziliox, M.; Aimoto, S. Implications of threonine hydrogen bonding in the glycophorin A transmembrane helix dimer. Biophys. J. 2002, 82, 2476-2486.

34. Marrink, S.J.; Tieleman, D.P. Perspective on the Martini model. Chem. Soc. Rev. 2013, 42, 6801-6822.

35. Anbazhagan, V.; Schneider, D. The membrane environment modulates self-association of the human GpA TM domain—Implications for membrane protein folding and transmembrane signaling. Biochim. Biophys. Acta 2010, 1798, 1899-1907.

36. Dony, N.; Crowet, J.M.; Joris, B.; Brasseur, R.; Lins, L. SAHBNET, an accessible surface-based elastic network: An application to membrane protein. Int. J. Mol. Sci. 2013, 14, 11510-11526.

37. Globisch, C.; Krishnamani, V.; Deserno, M.; Peter, C. Optimization of an elastic network augmented coarse grained model to study CCMV capsid deformation. PLoS One 2013, 8, e60582.

38. North, B.; Cristian, L.; Fu Stowell, X.; Lear, J.D.; Saven, J.G.; Degrado, W.F. Characterization of a membrane protein folding motif, the Ser zipper, using designed peptides. J. Mol. Biol. 2006, 359, 930-939.

39. Cuthbertson, J.M.; Bond, P.J.; Sansom, M.S.P. Transmembrane helix-helix interactions: Comparative simulations of the glycophorin a dimer. Biochemistry 2006, 45, 14298-14310.

40. Polyansky, A.A.; Volynsky, P.E.; Efremov, R.G. Multistate organization of transmembrane helical protein dimers governed by the host membrane. J. Am. Chem. Soc. 2012, 134, 14390-14400.

41. Doura, A.K.; Kobus, F.J.; Dubrovsky, L.; Hibbard, E.; Fleming, K.G. Sequence context modulates the stability of a GxxxG-mediated transmembrane helix-helix dimer. J Mol Biol. 2004, 341, 991-998.

42. Manni, S.; Mineev, K.S.; Usmanova, D.; Lyukmanova, E.N.; Shulepko, M.A.; Kirpichnikov, M.P.; Winter, J.; Matkovic, M.; Deupi, X.; Arseniev, A.S.; et al. Structural and Functional Characterization of Alternative Transmembrane Domain Conformations in VEGF Receptor 2 Activation. Structure 2014, 22, 1077-1089.

43. Bocharov, E.V.; Mayzel, M.L.; Volynsky, P.E.; Goncharuk, M.V.; Ermolyuk, Y.S.; Schulga, A.A.; Artemenko, E.O.; Efremov, R.G.; Arseniev, A.S. Spatial structure and pH-dependent conformational diversity of dimeric transmembrane domain of the receptor tyrosine kinase EphA1. J. Biol. Chem. 2008, 283, 29385-29395.

44. Bocharov, E.V.; Mayzel, M.L.; Volynsky, P.E.; Mineev, K.S.; Tkach, E.N.; Ermolyuk, Y.S.; Schulga, A.A.; Efremov, R.G.; Arseniev, A.S. Left-handed dimer of EphA2 transmembrane domain: Helix packing diversity among receptor tyrosine kinases. Biophys. J. 2010, 98, 881-889.

45. Tang, T.C.; Hu, Y.; Kienlen-Campard, P.; El Haylani, L.; Decock, M.; van Hees, J.; Fu, Z.; Octave, J.N.; Constantinescu, S.N.; Smith, S.O. Conformational changes induced by the A21G Flemish mutation in the amyloid precursor protein lead to increased Ap production. Structure 2014, 22, 387-396.

46. Janosi, L.; Prakash, A.; Doxastakis, M. Lipid-modulated sequence-specific association of glycophorin A in membranes. Biophys. J. 2010, 99, 284-292.

47. Petrache, H.I.; Grossfield, A.; MacKenzie, K.R.; Engelman, D.M.; Woolf, T.B. Modulation of glycophorin A transmembrane helix interactions by lipid bilayers: Molecular dynamics calculations. J. Mol. Biol. 2000, 302, 727-746.

48. Orzáez, M.; Pérez-Payá, E.; Mingarro, I. Influence of the C-terminus of the glycophorin A transmembrane fragment on the dimerization process. Protein Sci. 2000, 9, 1246-1253.

49. Lemmon, M.A.; Flanagan, J.M.; Hunt, J.F.; Adair, B.D.; Bormann, B.J.; Dempsey, C.E.; Engelman, D.M. Glycophorin A dimerization is driven by specific interactions between transmembrane alpha-helices. J. Biol. Chem. 1992, 267, 7683-7689.

50. Gerber, D.; Shai, Y. In vivo detection of hetero-association of glycophorin-A and its mutants within the membrane. J. Biol. Chem. 2001, 276, 31229-31232.

51. Brosig, B.; Langosch, D. The dimerization motif of the glycophorin A transmembrane segment in membranes: Importance of glycine residues. Protein Sci. 1998, 7, 1052-1056.

52. Lindberg, L.; Santos, A.X.; Riezman, H.; Olsson, L.; Bettiga, M. Lipidomic profiling of Saccharomyces cerevisiae and Zygosaccharomyces bailii reveals critical changes in lipid composition in response to acetic acid stress. PLoS One 2013, 8, e73936.

53. Ghorbal, S.K.B.; Chatti, A.; Sethom, M.M.; Maalej, L.; Mihoub, M.; Kefacha, S.; Feki, M.; Landoulsi, A.; Hassen, A. Changes in membrane fatty acid composition of Pseudomonas aeruginosa in response to UV-C radiations. Curr. Microbiol. 2013, 67, 112-117.

54. Oberg, T.S.; Ward, R.E.; Steele, J.L.; Broadbent, J.R. Genetic and physiological responses of Bifidobacterium animalis subsp. lactis to hydrogen peroxide stress. J. Bacteriol. 2013, 195, 3743-3751.

55. Shimura, Y.; Shiraiwa, Y.; Suzuki, I. Characterization of the subdomains in the ^-terminal region of histidine kinase Hik33 in the cyanobacterium Synechocystis sp. PCC 6803. Plant Cell Physiol. 2012, 53, 1255-1266.

56. Mikami, K.; Kanesaki, Y.; Suzuki, I.; Murata, N. The histidine kinase Hik33 perceives osmotic stress and cold stress in Synechocystis sp PCC 6803. Mol. Microbiol. 2002, 46, 905-915.

57. Mironov, K.S.; Sidorov, R.A.; Trofimova, M.S.; Bedbenov, V.S.; Tsydendambaev, V.D.; Allakhverdiev, S.I.; Los, D.A. Light-dependent cold-induced fatty acid unsaturation, changes in membrane fluidity, and alterations in gene expression in Synechocystis. Biochim. Biophys. Acta 2012, 1817, 1352-1359.

58. Cybulski, L.E.; Martín, M.; Mansilla, M.C.; Fernández, A.; de Mendoza, D. Membrane thickness cue for cold sensing in a bacterium. Curr. Biol. 2010, 20, 1539-1544.

59. Andersen, O.S.; Koeppe, R.E. Bilayer thickness and membrane protein function: An energetic perspective. Annu. Rev. Biophys. Biomol. Struct. 2007, 36, 107-130.

60. Carruthers, A.; Melchior, D.L. Human erythrocyte hexose transporter activity is governed by bilayer lipid composition in reconstituted vesicles. Biochemistry 1984, 23, 6901-6911.

61. Hess, B.; Kutzner, C.; van der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for Highly Efficient, Load-Balanced, and Scalable Molecular Simulation. J. Chem. Theory Comput. 2008, 4, 435-447.

62. Berendsen, H.J.C.; Postma, J.P.M.; van Gunsteren, W.F.; DiNola, A.; Haak, J.R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684. doi:10.1063/1.448118.

63. Marrink, S.J.; Risselada, H.J.; Yefimov, S.; Tieleman, D.P.; de Vries, AH. The MARTINI force field: Coarse grained model for biomolecular simulations. J. Phys. Chem. B 2007, 111, 7812-7824.

64. Krieger, E.; Vriend, G. Models@Home: Distributed computing in bioinformatics using a screensaver based approach. Bioinformatics 2002, 18, 315-318.

65. Jones, E.; Oliphant, T.; Peterson, P. SciPy: Open Source Scientific Tools for Python. 2001. Available online: http://www.scipy.org/ (access on 22 May 2014).

66. Dice, L.R. Measures of the Amount of Ecologic Association between Species. Ecology 1945, 26, 297-302.

67. Wassenaar, T.A.; Pluhackova, K.; Böckmann, R.A.; Marrink, S.J.; Tieleman, D.P. Going Backward: A Flexible Geometric Approach to Reverse Transformation from Coarse Grained to Atomistic Models. J. Chem. Theory Comput. 2014, 10, 676-690.

© 2014 by the authors; licensee MDPI, Basel, Switzerland. This article is an open access article

distributed under the terms and conditions of the Creative Commons Attribution license

(http://creativecommons.org/licenses/by/3.0/).

Supplementary Information

Table S1. The 5 minimal distances between beads from different monomers of homodimeric bitopic a-helices with known NMR structure, which was converted into a CG representation for the analysis.

2ka2 2k9y 1afo

helix A helix B distance helix A helix B distance helix A helix B distance

BBPHE BBSER 3.10488 BBLEU BBALA 3.2457 BBGLU BBPRO 2.74635

BBPHE BBSER 3.18935 BBLEU BBALA 3.25337 BBILE BBLYS 2.86863

BBLEU BBPRO 3.2447 BBPHE BBILE 3.25935 BBGLU BBPRO 3.29006

SC1SER SC2HIS 3.31513 BBPHE BBILE 3.27329 BBILE BBLYS 3.37511

BBGLY SC1SER 3.31833 SC1SER SC1LEU 3.56146 BBGLU SC1PRO 3.49815

2k1k 2ka1 2l2t

helix A helix B distance helix A helix B distance helix A helix B distance

BBLEU BBTHR 3.64333 BBPHE BBSER 3.18135 BBARG BBLYS 3.15629

BBLEU BBTHR 3.64489 BBLEU BBPRO 3.2233 BBARG BBLYS 3.16398

BBARG BBSER 3.67019 BBLEU BBPRO 3.25696 BBARG BBTHR 3.6572

BBARG BBSER 3.67661 BBPHE BBSER 3.31495 BBARG BBTHR 3.72581

BBVAL BBARG 4.04698 BBARG BBLEU 3.35745 SC3HIS SC1GLN 3.78445

2l9u 2loh 2l6w

helix A helix B distance helix A helix B distance helix A helix B distance

SC3PHE SC3PHE 3.13432 BBSER BBGLN 3.08917 BBGLN BBLYS 3.12325

BBTHR BBHIS 3.20186 BBASN BBLYS 3.12417 BBGLN BBLYS 3.28172

BBTHR BBHIS 3.21512 BBVAL BBGLY 3.3569 BBPRO BBPHE 3.49944

SC2PHE SC2PHE 4.06314 BBGLY BBSER 3.3671 BBPRO BBPHE 3.54458

SC3PHE SC2PHE 4.11281 BBVAL BBGLY 3.3774 BBGLN SC 1 LYS 4.2117

2j5d 2l34

helix A helix B distance helix A helix B distance

BBPHE BBLEU 3.19401 BBARG BBLEU 2.67973

BBPHE BBSER 3.32284 BBARG BBLEU 2.71681

BBLEU BBTHR 3.37156 BBLEU BBGLY 3.22939

SC1SER SC2HIS 3.38521 BBSER BBPRO 3.40705

BBLEU BBTHR 3.45504 BBSER SC1PRO 3.52303

Figure S1. Representative Structures for the most frequent clusters. Shown are the contact plot for the whole cluster (as shown in Figure 1), the contact plot of the representative structure and the CG-structure. In the CG-structure all amino acids which are typical interface for the interface are colored in red, other residues of molecule A/molecule B are colored in light grey/dark grey. All backbone beads are show as small ball and side chain beads are shown as bigger balls. For celerity the distance between the Glycine residues of the GxxxG motif is indicated for the NMR like structures.

Figure S2. Frequency of Transitions between frequent clusters. Each of the eleven plots displays the transitions from one cluster to all other frequent clusters listed on the X-axis for the three different lipids (7-axis). The darker the color the more frequent the transition from cluster x to y.

Figure S3. Number of contacts between amino acids. Number of contacts between amino acids of different helices (X-axis) is plotted against the frequency (7-axis) observed in the simulations of the GpA_TM NMR structure (25 mer) with the Martini 2.2 force field and different fatty acids. The red, dashed line corresponds to the cutoff used for the short self-assembly simulations: The mean of observed contacts per frame is 31.79 and the standard derivation 4.01 ([mean + 3a J = [19.76 J = 19).

£ H—

£ 150

100 50

°0 5 10 15 20 25 30 35 40 45

contacts

Copyright of International Journal of Molecular Sciences is the property of MDPI Publishing and its content may not be copied or emailed to multiple sites or posted to a listserv without the copyright holder's express written permission. However, users may print, download, or email articles for individual use.