Lipid-Dependent Alternating Access Mechanism of a Bacterial Multidrug ABC Exporter.

By undergoing conformational changes, active membrane transporters alternate between an inward-facing (IF) and an outward-facing (OF) state to transport their substrates across cellular membrane. The conformational landscape of membrane transporters, however, could be influenced by their environment, and the dependence of the alternating access mechanism on the lipid composition has not been understood at the molecular level. We have performed an extensive set of microsecond-level all-atom molecular dynamics (MD) simulations on bacterial ATP binding cassette (ABC) exporter Sav1866 in six different phosphocholine (PC) and phosphoethanolamine (PE) lipid membrane environments. This study mainly focuses on the energetically downhill OF-to-IF conformational transition of Sav1866 upon the ATP hydrolysis. We observe that the transporter undergoes large-scale conformational changes in the PE environment, particularly in the POPE lipids, resulting in an IF-occluded conformation, a transition that does not occur when the transporter is embedded in any of the PC lipid bilayers. We propose that the PE lipids facilitate the closing of the protein on the periplasmic side due to their highly polar headgroups that mediate the interaction of the two transmembrane (TM) bundles by a network of lipid-lipid and lipid-protein hydrogen bonds. POPE lipids in particular facilitate the closure of periplasmic gate by promoting a hinge formation in TM helices and an interbundle salt bridge formation. This study explains how the alternating access mechanism and the flippase activity in ABC exporters could be lipid-dependent.

Introduction

ABC transporters are a group of integral membrane proteins that can harness the energy stored in ATP to actively transport a variety of substrates including ions, amino acids, peptides, and therapeutic drugs across the membrane.[1] Members of the ABC transporter superfamily are found in both prokaryotes and eukaryotes.[2] Most of the 48 ABC transporters identified in humans are exporters,[3] and some such as P-glycoprotein (P-gp), ABCG2, and ABCC1 are also multidrug resistance (MDR) transporters.[4−6] In bacteria, numerous ABC exporters are responsible for active extrusion of antibiotics. In humans, transporters such as P-gp (also known as MDR1 and ABCB1) participate in tumor resistance to chemical treatment when overexpressed.[7,8] ABC exporters are either homo- or heterodimers[9−11] and are composed of two transmembrane domains (TMDs) that form the substrate translocation pathway, and two nucleotide binding domains (NBDs) that bind to the ATP.

Membrane transporters undergo large-scale conformational changes to function, and the widely accepted model that describes these conformational changes is the alternating access mechanism.[12] In line with this mechanism, an ABC transporter alternates between an ATP-bound OF and a nucleotide-free IF conformation (i.e., the resting state) in the transport cycle, while the occluded conformation is a necessary intermediate.[13−17] For exporters, the substrate enters into the translocation chamber when the transporter is in the IF state, while binding of ATPs to the NBDs triggers the IF → OF transition. This conformational transition occurs through an occluded intermediate to prevent the passive diffusion of substrate along its concentration gradient from the periplasmic/extracellular side to the cytoplasm. ATP hydrolysis and the release of ADPs from NBDs resets the transporter back to the IF state.[18]

Lipid composition is known to significantly influence the structure and function of membrane transporters[19−22] including those in the ABC superfamily.[23−25] In addition, some transporters (including various ABC transporters[26−28]) act as flippases by selectively translocating lipids from one bilayer leaflet to the other.[29] Several studies have reported on the lipid dependence of the function of various ABC transporters such as ABCR,[30] MsbA,[31] HorA,[23] ABCA1,[32] MalFGK2,[24] and TmrAB.[25] However, a detailed molecular-level study of how different lipids influence the activity of ABC transporters is missing. In this study, we have probed the lipid-dependent nature of alternating access mechanism in a bacterial multidrug ABC exporter, namely, Sav1866, using microsecond-level all-atom MD simulations. Sav1866 is a homodimer, consisting of two identical TMDs and two identical NBDs. The two available crystal structures of Sav1866 are both in the OF state.[33,34] We have modeled Sav1866 based on the one with the higher resolution of 3 Å (PDB, 2HYD; Figure A) for this study.[33]

Figure 1

Conformational dynamics of Sav1866 in POPE and POPC lipid environments. (A) Crystal structure of Sav1866 in the OF state: the monomers A and B are colored green and orange, respectively, and the two ATP molecules are colored in magenta. β time series in POPE (B) and POPC (D) simulations along with the β distribution in the last half of the simulations (i.e., from 1.2 to 2.4 μs) for POPE (C) and POPC (E) simulations. (F) Projection of POPE and POPC trajectories onto the (PC1, PC2) space, i.e., the two first principal components of protein Cα atoms. Snapshots of the transporter conformation before and after the equilibration in the POPE (G) and POPC (H) bilayers. Water density in and around the transporter associated with a 10 Å thick cross section of simulation box during the first and last 0.5 μs of POPE1 (I) and POPC1 (J) simulations.

The MD simulation technique is a powerful tool to study the structural dynamics of transporters at an atomic level.[35−37] Many MD simulations of ABC transporters have been conducted including Sav1866,[38−42] MsbA,[43] P-gp,[44−47] ABCB10,[48] McjD,[49,50] TM287/288,[51] and TAP;[52] however, they have been mostly restricted to nanosecond-level time scales. Studying the conformational changes of membrane transporters via all-atom MD is computationally demanding as such conformational changes occur on the order of microseconds to seconds. Very recently, a study by Goddeke et al.[35] has used microsecond-level all-atom MD simulations of a heterodimeric ABC exporter TM287/288; however, this study differs from ours in two major aspects: (1) only one lipid type (POPC) is used in the aforementioned study, and (2) the study focuses on the IF → OF transition rather than the OF → IF transition, which is the subject of our study. Therefore, our study provides a complementary picture, which particularly focuses on the lipid-dependent behavior of ABC transporters.

To our knowledge, this is the first comprehensive computational study that gives insight into the influence of lipid environment on the alternating access mechanism of ABC exporters. Based on an extensive set of microsecond-level Sav1866 simulations performed here we observe that the apo protein reproducibly adopts an IF-occluded (IFocc) conformation in the POPE (and some of the other phosphoethanolamine, PE) environments, but not in any of the PC environments. With the wealth of information provided by these simulations (approximately about 50 μs), we explain why the PE lipids favor the occlusion of the nucleotide-free transporter on the periplasmic side, and how these lipids promote the flippase activity of Sav1866. We note that, here, we have only studied the energetically downhill OF → IF conformational transition of Sav1866 upon ATP hydrolysis. The energetically uphill conformational transitions of ABC exporters, which are associated with much longer time scales, can be studied using enhanced sampling techniques.[53−56]

Overall, our computational findings are in line with and could help understand various experimental studies that, due to the limitations of experimental techniques, could not be fully understood. For instance, we observe that the Sav1866 transporter shows a higher dynamic activity in the PE membrane environment as compared to PC, consistent with the conclusions made based on experimental studies of the photoreceptor ABC transporter ABCR.[30] Our observations are also in line with experimental studies of ABC importer MolBC that shows a lipid-specific trend in its activity and conformational changes.[57] More specifically, our results explain the loss of substrate transport in Sav1866 homologue HorA in the PC environment as shown by Gustot et al., who clearly show that while the ATPase activity is boosted in the PC environment (as compared to PE), the transport activity is lost.[23] Employing infrared spectroscopy, Gustot et al. also showed that the transmembrane helices of HorA adopt substantially different orientations in the PE and PC environments,[23] which is also consistent with our findings. Our study provides an explanation for the experimental observations mentioned above and more generally provides a nuanced picture of the lipid-dependent behavior of membrane transporters. It also calls for a more careful interpretation of studies conducted in non-native environments for membrane transporters. For instance, bacterial and mammalian membrane proteins that function in membranes of substantially different compositions may not show a similar behavior when reconstituted in the same membrane as shown by Moeller et al. for bacterial MsbA and mammalian P-gp, both of which are ABC exporters but display a divergence in their conformational landscapes when embedded in the same membrane-mimicking environment.[58]

Results and Discussion

The nucleotide-free apo model of Sav1866 was generated by removing the nucleotides from the crystal structure of the nucleotide-bound OF state of Sav1866 (PDB, 2HYD; Figure A)[33] and embedded in six different lipid environments, namely, POPE, DOPE, DPPE, POPC, DOPC, and DPPC lipid bilayers. We performed three independent sets of 2.4 μs long MD simulations based on each model, along with one negative control simulation with the Mg-ATP bound OF state embedded in the POPE lipid bilayer, which shows no major conformational changes as expected. We observed that the ATP-removed OF state behaves strikingly differently in the PE and PC lipid environments. To simplify the discussion, we start with the POPE and POPC simulations that consistently show a substantially different behavior in all three repeats. We will then extend our discussion to other types of PE and PC lipids. In addition to protein conformational changes, we have also observed indications for the flippase activity of Sav1866 in some of the PE simulations that will be discussed prior to our concluding remarks. The methodological details can be found in the Methods section.

Sav1866 Adopts an IFocc Conformation in the POPE but not POPC Environment

To monitor the global protein conformational changes associated with the opening/closing of the periplasmic side of the protein, we measured an angle between the two bundles of TM helices denoted by β (see the Methods section) (Figure A). Roughly speaking, if β is less than 15°, the TMDs may be considered closed on the periplasmic side. In all three simulations of apo protein in the POPE lipids (Figure B), the transporter closes on the periplasmic side within the first half of the simulation (i.e., within 1.2 μs) and remains closed for the rest of the simulation. Comparing to the negative control (the Mg-ATP bound state in the POPE lipids, Figure B), there is clear indication that the closing of the periplasmic side of the TMDs in the POPE environment is due to the removal of the ATPs (mimicking the ATP hydrolysis and ADP release). On the other hand, in none of the POPC simulations (Figure D), the transporter adopts a closed state on the periplasmic side.

The β distributions based on the last 1.2 μs of the simulations (Figure C,E) also reveal a clear difference in TMD dynamics between the POPE and POPC environments. The periplasmic gate stays open throughout the 3  ×  2.4 μs aggregate simulations of Sav1866 in the POPC environment. This observation about the dynamics of a bacterial ABC exporter in the POPC environment is consistent with that made by Goddeke et al. with regards to TM287/288 that does not undergo any major conformational changes in any of the 30 500 ns long MD simulations carried in the presence of the POPC lipid bilayer, when the ATPs were removed.[35] However, since the same starting point in the POPE environment consistently results in a conformational transition, we can conclude that these conformational changes are lipid-dependent. It is also interesting to note that Goddeke et al. have observed major conformational changes occurring in the POPC environment in the reverse direction (at least in some of the simulations) when the initial point is the ATP-docked IF state of TM287/288. We have not studied the reverse direction here since the IF state of Sav1866 is not available. We may speculate the reverse reaction for the ATP-docked Sav1866 can occur in the POPC environment as well; however, it is not clear whether or not this would occur at a rate similar to that in the POPE environment. What is clear from our simulation data is that the behavior of the protein in the two lipid environments is substantially different.

It is important to note that although the protein stays open to the periplasm in the POPC environment, other conformational changes take place even in this environment. This can be seen by monitoring the root mean square deviation (RMSD) of the protein embedded in both environments with respect to the crystal structure (Figure S1A,B). This analysis shows a consistent increase in the RMSD of apo protein in all POPE (∼5–6 Å) and POPC (∼4–5 Å) simulations as compared to that of the Mg-ATP bound protein, i.e., the negative control (RMSD ≈ 3 Å). To compare the global protein conformational changes in the two environments, we have also performed principal component analysis (PCA) using a combination of all six POPE and POPC apo Sav1866 trajectories. When projected onto the first two principal components (PC1 and PC2), it is evident that the transporter occupies two well-separated regions of the conformational landscape in the POPE and POPC environments (Figure F). The contribution of these two PCs to the variance is ∼41% (Figure S2A), indicating that PC1 and PC2 capture a significant portion of the conformational changes.

We note that the occluded conformation observed in the POPE simulations is an intermediate and not the end state. The two NBDs are expected to completely dissociate upon the ATP hydrolysis and ADP release. One may assume that the heat released from the ATP hydrolysis may facilitate the process of NBD dissociation; however, even without this heat, the NBDs will eventually dissociate, leading to the complete opening of the periplasmic gate of Sav1866. In other words, since the IF state is thermodynamically favored in the absence of the nucleotides, the transition to the IF state is expected to occur even without the help of the heat from ATP hydrolysis. The initial stages of the NBD dissociation process can be already seen in our nucleotide-free simulations, regardless of the lipid types (made evident from an increase in the distance between the two NBDs; Figure S3). Since the transporter adopts a closed state on the periplasmic side in the POPE environment and due to the slight opening of the cytoplasmic side, we conclude that the transporter consistently adopts an IFocc conformation in the POPE environment within microseconds (Figure G). The complete transition to the IF state may occur on a longer time scale; however, the ATP hydrolysis can speed up the process. On the other hand, the protein does not make a transition to the IFocc state in the POPC environment within microseconds (Figure H). The distinct behavior of Sav1866 in POPE and POPC lipids is certainly indicative of substantially different kinetics. One may expect that the difference in the behavior of Sav1866 in different environments results in different thermodynamic properties as well.

A strong support for our claim comes out of the experimental observations made by Gustot et al.,[23] as discussed in the Introduction. The loss of transport activity in the PC environment while preserving the ATPase activity and the differential orientation of TM helices observed in the PE and PC environments for Sav1866 homologue HorA are indicative of the lipid-dependent behavior at the thermodynamic level. We note that while our computational results qualitatively agree with the experimental data reported by Gustot et al.[23] in terms of the differential orientation of TM helices, the comparison is not straightforward. This is due to the fact that our simulations are not long enough to allow the protein reach its resting state. In the POPE lipids, we expect this state to be the IF conformation, potentially with a large opening on the cytoplasmic side, similar to those observed in other ABC exporters.[59−62] In the POPC lipids, however, there is no indication of any TM conformational change toward the IF state in our simulations, which is in line with the experimentally observed loss of transport function in PC-reconstituted HorA.[23] We speculate that the distinct TM orientations observed experimentally for HorA in PE and PC lipids could be due to the difference in the functionality of the resting state in the two environments, where the apo protein adopts a wide open IF state in the PE environment and an OF state (or a combination of OF and other states) in the PC environment.

Our computational data on Sav1866 may also be used to shed light on biophysical data on ABC exporters such as double electron–electron resonance (DEER) measurements of BmrCD.[63] Although these measurements seem to be able to capture the protein in different functional states even in non-native environments, e.g., in detergent micelles or in non-PE-containing lipids, a closer look at the data shows that the extracted DEER distance distributions associated with the TMDs of BmrCD seem to be quite similar in the presence and absence of different nucleotides in these environments, an observation that may be explained by the decoupling of ATP hydrolysis and TMD conformational changes in certain environments as observed by Gustot et al. for HorA.[23] One exception is when vanadate is added along with the ATP to trap the high-energy transition state of ATP hydrolysis.[63] The focus of our work, however, is not on the ATP hydrolysis, and reconciling our simulations with these experimental data may require involving ATP hydrolysis in our simulations. Finally, we note that our results are not consistent with computational observations made by Xu et al.,[42] who report the closure of the periplasmic gate of Sav1866 in a very short simulation (100 ns) in the POPC environment upon the removal of the nucleotides. Given the fact that both studies use the exact same force field and crystal structure (but different temperatures; i.e., 310 K here and 323.15 K in the study by Xu et al.[42]), it is surprising to see that the periplasmic closure has been observed in the only simulation of apo protein performed by Xu et al.,[42] which is only 100 ns, while none of the three 2.4 μs simulations performed by us in the same environment (POPC) has resulted in the closure of the periplasmic gate. The difference in the temperature is very unlikely to explain such a significant difference, particularly when we consider the computational observations of Goddeke et al.[35] that use a temperature of 375 K to study ABC exporter TM287/288. They show that, upon the removal of nucleotides, the apo protein (starting from different conformations) does not undergo any major conformational changes in any of the 30 (500 ns long) MD simulations performed in the POPC environment.[35] The very rapid closure of the periplasmic gate (within a few nanoseconds) in the apo Sav1866 simulation reported by Xu et al.[42] is clearly inconsistent with both our results and those reported by Goddeke et al.[35]

Water Density Maps Confirm the Adoption of an IFocc State in the POPE Lipids

A 10 Å cross section of the simulation box centered around the protein transport lumen was used to build the water density map and determine the water accessibility of the lumen in the POPE and POPC lipid environments, results of which are shown in Figure I,J for one of the POPE and POPC trajectories (POPE1 and POPC1). A continuous water flow between the bulk water in the periplasm and the substrate translocation chamber in the protein is observed in the first 0.5 μs of the simulations in both POPE and POPC environments. The periplasmic gate, however, completely closes in the POPE simulations disconnecting the bulk water access to the protein transport lumen at the periplasmic side of the protein. The water density maps based on the last 0.5 μs of POPE1 and POPC1 simulations are shown in Figure I,J, illustrating how the periplasmic gate closes (stays open) in the POPE (POPC) environment. In both environments, however, the cytoplasmic gate stays closed (Figure I,J), which indicates that the protein conformation becomes occluded in the POPE environment and stays in the OF state in the POPC simulations. The water density in the NBD region shows some hydration over the course of the simulation in both POPE and POPC environments; however, the NBD hydration is more pronounced in the POPE environment.

PE Headgroups Mediate Closing of the Sav1866 Periplasmic Gate

The replacement of the primary ammonium (in PE) by quaternary ammonium (in PC) reduces the capability of the lipid headgroups to be involved in hydrogen bonding, which we believe is a major factor contributing to the discrepancy observed in the conformational transitions of Sav1866 in POPE and POPC lipid environments. The PE headgroups are capable of forming hydrogen bonds with more types of amino acid side chains as compared to PC due to the presence of the primary ammonium. Hakizimana et al. show that the functioning of secondary multidrug transporter LmrP is specifically dependent on hydrogen bonding capability of the PE headgroups rather than the bulk properties of the lipids.[22] In addition, the PE headgroups, unlike the PC headgroups, can form interlipid hydrogen bonds among themselves.[64] Here we provide evidence that the lipid–protein and lipid–lipid interaction networks present in the PE lipid environment facilitate the closing of the Sav1866 transporter, while the absence of these interactions along with the sheer presence of the PC lipids in the periplasmic gate of Sav1866 keep the transporter open to the periplasm. By monitoring the outer leaflet lipids around the protein, it is evident that all lipids leave the periplasmic gate in the POPE simulations (Figure A), while many lipids stay in this vicinity in the POPC simulations (Figure B). The presence of the lipids between the two TMD bundles that form the periplasmic gate could potentially prevents the protein from closing by blocking the interactions between the two bundles (Figure B). We hypothesize that the primary ammonium in PE promotes the interactions of the two TMD bundles in a stepwise manner illustrated in Figure C–E and discussed below.

Figure 2

Lipid-mediated periplasmic gate closure. Outer leaflet lipid headgroup occupancy isosurfaces in and around the protein based on POPE1 (A) and POPC1 (B) simulations. (C–E) Stepwise lipid-mediated periplasmic gate closing mechanism of Sav1866 in POPE lipids. TM helices h1 (cyan) and h6 (magenta), which are directly involved in the mechanism, are shown in cartoon representation. Periplasmic residues K38, D42, T276, and T279 are shown in thick licorice representation, while the lipid headgroups and tails are shown in thin ball-and-stick and licorice representations, respectively. Blue dashed lines represent the interbundle, lipid–protein, and lipid–lipid hydrogen bonds. (F) Maximal hydrogen bond occupancy (in percentage) of any single lipid with residues of interest in POPE and POPC simulations. (G) Average number of protein–lipid hydrogen bonds with an occupancy greater than 5% in POPE and POPC simulations. Number of interlipid hydrogen bonds as a function of time formed between the lipids positioned within 8 Å of the periplasmic residues K38, D42, T276, and T279 in POPE1 (H) and any of the POPC (I) simulations. Number of lipid–protein hydrogen bonds as a function of time formed between the lipids and the periplasmic residues K38, D42, T276, and T279 in POPE1 (J) and POPC1 (K) simulations. K38-T279 and D42-T276 hydrogen bond donor–acceptor distance time series for monomer-A in POPE1 (L) and POPC1 (M) simulations. (N, O) Same as parts L and M for monomer-B.

As shown in Figure C, the POPE lipids can form hydrogen bonds with residues K38, D42, T276, and T279 located at the periplasmic side of the transporter. Figure F,G compares the occupancy of these lipid–protein hydrogen bonds for the POPE and POPC lipids. On each side of the periplasmic opening, one can identify K38 and D42 from one TMD bundle (h1, specifically) and T276 and T279 from the opposite TMD bundle (h6, specifically) as residues with side chains interacting with the PE lipid headgroups. D42, T276, and T279 barely form any hydrogen bonds with the POPC lipids. K38 has a comparable interaction with both PE and PC headgroups (due to the fact that K38 interacts with the phosphate group, which is common in both PE and PC headgroups).

When the protein is completely open to the periplasmic side (β > 20°), one may identify two or more lipids bridging the two bundles through a combination of lipid–protein and lipid–lipid interactions (Figure C). The movement of these lipids can potentially pull the h1 and h6 helices (and their corresponding bundles) closer. Once the two bundles are closer (15°< β < 20°), only one lipid molecule is enough to bridge the two bundles (Figure D). Figure H–K shows the number of lipid–lipid and lipid–protein hydrogen bonds between the POPE/POPC headgroups and the 4 periplasmic residues mentioned above. The PC headgroups cannot form hydrogen bonds among themselves (Figure I) as mentioned above while the PE headgroups can form such hydrogen bonds even within the vicinity of the protein (Figure H). Note that the hydrogen bonds quantified here belong to only those lipids that are within 8 Å of the 4 residues mentioned above. In addition to lipid–lipid hydrogen bonding, lipid–protein hydrogen bonding also shows a differential behavior in the POPE and POPC environments (Figure J,K). Note that these protein–lipid and lipid–lipid hydrogen bonds are transient in nature, facilitating the closure of the two TMD bundles without being necessary to keep them together. Negatively charged D42 and positively charged K38 side chains attract the PE lipids by interacting with the positively charged primary ammonium and negatively charged phosphate groups of the PE headgroups. Once the lipids are close enough to T276 and T279, they also interact with the PE headgroups. The presence of these hydrogen bonds is probably necessary for bringing the two TMD bundles close to each other such that they can form protein–protein hydrogen bonds and stabilize the occluded conformation on the periplasmic side. It is also likely that the interaction of PE headgroups with these particular residues on the surface of the protein prevents the lipids from entering the vicinity of the protein periplasmic opening (i.e., what occurs in the POPC lipids preventing the closure of the periplasmic gate).

Once the periplasmic gate closes, two interbundle hydrogen bonds (K38-T279 and D42-T276) form on each side of the periplasmic gate (Figure E). Figure E,L,N illustrates the formation of these interbundle hydrogen bonds in one of the POPE simulations (POPE1), where the closure of TMDs take place in both monomers within the first 0.8 μs of simulations, resulting in the formation of hydrogen bonds between K38 and T279 as well as between D42 and T276. Such hydrogen bonds are not formed in the POPC lipid environment (Figure M,O).

Lipid Tails Play a Role in the Closure of the Periplasmic Gate

Our lipid–protein interaction analysis identifies a hydrophobic lipid binding pocket (in each monomer) between the helices h3, h4, and h6, where the lipid tails from the inner leaflet interact transiently but frequently with various amino acid side chains (Figure A–C). In particular, several h6 residues, such as L297, F301, L304, F208, and M311, interact with the lipid tails as shown and quantified in Figure C,D. We hypothesize that, upon the removal of nucleotides (or the ATP hydrolysis) and subsequent changes in the protein intramolecular interactions, the lipid–protein intermolecular forces become unbalanced. The hydrophobic lipid tails, for instance, interact with the TM helices, particularly h6, slowly pushing the two TMD bundles toward each other. This process is accompanied by the protein–lipid and lipid–lipid interactions involving the outer leaflet PE headgroups discussed above and continues until the system reaches a new equilibrium. Since the outer leaflet PC lipids are preventing the periplasmic gate from closing (rather than facilitating it), the transporter may not be able to close on the periplasmic side despite the presence of inner leaflet lipid–protein interactions. Note that since POPE and POPC lipids have the same lipid tail type, both lipids can interact with the identified hydrophobic pocket; however, the POPE lipids are interacting with a greater probability (Figure D). Of relevance to this observation is the presence of a basic amino acid side chain (R186) in close proximity to the hydrophobic pocket, which can attract the lipid headgroups (PE and PC alike) through hydrogen bonding with their negatively charged phosphate groups (Figure C,E). A greater average hydrogen bond occupancy is observed with the POPC lipids (at least in one of the monomers), as shown in Figure E. The differences observed between POPE and POPC lipids in terms of lipid tail protein interactions and lipid headgroup protein interactions (with the hydrophobic pocket and R186 residue, respectively) could be simply due to the difference in the conformation of Sav1866 in POPE and POPC environments.

Figure 3

Interactions between protein and lipid tails. Hydrophobic lipid binding pocket in monomer-A (A) and monomer-B (B). (C) TM helices and important residues within the lipid binding pocket of monomer-A. h3, h4, and h6 are colored blue, gray, and magenta, respectively. Key residues of the h6 (L297, F301, L304, and M311) interacting with the lipid tails are represented by orange spheres, and lipids binding in the pocket are shown in the ball-and-stick/licorice representation. Residue R186, which interacts with lipid headgroups, is represented by red spheres. (D) Average contact frequency of POPE and POPC lipids with h6 residues based on all three simulation sets. The contact was defined based on a cutoff distance of 3.0 Å. (E) Average hydrogen bond occupancy (%) of R186 from with the POPE and POPC lipids based on all three simulation sets.

TMD Hinge Formation and Interbundle Salt Bridge Formation Promote the Periplasmic Closure

In addition to the protein–protein interactions discussed above, we have also identified a particular interbundle salt bridge (i.e., R81-D145) between TM helices h2 and h3 (of the same monomer), near a hinge region in h3, in POPE simulations. Neither the hinge region, which is formed through disturbing the secondary structure of residues 138–141, nor the R81-D145 salt bridge exists in the Sav1866 crystal structure (Figure A,D). The formation of R81-D145 salt bridge in POPE simulations (Figure B,E,G) is accompanied by the weakening of the intrabundle R296-D145 salt bridge (between helices h6 and h3 of the same monomer) that exists in the crystal structure (Figure C,F,H) and the formation of the hinge region in h3 (as the result of the loss of helical content), as shown in Figure B,E,I. The R81-D145 salt bridge formation, however, does not completely abolish the R296-D145 salt bridge (Figure H).

Figure 4

Formation of h3 hinge and R81-D145 salt bridge in POPE lipid membrane. (A) Crystal structure of Sav1866 with the region of interest highlighted by a red circle. (B, C) Two perpendicular views of the highlighted region in part A. Salt bridge forming residues R81, D145, and R296 and their corresponding helices h2, h3, and h6 are highlighted. (D–F) Same as parts A–C using a representative structure of the transporter resulting from POPE simulations. R81-D145 (G) and R296-D145 (H) salt bridge distance time series in POPE simulations. (I) Helical content of the residues 138–141 in monomer-A in POPE simulations. (J–L) Same as parts G–I in POPC simulations.

None of the events mentioned above (i.e., the h3 hinge formation and the R296/R81-D145 salt bridge swapping) is observed in any of the POPC simulations (the data are only shown for monomer-A in Figure J–L). Interestingly, the significant loss of helical content in h3 and formation of the hinge region as well as the interbundle R81-D145 salt bridge formation, which are observed in all three POPE simulations, only occur in one of the monomers (monomer-A) in the three POPE simulations performed. One may conclude that the h3 hinge formation and R81-D145 salt bridge formation (1) are correlated and (2) are not required for periplasmic closure. However, both events are likely to promote the periplasmic closure even if they happen only in one monomer. Note that the R81-D145 interbundle salt bridge is located close to the middle of the membrane right above the formed hinge region, coupling the two wings of each monomer that form the periplasmic gate (Figure A,D). The R81-D145 salt bridge brings the two bundles closer to each other while the hinge region allows for the easier bending of h3. This eases the movement of the TMD bundle involving h3 of monomer-A (which also includes h6 in monomer-A) toward the opposite bundle involving h2 of monomer-A (which also includes h1 in monomer-A). As a result, the periplasmic sides of h1 and h6 in monomer-A will be positioned closer to each other, and K38-T279/D42-T276 interbundle hydrogen bonds can form, resulting in the closure of the periplasmic gate.

Conformational Changes of Sav1866 in Non-PO Lipids

While the POPE and POPC simulations have a distinct behavior in a very consistent manner, the distinction between the PE and PC lipids with non-PO tails seems to be less clear. Figure A summarizes the β measurements based on the last 1.2 μs of all PO and non-PO simulations, excluding the negative control (also see Figure S4). One may conclude that none of the 9 PC simulations (regardless of the tail type) transitions to the occluded state. However, among the non-PO PE (DOPE and DPPE) simulations, we observe the closure of the periplasmic gate only in one out of three repeats of each type, including DOPE2 (Figure S4A,B) and DPPE3 (Figure S4E,F) simulations. Also, see Figure S5A,C and Figure S5B,D, to examine the formation of interbundle hydrogen bonds among periplasmic residues in the DOPE2 and DPPE3 simulations, respectively. We do not observe the formation of h3 hinge as well as the R81-D145 salt bridge in any of the non-PO lipid simulations, even in DOPE2 and DPPE3, where the periplasmic gate closes. This observation further confirms that these two events are not necessary for the closure of periplasmic gate. However, the lack of the h3 hinge and the R81-D145 salt bridge could potentially explain why the gate closure was only observed in one out of the three repeats of DOPE and DPPE simulations. Therefore, we hypothesize that the PO lipid tails promote the closure of the periplasmic gate in the PE lipids by interacting with the TMDs in the middle of the membrane and promoting the formation of the h3 hinge and the R81-D145 salt bridge.

Figure 5

Sav1866 in non-PO environments. (A) Average β in all PO and non-PO lipid environments based on the last 0.5 μs of the simulations. Projection of DPPE/DPPC (B) as well as DOPE/DOPC (C) trajectories onto their corresponding (PC1,PC2) space. (D) Projection of all PE and PC simulation trajectories (excluding the negative control) onto the (PC1,PC2) space. Representative snapshots of side views of the DOPE (third set, DOPE3; E) and DPPE (first set, DPPE1; F) lipids blocking the protein from closing through interacting with charged residues of the substrate translocation chamber. Positively charged residues (R81, R295, and R296) of the translocation chamber are colored blue, and negatively charged ones (D145 and E288) are colored red. All non-carbon and non-hydrogen atoms of lipids are represented as spheres and carbons as ball and sticks. Oxygen, nitrogen, and phosphorus atoms of lipids are colored red, blue, and gold, respectively. Hydrogen bond interaction occupancies (in percentage) of respective DOPE (G) and DPPE (H) lipids (shown in panels E and F) with the charged residues of the translocation chamber.

Although the periplasmic closure occurs only in 2 out of 12 non-PO simulations, the NBD-associated conformational changes that are the direct consequence of the removal of the nucleotides occur in all 12 simulations. Figures S1 and S3 show the time series of RMSD and NBD distances for all PO and non-PO simulations performed. Comparing to the negative control, all apo simulations show a meaningful increase in terms of both the RMSD from the crystal structure and the NBD distance. PCA analysis was also used to characterize the conformational variation of Sav1866 in PO and non-PO lipid environments. When all six DP trajectories were used in the PCA analysis, the PE and PC trajectories were clustered differently in the (PC1,PC2) space (Figure B). A similar behavior was observed for the DO simulations (Figure C). The behavior of the transporter in PO (Figure F) and DO bilayers is somewhat similar, in that all the triplicate simulations of POPC and DOPC bilayers are tightly clustered compared to those of POPE and DOPE. Interestingly, this is the opposite of what is observed for DPPE and DPPC simulations (Figure B). Generally, by comparing Figure B,C to Figure F, one can observe that the DPPC and DOPE simulations are associated with larger variances as compared to other simulations systems. When all the trajectories (excluding the control simulation) were used in a single PCA analysis and projected onto the (PC1,PC2) space, we could also identify that the DPPC and DOPE simulations are the ones with a much more varying behavior (Figure D). Note that the accumulated contribution of these two PCs to the total variance is ∼35% (Figure S2B), and thus they do not necessarily represent all conformational changes. However, larger variation of transporter in DPPC and DOPE simulations is evident from the RMSD calculations as well (Figure S1). Overall, PCA analysis demonstrates that not only headgroups but also lipid tails are playing a role in structural dynamics of Sav1866. Figure S6 shows more analyses with respect to membrane properties of all six types of lipids used in this study.

We would like to note that while the distinction between the behavior of Sav1866 in POPE and POPC lipids is quite clear, and there is virtually no overlap between the β distributions in the two membranes (see Figure S7A), the distinction between the non-PO PE and PC lipids is less clear (see Figure S7B). We argue that there is a statistically meaningful distinction between the two as discussed in Figure S7; however, it is desired to explore the distinction more extensively either using more/longer MD simulations (which is beyond the current computational capabilities) or employing enhanced sampling techniques. We have recently illustrated the difficulties of employing unbiased MD simulations for the study of large-scale conformational changes in membrane transporters in a statistically sound way.[37] We have also shown how enhanced sampling techniques can be employed to provide more reliable information on large-scale conformational changes of membrane transporters.[53−55]

Finally, we note that, in addition to lipid composition, the lipid phase could play a role in the distinctive behavior of PE vs PC lipids. In this study, DPPE and DPPC lipids were simulated in the gel phase as determined by the order parameter calculations (Figure S8). Our simulations suggest that the differential behavior of Sav1866 in PE and PC lipids is retained even in the gel phase. Although the periplasmic gate closure occurs quite early (t ≈ 100 ns) in the DPPE3 simulation (see Figure S4E), the lipids transition into the gel phase much more quickly as shown in Figure S9 (within a few nanoseconds), and the periplasmic gate closure occurs while the lipids are in the gel phase. It is interesting to note that not only is the Sav1866 homologue P-gp known to function in the gel phase, but also it is reported to transport colchicine with a higher rate in the gel phase than in the fluid phase.[65] It is thus hypothesized that rigidity of the gel phase may not be a barrier for conformational changes of P-gp.[65] We note that, given the limited number of simulations in this study, further investigation is necessary to generalize these findings as discussed above.

Sav1866 is a Flippase

Several charged side chains are located in or close to the translocation chamber of Sav1866, some of which are shown in Figure E,F. Interactions between these charged residues and the lipid headgroups could result in binding of the lipids in a way that blocks the closure of the periplasmic gate as is the case in the examples shown in Figure E,F and quantified in Figure G,H for two non-PO PE lipid environments. These interactions, however, do not necessarily result in blocking of the periplasmic gate closure. If the lipids flip and align themselves along the membrane normal, the periplasmic gate can still close. This is exactly what we have observed in some of the PE simulations.

Flippases aid the movement of lipids between the two leaflets of a lipid bilayer. Several ABC transporters[28] such as TmrAB[25] and P-gp[66] have been reported to have flippase activity. We have observed initial stages of the flippase activity of Sav1866 in the PE lipid bilayers, which occur on relatively short time scales. In particular, we have observed a flippase-like activity in one of the POPE (POPE2) and two of the DOPE (DOPE2 and DOPE3) simulations (Figure A–C). We did not see this behavior in any of the PC or DPPE simulations.

Figure 6

Flippase activity of Sav1866. Side views of POPE (A) and DOPE (B) lipids binding parallel to the membrane normal in the substrate translocation chamber. In these two cases the periplasmic gate is closed. (C) Representative snapshot from the third set of DOPE simulations, where the lipid is binding in the translocation chamber. In this case, the periplasmic gate stays open. (D) Tilt angle of the lipids that are shown in parts A–C as a function of time. (E) Distance between the periplasmic residues K38 and T279 as a function of simulation time. (F–H) Hydrogen bond occupancies (in percentage) of the same lipids shown in parts A–C with the charged residues (E288, R295, and R296, respectively) lining the translocation chamber.

The translocation chamber of the transporter contains charged residues, some of which are close to the center of the lipid bilayer. These residues include E288 (slightly above the center, toward the periplasm) and R295/R296 (slightly below the center, toward the cytoplasm) (Figure A). PE lipid headgroups in the outer leaflet are attracted to these charged residues; they first rotate ∼90° to enter the substrate translocation chamber from the center of the membrane, and then flip toward the cytoplasm (rotate another 90°) such that the lipid headgroup faces the cytoplasmic side, and the lipid tail faces the periplasmic side. We propose that E288 acts as the main driving force to attract the PE lipids to the translocation chamber of the protein, although lipids interact with R295/R296 as well. Since E288 is slightly above the center, the lipid headgroups of the outer leaflet have a greater chance of interacting with E288 as compared to R295/R296. Unlike PC, PE headgroups are capable of hydrogen bonding with E288 due to the presence of the small and positively charged primary ammonium in their headgroups. The PE headgroups are thus more likely to enter the chamber through interactions with E288. They will then move further down by approaching the positively charged R295/R296 located below the membrane center (which interacts with the negatively charged phosphate group of the lipids). This completes the ∼180° rotation of the lipids that interact with these residues (Figure D).

The protein may or may not close on the periplasmic side after the lipids enter the chamber (compare Figure A,B to Figure C). We have seen the periplasmic closure in POPE2 and DOPE2 simulations, which occur after the lipid flipping, made evident by comparing the timelines of lipid flipping and periplasmic closure, quantified using the lipid tilt angle and the K38-T279 interbundle hydrogen bonding, respectively (Figure D,E). In these two cases, interactions between the lipid headgroups and E288 are strong (Figure F), but there are also interactions between the lipid headgroups and R295/R296 (Figure G,H). Even though the lipid molecules diffuse into the protein translocation chamber, the transporter is still able to close on the periplasmic side after the lipid is flipped and oriented along the membrane normal toward the cytoplasm. In the DOPE3 simulation, although a lipid is flipped and oriented along the membrane normal (Figure C,D), protein does not close within the simulation time (Figure C,E). This is due to the presence of other lipids located between the two wings of the transporter and interacting with the charged residues of the translocation chamber (Figure E,G).

The flippase activity and the periplasmic closure are two interrelated processes; they both occur (for Sav1866) only in PE lipids (at least within the microsecond time scale). The PC lipids, due to lack of primary ammonium in their headgroup, cannot interact with E288; therefore, they do not enter the chamber from the center and do not flip. We note that the ABC exporters are known to transport the lipids from the cytoplasmic side to the periplasmic side. However, we have observed the flippase activity in the opposite direction. This is perhaps evidence that the flippase activity is a bidirectional process even in ABC exporters, and it is the coupling with the ATP binding/hydrolysis that gives directional preference, thermodynamically.

Conclusions

A detailed study of the structural dynamics of bacterial ABC exporter Sav1866 using microsecond-level all-atom equilibrium MD simulations gives insight into lipid-dependent conformational dynamics of ABC transporters and provides evidence that the function and mechanism of ABC transporters are lipid-specific. We observed that none of the nine simulations of Sav1866 in the PC lipid environments could induce a state transition, at least within microseconds. However, in all POPE and one-third of other PE lipid environments, the protein underwent large-scale conformational changes resulting in a transition from the OF to the IFocc state.

A hydrophobic lipid binding pocket was identified in each monomer of Sav1866 on the surface of the protein surrounded by TM helices h3, h4, and h6. We propose that the constant bombardment of these binding pockets and the TM helices involved (particularly h6, which participates in the periplasmic gate closure) by lipid tails help to close the protein on the periplasmic side by pushing the two TMD bundles involved in the periplasmic gate closer to each other. The lipid types, particularly the lipid headgroups, however, play a crucial role in either facilitating or preventing the closure of the periplasmic gate. The PE lipids help to bridge the two TMD bundles by promoting the periplasmic interbundle hydrogen bonding (K38-T279 and D42-T276) through protein–lipid and lipid–lipid transient hydrogen bonding, which is absent in the PC lipids due to the lack of a small and positively charged primary ammonium that is present in the PE lipids. On the other hand, the PC lipids prevent the closure by blocking the interbundle protein–protein interactions on the periplasmic side.

We also observed that the POPE lipids influence the protein dynamics in a different manner as compared to the other types of PE lipids. Only the POPE lipids induce the formation of a hinge region in h3 (near the residues 138–141) and an interbundle salt bridge (R81-D145) right above the hinge (close to the membrane center), both of which promote the periplasmic gate closure.

Finally, we observed the partial flippase activity of the transporter in some of the PE but none of the PC simulations, further suggesting a lipid-specific functioning for ABC transporters. Lipids from the outer leaflet of the membrane enter the substrate translocation chamber attracted by the charged residues lining the chamber (particularly E288). The PE lipid headgroups enter the chamber from a site slightly above the membrane center and flip toward the cytoplasm to interact with other charged residues that are slightly below the membrane center. The flippase activity is coupled to the closure of the protein on the periplasmic side.

In conclusion, Sav1866 seems to be tailor-made to function in a POPE dominated membrane. Our study provides evidence and explanation for the alternating access mechanism; however, it calls into question the assumption that the membrane transporters always follow this mechanism even outside their native environment.

Methods

Characterizing the energetically downhill OF → IF conformational transition of an ABC transporter is computationally challenging as it demands MD simulations at the μs or longer time scales. In the current study, the 128-node Anton 2 machine at the Pittsburgh Supercomputing Center (PSC) was employed to achieve μs time scales. The OF state of the Sav1866 transporter (PDB, 2HYD)[33] was modeled in six different PE and PC lipid environments, namely, POPC (1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine), DOPC (1,2-dioleoyl-sn-glycero-3-phosphocholine), DPPC (1,2-dipalmitoyl-sn-glycero-3-phosphocholine), POPE (1-palmitoyl-2-oleoyl-sn-glycero-3-phosphoethanolamine), DOPE (1,2-dioleoyl-sn-glycero-3-phosphoethanolamine), and DPPE (1,2-dipalmitoyl-sn-glycero-3-phosphoethanolamine). The two nucleotide molecules present in the Sav1866 crystal structure (PDB, 2HYD(33)) were removed to obtain a nucleotide-free (apo) state of the transporter, mimicking the state of the system after the ATP hydrolysis, and before the start of the energetically downhill conformational transition. In doing so, we ignore (1) the energy released from the ATP hydrolysis and (2) the role of the ADP.

CHARMM-GUI[67,68] simulation input generator was used for building each lipid–protein system. Each simulation system consists of one protein, ≈360 lipids, 0.15 M NaCl, and ≈45 000 TIP3P[69] water molecules. Overall, the total size of the system was ≈134 × 142 × 174 Å3. We also generated a “negative control” with 2 Mg-ATPs bound to the NBDs of Sav1866 embedded in a POPE lipid bilayer, similar to that made for the apo protein.

CHARMM36 all-atom additive force field parameters[70,71] were used to model the systems. Preliminary MD simulations were performed using NAMD 2.10[72] prior to Anton 2 production runs. Each system was simulated with periodic boundary conditions at 310 K temperature using a Langevin integrator with a time step of 2 fs and collision frequency of 0.5/ps. The DPPE and DPPC lipids are in the gel phase at this temperature, while the remaining four lipid bilayers (i.e., the POPE, POPC, DOPE, and DOPC) are in the fluid phase. A 1 atm pressure was maintained using the Nosé–Hoover Langevin piston method.[73,74] The cutoff for the nonbonded interactions was 12 Å (with a smoothing function, starting at 10 Å), and the particle mesh Ewald (PME) method[75] was used to compute the long-range electrostatic interactions. Prior to equilibration, each system was energy minimized for 10 000 steps using conjugate gradient algorithm[76] and further relaxed using a multistep restraining procedure explained elsewhere[67] for ∼1 ns. The initial relaxation was done in NVT, followed by a 5 ns equilibrium simulation in the NPT ensemble.

Three different structures from each 5 ns preliminary simulation were obtained for the 2.4 μs long production runs on Anton 2 (i.e., three sets of simulations for each lipid environment). For the control simulation, only one production run was conducted. Each production run was carried out using a 2.5 fs time step. The pressure was maintained at 1 atm semi-isotropically, using the MTK barostat, while the temperature was maintained at 310 K, using the Nosé–Hoover thermostat.[73,74] The long-range electrostatic interactions were computed using the fast Fourier transform (FFT)[77] method implemented on Anton 2. Overall, we performed 19 simulations on Anton 2, each for 2.4 μs, with an aggregate simulation time of ∼46 μs. The configurations were collected every 2.4 ns, generating 1000 configurations for each 2.4 μs trajectory. Visual molecular dynamics (VMD)[78] was used for visualization, making molecular snapshots, and some of the analyses.

In each monomer, TM helices h1, h2, h3, h4, h5, and h6 consist of residues 1–44, 52–106, 117–159, 161–216, 218–272, and 277–319, respectively. Each NBD consists of residues 321–578 of each monomer. Bundle-1 is composed of residues 1–107 of monomer-A (TM helices h1 and h2) and residues 116–320 of monomer-B (TM helices h3–h6), whereas bundle-2 is composed of residues 1–107 of monomer-B and 116–320 of monomer-A. Interbundle angle β was calculated as the angle between the third principal axes of the two TMD bundles. For calculating hydrogen bonds, salt bridges, secondary structure, membrane properties (such as thickness, interdigitation, deuterium order parameter (SCD), and tilt angle), and water and lipid densities we used VMD and its various plug-ins.[78−80] The density isosurface calculations were used to estimate the probability of finding water and lipid molecules around and within the protein. Cross sections of such density isosurfaces of water were used to determine the accessibility of water within/around the protein in different lipid environments. The cross sections were taken along the yz plane, and a slice of about 10 Å was considered along the x axis. For generating the lipid density profiles, only lipids in the outer leaflet within 5 Å of the protein were considered. In the hydrogen bond analysis, the cutoff distance and angles considered were 3.5 Å and 30°, respectively. Salt bridge distances were calculated as the minimum distance between donor and acceptor atoms. Principal component analysis (PCA) was carried out using PRODY software.[81] Only Cα atoms were considered for the PCA calculations. Distances between the two NBDs were calculated as the distance between the centers of mass of their Cα atoms.

Area per lipid (APL) was calculated using GridMAT-MD software,[82] which follows the polygon-based tessellation approach, and it uses the phosphorus atom of the headgroup of each lipid as a reference point. Membrane thickness was measured as the distance between the average positions of the phosphorus atoms of the two leaflets. Width of the overlap region of the lipid acyl chains belonging to the two leaflets was measured using all heavy atoms. The cutoff distance for the contacts was set to 4 Å. The lipid tilt angle was measured as the change of vector connecting the phosphorus atom in the headgroup to the terminal carbon of the acyl chains. The SCD for each methyl group of the lipid tails was calculated by averaging 3 cos2θ – 1 over all lipid molecules of single/multiple snapshots of the simulation system for instantaneous/time-averaged SCD values. Here θ is the angle between each C → H vector and the membrane normal.