Multiscale simulation of surfactant–aquaporin complex formation and water permeability

Abstract

Multiscale simulation has been conducted for the formation of a surfactant–protein complex that uses sodium dodecyl sulfate (SDS), a negatively charged surfactant, and aquaporin Z (AqpZ), a membrane protein that facilitates water transport across lipid membranes. A detailed analysis of the molecular driving forces of the self-assembly at different pH values reveals distinctive contributions of electrostatic and hydrophobic interactions to the complex structure and formation kinetics. The electrostatic interactions become more significant at low pH and are responsible for the formation of larger complexes. A comparison of the protein conformation in the SDS complex with that in the lipid bilayer of palmitoyloleoyl phosphatidyl ethanolamine (POPE) shows that the SDS molecules have only marginal effects on AqpZ conformation including the water channel structure. The simulation indicates that AqpZ preserves its secondary structures after being bound with SDS molecules, while the arrangement of the helical structures leads to a coiled-coil to single helix transition as suggested by experiments. AqpZ may lose water permeability either due to the blockage of the water channel by individual SDS molecules or due to the attachment of micelle-like structures at the hydrophilic ends of the water pore. Reconstitution of the AqpZ complex into a POPE bilayer shows that the membrane protein regains its activity after the complete removal of SDS molecules from the protein pores. The molecular insights gained from multiscale simulation will be helpful for future development of AqpZ-embedded membranes.

---

Introduction

Aquaporins are homotetrameric membrane proteins¹ prevalent in a wide range of organisms including bacteria, archaea, plants, insects and mammals.² These proteins contain water channels that are essential for the cellular regulation of osmotic pressure in response to the fluctuations of water content in the extracellular environment. Aquaporins facilitate water permeation across biological membranes at an impressive rate of 3 × 10⁹ molecules per second per channel.^3,4 The high permeability and near perfect selectivity make aquaporin an ideal candidate for water purification and biomedical applications. Over the past decades both theoretical and experimental approaches have been actively pursued for extraction of aquaporins from bacteria cells and their incorporation into artificial membranes⁵ such as supported lipid bilayers^6,7 and block copolymer lamellae.^8,9 From a practical perspective, one major challenge for the fabrication of the biomimetic membranes is to reconstruct aquaporin proteins into an artificial membrane that preserves water permeability and selectivity.⁹ For that purpose, the membrane proteins must be extracted from the bacteria cells with suitable surfactants.¹⁰ To ensure the protein integrity and functionality, practical applications require a good understanding of the surfactant–protein complex formation and the influences of surfactant molecules on the protein structure and water permittivity.

A number of experimental procedures have been applied to studying the interactions of aquaporin proteins with surfactant molecules. Because membrane proteins and surfactant molecules have distinctly different electron densities, their complex formation makes the analysis of SANS or SAXS data extremely difficult.¹¹ Recently, Berthaud et al. reported that, at adequate conditions, n-dodecyl-β-D-maltopyranoside molecules and aquaporins form an elliptical, toroid-shaped corona in an aqueous solution. The shape of the protein–surfactant complex was deduced from complementary results from size exclusion chromatography, SAXS, refractometry experiments and coarse-grained modelling.¹² Because this work accounts for only surfactants adsorbed at the hydrophobic trans-membrane regions of aquaporin proteins, it provides no information on the microscopic structures of the entire aquaporin–surfactant complexes, in particular those pertinent to the hydrophilic domains that directly affect the functionality of the water channels. Little is known on how the surfactant affects molecular level events underpinning the kinetics of water transport through aquaporin channels. Direct inspections of the structure and functionality of membrane proteins surrounded by surfactant molecules are experimentally challenging. For example, the strong binding of sodium dodecyl sulfate (SDS), a negatively charged surfactant, with membrane proteins changes the extinction coefficient of the peptide carbonyl group, causing difficulties for the interpretation of the secondary structure from circular dichroism (CD) data.¹³

In the present work, we use all-atom and coarse-grained molecular dynamic simulations to investigate the interactions of SDS surfactant molecules with aquaporin Z (AqpZ) and the effects of the surfactant on the protein stability and water permeability at two different pH values.^14–16 The objective is to identify the molecular driving forces of the surfactant–protein complex formation and examine the effects of SDS on the aquaporin structure and functionality. SDS has been used for the extraction of AqpZ from recombinant E. coli and in subsequent reconstitution of the surfactant–protein complex into a lipid bilayer. It is known that the surfactant may change the secondary or tertiary structures of certain membrane proteins such as diacylglycerol kinase,¹⁷ bacteriorhodopsin^18,19 and rhodopsin.²⁰ However, the effect of SDS on the AqpZ structure and water permeability remains inconclusive. Previous research indicates that AqpZ can form stable tetramers in SDS solutions. While CD spectroscopy confirms the stability of the helical structure at high SDS concentrations, it is not clear whether the tertiary structure was also stable under similar conditions.²¹ By examining SDS–AqpZ interaction in a random mixture of the protein and the surfactant at pH 4 and 7, we will attain a better understanding on the effects of SDS on both singular helix and the tertiary conformation of the protein. The simulation will also allow us to investigate the effect of SDS binding on water permeability.

Materials and methods

Molecular models

Protein model. All simulations conducted in this work are concerned with aquaporin Z (AqpZ), a common membrane protein from Escherichia coli. The all-atom structure of this protein is available from the RCSB Protein Data Bank (PDB code: 1RC2). We use the united atom force field Gromacs ffG43a1 (ref. 22) to describe the protein interacting with water and surfactant molecules. To determine the charge distribution of AqpZ at different pH conditions (pH = 4 and 7), we obtained the internal pK_a values of the titratable amino acids (i.e., lysine, arginine, histidine, glutamic acid, and aspartic acid) from the H++ web tool (http://biophysics.cs.vt.edu). A brief summary of the protein net charges and the total numbers of atoms at different pH values are given in ESI (Table S1†). For efficient sampling, we have also developed a coarse-grained model for AqpZ on the basis of the MARTINI force field.²³ Briefly, each amino acid is represented by one backbone site and up to 4 interaction sites are used for each side chain. The coarse-graining procedure is consistent with the standard rules for the MARTINI force field.

SDS, POPE and water models. The atomic parameters for SDS were determined from existing atom types and force constants in the Gromacs ffG43a1 force field. The isoelectric point for SDS is 1.0, indicating that the surfactant molecules are in the deprotonated state at both pH = 4 and 7. The atomic charges for the SDS molecules were adopted from the literature.²⁴ It has been shown by Sammalkorpi et al.¹⁶ that the SDS model is able to reproduce the experimental results for the SDS micelle structure, the degrees of ionization and hydration. The parameters for palmitoyloleoylphosphatidylethanolamine (POPE) were adopted from Berger and co-workers.²⁵ The simple point charge (SPC) model was used for water molecules.

The MARTINI model for SDS includes one Qa site for the hydrophilic head group and three C1 sites, each representing four CH₂/CH₃ groups, for the hydrophobic tail.²³ The water cluster, Na⁺ and Cl⁻ ions are modelled as single polar sites, P4, Qd and Qa, respectively. We assume that all polar sites have the mass of 72 a.m.u.

Simulation methods. We used the GROMACS 4.5.3 package for all molecular dynamics (MD) simulations described in this work.²⁶ In all-atom simulations, a typical setup includes an AqpZ tetramer placed at the centre of a 10 nm cubic box. Fig. 1(a) shows a representative snapshot of the starting configuration. In the initial configuration, 300 SDS molecules are placed randomly inside the simulation cell. The molar ratio of SDS to AqpZ conforms to the typical range of SDS concentration used in experiments.²⁷ At pH = 7, the AqpZ tetramer and SDS molecules were solvated with 23 403 SPC water molecules, 390 Na⁺ ions, and 90 Cl⁻ ions. At pH = 4, the number of water molecules, sodium ions, and chloride ions were 23 441, 350, and 90, respectively. The slight difference in the particle numbers reflects the change in protonation of amino acid residues at different pH values. In both cases, the salt concentration amounted to 0.15 M in the bulk, and there was no net charge for the simulation box. All-atom MD simulations were performed with periodic boundary conditions in the NPT ensemble. The Particle-Mesh Ewald method was used to account for the long-range electrostatic interactions.²⁸ The covalent bonds were constrained by LINear Constraint Solver²⁹ (LINCS). In all simulations, the temperature and pressure were fixed at 310 K and 1.0 bar, respectively. The Nose–Hoover thermostat³⁰ was used for temperature coupling (with coupling constant of 0.2 ps), and the Parrinello–Rahman method³¹ was used for isotropic pressure coupling (with coupling constant of 5 ps). In all-atom simulations, the systems were first subjected to 5000 steps of steepest-descent energy minimization. The integration step was 0.002 ps. To investigate the assemble dynamics of SDS molecules, 100 ns long MD simulation were applied to small ions, water and surfactant molecules without changing the AqpZ position and conformation.

Fig. 1 The initial random configurations of AqpZ and SDS molecules for all-atom (a) and coarse-grained simulations (b). For clarity, the water molecules are not shown.

To study the complex structure and the configurations of SDS molecules around the AqpZ protein, we carried out five rounds of simulated annealing with the temperature decreased from 450 K to 250 K within 10 ns. The positions of AqpZ atoms were fixed during the simulated annealing processes. It has been shown³² that simulated annealing can speed up the system to reach equilibrium by 10 times without compromising the accuracy. After simulated annealing, MD simulation were carried out for 40 ns at 310 K so that the configuration of AqpZ can be fully relaxed in the presence of water and SDS molecules.

Fig. 1(b) shows a snapshot of the starting configuration of the protein–surfactant complex in the coarse-grained molecular dynamics (CGMD) simulation. The system includes an AqpZ tetramer, 300 SDS molecules, 60 338 water molecule, 1023 Na⁺ and 723 Cl⁻ ions. NPT ensemble MD simulation was carried out in a cubic simulation box, 20.0 nm for the side length, with the periodic boundary conditions in all directions. As in all-atom MD simulations, the particle-mesh Ewald method was used to deal with the electrostatic interactions.³³ The temperature was again fixed at 310 K, while the pressure was fixed at 1.0 bar, both implemented with the Berendsen method with a coupling coefficient of 2.0 ps. After initial energy minimization, we ran coarse-grained molecular dynamics for 500 ns with a 0.025 ps integration step.

To investigate water permeability through AqpZ, we incorporated the protein tetramer, either with or without SDS molecules, into a preformed POPE lipid bilayer according to the protocol proposed by Wolf et al.³⁴ Subsequently, we conducted energy minimization and 1.0 ns MD dynamics with the protein configuration fixed. After the removal of the constraint for protein configuration, the system runs for 50 ns MD steps to study water permeation through the aquaporin channel.

Analysis method. The simulation results were analysed with the tools provided by the GROMACS (version 4.5.3) software package plus home-written codes. The snapshots were drawn with VMD (version 1.9). The pore diameter profiles of AqpZ tetramer were analysed using the HOLE software.³⁵ These profiles were obtained by averaging over 100 snapshots, taken every 10 ps during the last 1.0 ns of the simulation. To monitor the structural changes of AqpZ tetramer, we followed the RMSD of the AqpZ backbone, namely all atoms in the main chain. The stability of the protein pores was evaluated in terms of all atoms from amino-acid residues lying in the water channel, i.e., Phe43, Gly60, His61, Phe62, Asn63, Pro64, Ala65, His174, Thr183, Ser184, Val185, Asn186, Pro187, Ala188 and Arg189. A cylindrical pore with a radius of 0.8 nm was used to evaluate the statistics of water permeation events during the last 20 ns of MD simulations with AqpZ reconstituted in the POPE bilayer. The upper and lower ends of the cylinder coincide with the positions of ASN63 and ARG189 residues, which are located at the entrance and exit regions of the AqpZ pore, respectively. If a water molecule enters from either end of the pore and leaves from the opposite side, we consider this water molecule having successfully penetrated through the channel and count it as a permeation event.

Results and discussion

The dynamics of SDS assembling around AqpZ

To investigate the effect of pH (or equivalently electrostatic interactions) on the dynamics of the complex formation, we conducted 100 ns all-atom NPT simulations (T = 310 K, P = 1.0 bar) at two different pH values (i.e., pH = 4 and 7). The MD simulation was carried out right after energy minimization by simulated annealing with a fixed AqpZ configuration. Approximately, the dynamics of complex formation between SDS and AqpZ molecules can be measured by monitoring the number of neighbouring surfactant molecules, i.e., those with the centre of mass to the closest atom at the AqpZ surface less than 1.0 nm. The distribution of surfactant molecules around the protein can be directly monitored with MD simulation.

Fig. 2 presents the variation of the number of neighbouring SDS molecules around the AqpZ surface with time within the first 10 ns of the MD simulation. At pH = 7, the number of neighbouring SDS molecules increases from 70 to 100 during the first nanosecond of simulation, and reaches an equilibrium value (∼110) in a few nanoseconds. At pH = 4, the protein has a net charge of +10 e per monomer and attracts more SDS molecules at equilibrium (∼150). The larger protein–SDS complex can be explained by the increased electrostatic attraction between the protein and surfactant molecules. The stronger electrostatic attraction also accelerates the dynamics of self-assembly.

Fig. 2 The number of SDS molecules with the centre of mass less than 1.0 nm from the protein surface in the first 10 ns of all-atom MD simulation for the surfactant–AqpZ complex formation.

To analyse the driving forces of protein–surfactant complex formation, we monitored the motions of the hydrophobic and hydrophilic segments of the surfactants moving toward the protein. Fig. 3 gives the evolution profiles for the positions of SDS heads and tails to the centre of mass of AqpZ tetramer within the first 10 ns of MD simulation. At the early stage of the assembling process, SDS heads are closer to AqpZ than SDS tails on average, indicating that the SDS molecules prefer an orientation with their heads pointing toward AqpZ. When SDS molecules are far from the AqpZ surface, the long-range electrostatic force is stronger than the van der Waals force, leading to the preferential attraction of negatively charged SDS heads. Near equilibrium configurations, however, the SDS molecules are so close to AqpZ that their orientations are dominated by the short-range hydrophobic interactions. While the mean distance of surfactant heads near the surface varies little with time, significant changes can be identified in the SDS tail positions, suggesting that the dynamics of the protein–surfactant complex formation is mainly determined by the hydrophobic effects.

Fig. 3 Evolution of the distance between SDS heads or tails and the AqpZ tetramer at pH = 4 and 7. Here the distance is defined as the centre of mass for all SDS tails/heads to that of the tetramer.

Fig. 4 shows the average electrostatic (EC) and the Lenard-Jones (ELJ) binding potentials between AqpZ tetramer and SDS molecules at different pH values. Here the binding energy is defined as the change in the total non-bonded average energy between AqpZ and all SDS molecules in the system at the beginning and the final 1 ns of MD simulation. Both the electrostatic and LJ energies are enhanced at pH = 4 because, as shown in Fig. 2, the complex contains more SDS molecules when AqpZ carries a net positive charge. Not surprisingly, the binding energy is dominated by LJ interactions regardless of the protein charge, suggesting the importance of hydrophobic interactions in complex formation.

Fig. 4 The binding energy between AqpZ and SDS molecules. Here EC and ELJ stand for, respectively, the electrostatic and LJ contributions to the change in the total non-bonded energy of the system during the 40 ns MD simulation.

In addition to the overall distance between the protein and the surfactant molecules, we have calculated the radial distribution functions (RDF) of SDS heads and tails from the surface of AqpZ tetramer. Here RDF is defined as the reduced SDS head/tail density (i.e., the local density of the centre of the mass/average SDS density) as a function of the distance from the mass centre of the protein tetramer. In Fig. 5(a) and (b), each line corresponds to 1 ns sampling interval. While both SDS tails and heads are moving toward AqpZ during the first 30 ns of MD simulation for the assembling process, the change in the local density of the SDS tails is more dramatic than that for the SDS heads, indicating again that the hydrophobic interaction between SDS and AqpZ dominates the dynamics of the complex formation.

Fig. 5 Evolution of radial distribution function (RDF) of SDS tails (a) and heads (b) in all-atom MD simulation of the assembling process. Each line represents 1 ns interval.

According to Fig. 2 and 5, the dynamics of surfactant–protein complex formation may be approximately divided into two stages. The first stage occurs within about the first 2 ns of MD simulation, in which SDS molecules move fast towards AqpZ, driven primarily by electrostatic interactions. In the second stage, which occurs between 2 ns and 30 ns, the SDS tails continue approaching AqpZ even though there is no appreciable movement of the mass centre of the entire SDS molecules. At this stage, the dynamics is dominated by the rearrangement of SDS molecules around the AqpZ surface. This stage is driven by the hydrophobic interactions between SDS and AqpZ, leading to a more compacted structure of the complex. As shown by the arrows in Fig. 5, the RDF curves do not change much after 30 ns, indicating that the time scale of the self-assembly is around 30 ns. This time scale is consistent with previous studies for the formation of SDS micelles in aqueous solutions.¹⁶ It should be noted that the time scale of protein structure evolution do not necessarily coincide with that pertinent to the self-assembly of SDS molecules, i.e. the protein configuration may change in a time scale much longer than that for the SDS assembling process.

The structure of SDS–AqpZ complex

We have analysed the equilibrium structure of the SDS–AqpZ complex using both all-atom and coarse-gained models. To obtain an equilibrated structure of SDS within the atomistic model, we conducted 10 ns simulated annealing before each production run. To ensure that the equilibrium state is reached in all-atom simulations, we compared the simulation results with those from a 500 ns long coarse-grained MD simulation. Fig. 6(a) and (b) show the RDFs of the SDS tails and heads around the geometry centre of the AqpZ tetramer at equilibrium from all-atom and coarse-grained simulations, respectively. We can see that the profiles from coarse-grained agree well with those from all-atom simulations, confirming the validity of simulated annealing involved in all-atom MD simulations.

Fig. 6 RDFs of SDS tails and heads from the geometry centre of the AqpZ tetramer for all-atom MD (a) and coarse-grained MD (b). The snapshots show three representative positions of the surfactant molecules in all-atom MD (c–e).

The RDF profiles of SDS tails shown in Fig. 6(a) allow us to identify three major clusters of peaks around 0.5 nm, 2.0 nm, and 3.5 nm. Similar distribution patterns are observed for the RDF of SDS heads, but at 1.0 nm, 2.5 nm and 4.0 nm, respectively. These patterns suggest that the distribution of SDS molecules around AqpZ can be classified into three different types. Recall that, approximately, the AqpZ tetramer has a cylindrical shape with the diameter and the height about 6 nm and 3 nm, respectively. Accordingly, the nearest peaks shown in Fig. 6(a) should be affiliated with the SDS molecules migrated inside the AqpZ channels. The snapshot shown in Fig. 6(c) confirms the existence of SDS molecules inside the water pore. Although there is only one SDS molecule inside the channel, its influence on water permeability is significant. As to be discussed later, these SDS molecules may occlude the permeation of water molecules.

The relatively weak peaks around 2–3 nm are affiliated with SDS molecules accumulated at the entrance of the AqpZ pore. This hypothesis is supported by a snapshot at the hydrophilic side of the AqpZ pore as shown in Fig. 6(d). In this case, some SDS molecules do not bind to the AqpZ surface. Instead, they associate with other SDS molecules to form a micelle-like structure, which then attaches to the AqpZ surface by hydrophilic interactions. The weak binding between the micelle and protein is evident in the wider and more fluctuated peaks between 2–3 nm in the RDFs of both SDS tails and heads shown in Fig. 6(a) and (b).

It is interesting to note that the farthest peaks are most distinctive in comparison with other peaks in the RDFs for both SDS heads and tails. These peaks correspond to surfactants adsorbed at the hydrophobic regions of AqpZ, namely the trans-membrane regions, as shown in Fig. 6(e). The large value of the farthest peak in the RDF of SDS tails implies a strong hydrophobic interaction between SDS tails and AqpZ. The assembly of SDS around the trans-membrane region of AqpZ is rather different from micellization, but resembles the random adsorption of the surfactant molecules on the protein surface. Recently, Berthaud et al. demonstrated a similar adsorption behaviour for complex formation between n-dodecyl-β-D-maltopyranoside molecules and aquaporin. The formation of an elliptical, toroid-shaped corona around the aquaporin surface has been validated with complimentary results from size exclusion chromatography, SAXS, and refractometry.¹² For SDS tails, the farthest peak appears around r = 3.5 nm; while for SDS heads, the farthest peak appears at 4 nm. The difference in these peak positions suggests that, on average, the hydrophobic tails are closer to the AqpZ surface, and that the hydrophilic heads stand up normal to the surface. The equilibrium structure is born out of the conformation changes of SDS molecules after being adsorbed at the AqpZ surface (see Fig. S1†). The self-organization of SDS molecules around the hydrophobic region of AqpZ resembles surfactant adsorption on carbon nanotubes (CNTs) as studied before by both experiments³⁶ and simulations.³⁷ Because both AqpZ and CNT are extremely hydrophobic at the side surface, it is not surprising that the SDS adsorption on the trans-membrane region of AqpZ is similar to surfactant adsorption at the CNT surface.

AqpZ structure and water permeability

In addition to surfactant distributions, we have examined the AqpZ structure and functionality. Table 1 gives the root-mean-square deviation (RMSD) of the AqpZ backbone and RMSD for the water channels in the surfactant aggregation for the last 1 ns (10 ns) of all-atom (coarse-grained) MD simulation. Also shown are the corresponding RMSD results for the protein in palmitoyloleoyl phosphatidylethanolamine (POPE) bilayer. Both all-atom and coarse-grained models indicate that there is no significant difference among the RMSD values for the AqpZ tetramer at different environments. In other words, the AqpZ structure in the surfactant solutions is similar to that in the native lipid environment. To further study the protein structure, Table 2 gives the percentage of the helical structure of AqpZ in its crystal state in comparison with those at different solution conditions in all-atom MD simulation. We find that the secondary structure does not show appreciable changes in three different cases, which is consistent with the previous experimental investigation.

Table 1 Structures of AqpZ in surfactant solutions and in POPE bilayer

	RMSD (nm)
	Backbone of AqpZ tetramer	Atoms of residues in water pore
In POPE bilayer	0.193 (0.004)	0.234 (0.010)
In pH 4 at 310 K	0.201 (0.005)	0.236 (0.011)
In pH 7 at 310 K	0.212 (0.007)	0.229 (0.012)
CG, in POPE at 310 K	0.524 (0.010)	0.558 (0.013)
CG, in SDS at 310 K	0.535 (0.012)	0.558 (0.015)

---

Table 2 AqpZ helical structure in crystal form and in different SDS solutions

	Percentage of helical structure	Standard deviation
Crystal	68.9%	—
In pH 4, at 310 K	68.2%	0.776%
In pH 7, at 310 K	68.7%	0.859%

---

Fig. 7 A and B present the evolutions of the angle between helix 1 and helix 2 of AqpZ monomer in all-atom MD and coarse-grained MD simulations, respectively. Both all-atom and coarse-grained simulations show a decrease of the angle in the first 30 ns, which corresponds to the alignment of AqpZ helices in the SDS environment. The SDS-induced transition of the helical arrangement in AqpZ corroborates the circular dichroism results reported by Hansen et. al.²¹ who demonstrated that SDS leads to transition from a coiled-coil helix to a single helix.

Fig. 7 Evolution of the angle between the axis of the first and second helices in AqpZ monomer from all-atom MD (a) and coarse-grained MD (b) simulations.

Fig. 8 gives the radial profile, the radius of water channel along Z direction, for AqpZ in the POPE bilayer and in different SDS solutions. We see minor deviations of the AqpZ structures in both natural lipid environment and SDS solutions from that in the crystal structure. The minor deviations indicate that the intact structures of water channels are maintained in the SDS solutions even though the arrangement of the helixes changes when embedded in the SDS aggregations.

Fig. 8 Radial profile of the water channel for ApqZ in POPE bilayer. Z stands for the position of the foremost surface atom along cylindrical direction of the water channel (see Fig. 1).

The structural integrity of AqpZ is essential to preserve its superior water permeability. Will the small changes in AqpZ structure by SDS affect water permeation? To address this question, we carried out three additional sets of MD simulations. In Case 1, AqpZ with all neighbouring SDS molecules was embedded into a POPE lipid bilayer. In Case 2, AqpZ with SDS molecules attached at the trans-membrane region was embedded into a POPE lipid bilayer. In Case 3, only AqpZ is inserted into the lipid bilayer. Fig. 9 gives a snapshot of the starting conformation for each case.

Fig. 9 Starting conformation for MD simulation of AqpZ embedded in POPE bilayer. The AqpZ in these MD simulations is extracted from AqpZ SDS assemble with different residual SDS. (a) AqpZ with all neighbouring SDS molecules; (b) AqpZ with SDS molecules attached at the trans-membrane region; (c) Only AqpZ in the assembly.

Fig. 10 shows typical configurations of the water chain inside each AqpZ channel. In Case 1, there is no water permeation through the AqpZ channels, indicating that the SDS molecules at the pore entrance or exit, as well as those penetrating into pore, prohibit the transport of water molecules. The snapshot shown in Fig. 10(a) indicates no single-filed water chain formed in the pore region of AqpZ, affirming the breakage of the water hopping process. In Case 2, two of the pores are totally open, as indicated by the persistent water permeation (Table 3). In Case 3, three water pores of the AqpZ tetramer are opened for water permeation, and the fourth one is permanently closed. The open or close status of each channel depends on a randomly determined orientation of amino acid residue Arg189, located at the restriction site of the pore region. The single-filed water chain forms only when the Arg189 side chain is parallel to the channel surface, namely, in its up state.³⁸ In both cases, a single file of water molecules is formed inside the open channel, as shown in Fig. 10(b) and (c). The simulation results suggest that, despite a minor alternation in the conformation of AqpZ caused by SDS binding, the AqpZ functionality is preserved after the protein is reconstituted into the lipid bilayer. While the SDS molecules surround the hydrophobic trans-membrane part of the AqpZ protein do not hinder the water permeation, those at the entrance and exit of the aquaporin pore block the water transport.

Fig. 10 Typical water structure in the pore region of AqpZ corresponding to the three cases shown in Fig. 9. Only one subunit is shown for each case; water, SDS and lipid molecules are not shown for clarity.

Table 3 Water permeation in the last 20 ns of the NPT-MD simulations

Cases	Permeation events
	Monomer A	Monomer B	Monomer C	Monomer D	Total
1	0	0	1	0	1
2	3	6	9	7	25
3	9	3	16	0	28

---

Conclusions

We used all-atom and coarse-grained molecular dynamics (MD) simulation to examine the interaction of a common anionic surfactant (SDS) with aquaporin (AqpZ) in an aqueous solution at neutral and acidic pH values. The simulation results indicate that electrostatic forces drive the formation of SDS and AqpZ complex at the beginning stage while the hydrophobic forces play a major role at the late stage. The analysis of the radial distribution functions (RDF) of the SDS heads and tails suggested that SDS molecules are distributed around AqpZ in three different modes, namely (1) monolayer formation at the hydrophobic surface of the trans-membrane region of the protein by random adsorption; (2) formation of micelle-like structures at the entrance of AqpZ channels, and (3) insertion of SDS molecules inside the water pore of the AqpZ.

The simulation suggested that, although the secondary structures of AqpZ remain intact in SDS solutions, the helix arrangement might change from a coiled-coil to a single helix configuration as suggested by earlier experiments. SDS molecules lead to only a minor change of the radius profile of the monomeric water pores. Reconstitution of SDS assembled AqpZ into POPE lipid bilayer showed that water permeation is unaffected when SDS molecules bound at the entrance or exit of the water pores are removed. In other words, SDS may influence the AqpZ activity only by clogging at the monomeric pores. SDS molecules bound at the hydrophobic region have negligible influence on the protein functionality. The simulation results provide molecular-level insights into the formation of the SDS–AqpZ complex useful for the design, synthesis and practical applications of AqpZ-incorporated membranes.

Acknowledgements

The authors are grateful to the Chinese National Natural Science Foundation (no. 21276138), 973 Project (no. 2013CB733604), Tsinghua University Fundamental Project (no. 20131089301) and the National Science Foundation of the US for the financial support of this project.

Notes and references

1. G. Ren, A. Cheng, V. Reddy, P. Melnyk and A. K. Mitra, J. Mol. Biol., 2000, 301, 369–387.
2. G. M. Preston, T. P. Carroll, W. B. Guggino and P. Agre, Science, 1992, 256, 385–387.
3. M. L. Zeidel, S. V. Ambudkar, B. L. Smith and P. Agre, Biochemistry, 1992, 31, 7436–7440.
4. M. J. Borgnia, D. Kozono, G. Calamita, P. C. Maloney and P. Agre, J. Mol. Biol., 1999, 291, 1169–1179.
5. C. H. Nielsen, Anal. Bioanal. Chem., 2009, 395, 697–718.
6. G. F. Sun, T. S. Chung, K. Jeyaseelan and A. Armugam, Rsc Adv, 2013, 3, 473–481.
7. G. F. Sun, T. S. Chung, K. Jeyaseelan and A. Armugam, Colloids Surf., B, 2013, 102, 466–471.
8. H. Wang, T.-S. Chung, Y. W. Tong, K. Jeyaseelan, A. Armugam, Z. Chen, M. Hong and W. Meier, Small, 2012, 8, 1185–1190.
9. M. Kumar, M. Grzelakowski, J. Zilles, M. Clark and W. Meier, Proc. Natl. Acad. Sci. U. S. A., 2007, 104, 20719–20724.
10. G. F. Sun, H. Zhou, Y. Li, K. Jeyaseelan, A. Armugam and T. S. Chung, Colloids Surf., B, 2012, 89, 283–288.
11. J. G. Grossmann, J. Appl. Crystallogr., 2007, 40, S217–S222.
12. A. Berthaud, J. Manzi, J. Perez and S. Mangenot, J. Am. Chem. Soc., 2012, 134, 10080–10088.
13. R. Renthal, Biochemistry, 2006, 45, 14559–14566.
14. S. Jalili and M. Akhavan, Biophys. Chem., 2011, 153, 179–186.
15. R. Braun, D. M. Engelman and K. Schulten, Biophys. J., 2004, 87, 754–763.
16. M. Sammalkorpi, M. Karttunen and M. Haataja, J. Phys. Chem. B, 2007, 111, 11722–11733.
17. F. W. Lau and J. U. Bowie, Biochemistry, 1997, 36, 5884–5892.
18. P. Curnow and P. J. Booth, Proc. Natl. Acad. Sci. U. S. A., 2007, 104, 18970–18975.
19. V. Krishnamani, B. G. Hegde, R. Langen and J. K. Lanyi, Biochemistry, 2012, 51, 1051–1060.
20. A. Dutta, T. Y. Kim, M. Moeller, J. N. Wu, U. Alexiev and J. Klein-Seetharaman, Biochemistry, 2010, 49, 6329–6340.
21. J. S. Hansen, A. Vararattanavech, I. Plasencia, P. Greisen Jr, J. Bomholt, J. Torres, J. Emnéus and C. Hélix-Nielsen, Biochim. Biophys. Acta, Biomembr., 2011, 1808, 2600–2607.
22. W. R. P. Scott, P. H. Hünenberger, I. G. Tironi, A. E. Mark, S. R. Billeter, J. Fennen, A. E. Torda, T. Huber, P. Krüger and W. F. van Gunsteren, J. Phys. Chem. A, 1999, 103, 3596–3607.
23. S. J. Marrink, H. J. Risselada, S. Yefimov, D. P. Tieleman and A. H. de Vries, J. Phys. Chem. B, 2007, 111, 7812–7824.
24. K. J. Schweighofer, U. Essmann and M. Berkowitz, J. Phys. Chem. B, 1997, 101, 3793–3799.
25. O. Berger, O. Edholm and F. Jähnig, Biophys. J., 1997, 72, 2002–2013.
26. H. J. C. Berendsen, D. v. d. Spoel and R. v. Drunen, Comput. Phys. Commun., 1995, 91, 43–56.
27. M. J. Borgnia and P. Agre, Proc. Natl. Acad. Sci. U. S. A., 2001, 98, 2888–2893.
28. U. Essmann, L. Perera, M. L. Berkowitz, T. Darden, H. Lee and L. G. Pedersen, J. Chem. Phys., 1995, 103, 8577–8593.
29. B. Hess, H. Bekker, H. J. C. Berendsen and J. G. E. M. Fraaije, J. Comput. Chem., 1997, 18, 1463–1472.
30. A. L. Cheng and K. M. Merz, J. Phys. Chem., 1996, 100, 1927–1937.
31. M. Parrinello and A. Rahman, J. Appl. Phys., 1981, 52, 7182–7190.
32. C. Yang, D. N. Lu and Z. Liu, Biochemistry, 2011, 50, 2585–2593.
33. H. Lee and R. G. Larson, J. Phys. Chem. B, 2008, 112, 7778–7784.
34. M. G. Wolf, M. Hoefling, C. Aponte-Santamaria, H. Grubmuller and G. Groenhof, J. Comput. Chem., 2010, 31, 2169–2174.
35. O. S. Smart, J. M. Goodfellow and B. A. Wallace, Biophys. J., 1993, 65, 2455–2460.
36. K. Yurekli, C. A. Mitchell and R. Krishnamoorti, J. Am. Chem. Soc., 2004, 126, 9902–9903.
37. E. J. Wallace and M. S. P. Sansom, Nano Lett., 2007, 7, 1923–1928.
38. L. Xin, H. B. Su, C. H. Nielsen, C. Y. Tang, J. Torres and Y. G. Mu, Biochim. Biophys. Acta, Biomembr., 2011, 1808, 1581–1586.

---

Footnote

† Electronic supplementary information (ESI) available: The number of atoms (N) and the net charge (Z) of each AqpZ monomer at different pH values (Table S1). Molecule model for SDS. Force-field parameters for the bonded and non-bonded interactions of SDS molecules (Table S2 to S5). Gauche defect analysis of SDS tails (Figure S1). See DOI: 10.1039/c4ra03759f

---