Article pubs.acs.org/JPCB
Computation and Simulation of the Structural Characteristics of the Kidney Urea Transporter and Behaviors of Urea Transport Zhe Wang,† Tao Yu,†,‡ Jian-Ping Sang,*,†,‡ Xian-Wu Zou,† Chengfei Yan,§ and Xiaoqin Zou*,§ †
Department of Physics, Wuhan University, Wuhan 430072, China Department of Physics, Jianghan University, Wuhan 430056, China § Dalton Cardiovascular Research Center, Department of Physics and Astronomy, Department of Biochemistry, and Informatics Institute, University of Missouri, Columbia, Missouri 65211, United States ‡
S Supporting Information *
ABSTRACT: Urea transporters are a family of membrane proteins that transport urea molecules across cell membranes and play important roles in a variety of physiological processes. Although the crystal structure of bacterial urea channel dvUT has been solved, there lacks an understanding of the dynamics of urea transport in dvUT. In this study, by using molecular dynamics simulations, Monte Carlo methods, and the adaptive biasing force approach, we built the equilibrium structure of dvUT, calculated the variation in the free energy of urea, determined the urea-binding sites of dvUT, gained insight into the microscopic process of urea transport, and studied the water permeability in dvUT including the analysis of a water chain in the pore. The strategy used in this work can be applied to studying transport behaviors of other membrane proteins.
monomer contains 12 α-helices. Six of them, Pa and T1a− T5a, form a hemicylindrical domain. The other six helices, Pb and T1b−T5b, form the other hemicylindrical domain in an opposite orientation. Highly conserved residues from six different helices of the protein, Pa and Pb, T3a and T3b, and T5a and T5b, are brought together to form the ureapermeation pore.12 Water permeability in dvUT has never been studied experimentally, though water transport was measured in ApUT.14 Mammalian UT17 was also found to be water permeable, although the conclusion remained controversial.13 Despite big progress achieved in characterizing the structure and function of dvUT, several important questions remain unanswered: Is dvUT water permeable? How many ureabinding sites are there and where are they in dvUT? How is a urea molecule transported in the pore of dvUT? These fundamental questions will be addressed in the present study.
1. INTRODUCTION Urea, a common organic molecule, is the end-product of mammalian metabolism of nitrogenous compounds and is concentrated in the mammalian kidney.1 Urea is also a common source of nitrogen for bacteria, fungi, and plants.2 Specific species of bacteria catabolize urea into ammonia and CO2, which is useful for buffering acid loads.3 Although urea is uncharged, it has a strong dipole moment, which prohibits it from penetrating the nonpolar lipid membrane. Therefore, to transport urea molecules rapidly through the membrane, a special kind of integral membrane protein is required.4 Several specific urea transporter proteins have been found in prokaryotes and eukaryotes.5 UT-A (six isoforms, encoded by alternative splicing of the SLC14A2 gene) and UT-B (two isoforms, encoded by the SLC14A1 gene) are most common in mammalian, avian, amphibian, teleost, elasmobrach, and bacterial species.6 UT-C is observed only in two species of fish.7 Their homologues in bacteria include the H+-gated channel,3,8 the ABC-type urea active transporter,9 the iondependent urea transporter,10 the urea transporter-like ApUT (Actinobacillus pleuropneumoniae urea transporter),11 and dvUT (Desulfovibrio vulgaris urea transporter).12 Several experimental studies5,13−15 have characterized the structures and functions of a few mammalian and bacterial UTs to understand the permeation mechanisms of why and how urea transporters facilitate the diffusion of urea. These channellike proteins were found to facilitate the diffusion of urea and urea analogues along their concentration gradients at rates ranging from 104 to 106 s−1 without costing extra energy.13,16 The high-resolution crystal structure of dvUT revealed that the transmembrane protein is a homotrimer.12 Each dvUT © 2015 American Chemical Society
2. MODELS AND METHODS Construction of the dvUT System. The crystal structure of dvUT (PDB entry code: 3K3F)12 was used as the starting conformation. All the detergent molecules and heavy metal atoms were removed using VMD (visual molecular dynamics).18 The four missing residues 164−167 were modeled using Modeller.19−21 Hydrogen atoms were added using the PSFGEN module of the VMD package. The N-terminal end was capped with -NH3+, and the C-terminal end was capped Received: November 11, 2014 Revised: February 14, 2015 Published: March 17, 2015 5124
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 1. Equilibrium structure of the urea transporter dvUT. (A) Three-dimensional equilibrium structure of the dvUT monomer obtained by running MD simulations. The urea pore is colored green. The neighboring solution is colored white. (B) A view from the direction parallel to the plane containing the oxygen ladders. The oxygen atoms (On) can form hydrogen bonds with the amide hydrogen atoms (Ahn) of urea. The hydroxyl arms (Han) of T294 and T130 can form hydrogen bonds with the carbonyl oxygen of urea. (C) A view from the direction perpendicular to the plane containing the oxygen ladders. There are six hydrophobic side chains, which destabilize water and urea. The inset shows the carbonyl oxygen (Cox) and the amide hydrogen (Ahn) atoms of the urea molecule. Oxygen, red; carbon, cyan; nitrogen, blue; hydrogen, white.
with an -COO−, leaving a net charge of +e for the whole protein. To mimic the in vivo membrane environment, the protein was embedded into a 1-palmitoyl-2-oleoyl-phosphatidylethanolamine (POPE) lipid bilayer (90 Å × 90 Å), with the membrane aligned along the Z-axis. The lipid molecules that overlapped with the protein were removed. The system was then solvated using TIP3P water molecules and ionized with 200 mM NaCl by random substitution of water molecules with Na+ or Cl− ions. The dimension of the system was 90 Å × 90 Å × 80 Å. The whole system consisted of 59 176 atoms, in which 5060 atoms belonged to dvUT, 19 375 atoms to the membrane, 34 698 atoms to TIP3P water molecules, and 43 atoms to the ions. The force-field parameters for lipid, protein, ions, and water were taken from the CHARMM27 force field.22,23 Molecular Dynamic (MD) Simulation. MD simulations were performed for the whole system using the NAMD (NAnoscale Molecular Dynamics) software.24 Simulations were carried out in the NPT ensemble. Pressure was kept at 1 atm using the Langevin piston method.25,26 Temperature was kept at 298 K by coupling to a Langevin thermostat.26 The periodic boundary condition was applied to all the dimensions. Electrostatic interactions were treated using the Particle Mesh Ewald algorithm.27 van der Waals interactions and electrostatic interactions were truncated at 10 Å and smoothened at 8 Å.28 To equilibrate the initial dvUT channel−membrane− solution system, a series of simulations were performed. First, the initial structure was compacted by running energy minimizations for 5000 steps. Then the energy-minimized structure was simulated for 1.35 ns with fixed protein atoms to enable good contact between the protein and lipid molecules. Finally, the refined system was equilibrated for about 16 ns by MD simulations. The time step was set to 1 fs. The root-meansquare deviation (RMSD) of all the heavy atoms in the protein backbone as a function of MD simulation time is shown in Figure S1 of Supporting Information. Steered Molecular Dynamics (SMD) Simulation. Under normal conditions, the time scale for a urea molecule to cross the pore is 10−6 s to 10−4 s, but most MD runs simulate no more than tens of nanoseconds. To investigate urea transport by dvUT, we used the technique of constant-velocity steered molecular dynamics.29 In the SMD simulation, an external force was applied to the carbon atom of the urea molecule to make the urea travel quickly through the pore. The initial
configuration was generated by placing a urea molecule at the entrance of the equilibrium structure of dvUT (i.e., at around Z = −12.0 Å, see the Equilibrium Structure of dvUT and UreaPermeation Pore. To achieve equilibrium, the initial configuration was minimized using 5000 steps followed by 2 ns of MD simulations in the NPT ensemble. A harmonic constraint (with a spring constant of 5.0 kcal/mol·Å2) was attached to the carbon atom of the urea molecule; this carbon atom was pulled along the Z-axis at a constant velocity (10−2 Å/ps) so that the urea molecule was able to pass through the channel in ∼2.5 ns. The corresponding force fluctuation is shown in Figure S2 of Supporting Information. Adaptive Biasing Force (ABF) Calculations. The ABF method has been widely used in conjunction with molecular dynamics simulations to evaluate the free energy of binding for many biochemically interesting systems.30−33 Variation in the free energy of binding, ΔG(Z), along the Z-axis was calculated using the ABF method,34 which relies upon the integration of the average force acting on the urea molecule along the Zdirection. In the NAMD implementation of ABF,35 this force was evaluated within the classical thermodynamic integration formalism.36 The derivative of the free energy, dG(Z)/dZ, was estimated locally throughout the simulation, thereby providing a continuous update of the biasing force. The computation in the present study covered a distance ranging from −13 Å to 11 Å along the Z-axis, for which the origin (Z = 0) was located at the center of mass of the Cα atoms in the protein. In our calculations, the distance between the projected position of the carbon atom in the urea molecule along the Z-axis (i.e., the pore axis) and the center of mass of the Cα atoms in the protein was chosen as the reaction coordinate. To enhance the sampling efficiency, the full reaction coordinate along the Zaxis was divided into 24 nonoverlapping windows and the length of each window was set to 1 Å. For each window, at least 5 ns of MD trajectory were generated to obtain relatively uniform sampling and to reduce statistical errors. This protocol led to a total simulation time of at least 130 ns for the system. Instantaneous values of the force were accrued in bins of 0.1 Å width, which had been used to model transport phenomena.35 The resulting data were then integrated to produce the profile of the variation in the free energy between −13 Å to 11 Å. The starting conformation of the system for each window was taken from a SMD conformation in the corresponding window. 5125
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 2. (A) Variation in the free energy of urea along the axis of dvUT pore. The black arrows indicate the urea-resident sites, which can be considered as the actual urea-binding sites. The red arrows indicate the positions of the ladder oxygen atoms. The gray dashed lines show the swinging range of the hydroxyl arms (Ha1 of T294 and Ha2 of T130). The blue arrows represent the positions of the centers of mass of the hydrophobic residue pairs, including F243/F190, L84/L247, and F27/F80, respectively. (B) The occupation probability of urea along the pore axis (black solid line).
Monte Carlo (MC) Simulation. The MC simulation is a powerful tool for characterizing random walks of particles. To be able to apply MC simulations to urea transport, we need to validate that the urea movement in dvUT can be regarded as diffusion of particles. As shown in our MD simulations, as a urea molecule travels along the pore (Z-axis), the carbonyl group of the urea is parallel to the X-axis most of the time. Therefore, the movement of the urea molecule through the pore can be regarded as a translational motion. Namely, the urea molecule walks along the pore axis (Z-axis) with its carbonyl group parallel to the X-axis as the diffusive motion of a particle in the urea pore. The specific protocol for our MC simulations is described as follows. The urea pore, i.e., the walking space, was divided into a grid of size 0.57 Å × 0.57 Å × 0.57 Å. The orifices of the pore at Z = −11 and 11 Å, respectively, corresponded to level 0 on the extracellular side and level 38 on the intracellular side of the grid, respectively. The walker moved on the grid from site i to its neighboring site j with the jumping probability Pij ∝ exp( −ΔGij/NAkBT )
To quantitatively characterize the pore and to perform permeation simulations, we need to define the locations of the two orifices of the urea pore. For convenience, the mass center of Cα atoms in the protein was set as the origin of the coordinate system (Z = 0). It can be seen from Figure 1A that the cross section of the pore increases dramatically at Z = −11 and 11 Å, where the pore becomes two wide vestibules. Therefore, the cross section at Z = −11 Å (11 Å) was defined as the orifice of the pore to extracellular (intracellular) solution. The more distant from the orifices, the narrower the pore. In the middle region, the pore is so narrow that only one urea molecule can pass through. Consequently, outside the pore the urea molecules are fully hydrated. They are partially dehydrated when they enter the pore and become fully dehydrated when they are in the middle of the pore. Similar dehydration behavior has been reported in detail for the transport of urea molecules in UT-B.39 We further examined our equilibrated structure of dvUT for the detailed composition of the urea pore. Figure 1B and 1C shows the stereoview of the residues lining the pore in two perpendicular directions. As shown in Figure 1B, along one side of the pore, E187 and V188 from the helix Pb and V25 and Q24 from the helix Pa form two linear arrays of three oxygen atoms, which are referred to as oxygen ladders. The ladder atoms O1 and O6 are the side-chain oxygen atoms of E187 and Q24, respectively. The ladder atoms O2, O3, O4, and O5 are the backbone oxygen atoms of E187, V188, V25, and Q24, respectively. These oxygen atoms can form hydrogen bonds with the amide hydrogen atoms (Ahn) of the urea molecule. O1, O2, and O3 form one oxygen ladder, and O4, O5, and O6 form the other ladder. Within each oxygen ladder, the distances between the adjacent oxygen atoms are 3.0−4.0 Å. Therefore, urea can move successively from one oxygen site to the next along the same oxygen ladder by forming and breaking hydrogen bonds. On the opposite side of the oxygen ladders on the pore, L293 and T294 from the helix T5b and T130 and L129 from the helix T5a provide hydroxyl or amide hydrogens that can potentially form hydrogen bonds with the carbonyl oxygen atoms of urea (denoted as Cox). Intriguingly, the distance between the two oxygen ladders is large, i.e., Z(O4) − Z(O3) = 7 Å. It is therefore difficult to move a urea molecule
(1)
where ΔGij = Gj − Gi, and Gi and Gj are the free energies of the walker positioning at the jumping-out-site i and jumping-ontosite j, respectively. NA is the Avogadro constant, and kB is the Boltzmann constant. In the present case, the absolute value of the free energy of the urea molecule in dvUT is unknown; ΔGij is therefore estimated as the difference between the variations in the free energy at site i and site j.
3. RESULTS AND DISCUSSION Equilibrium Structure of dvUT and Urea-Permeation Pore. Constructing an equilibrium conformation of dvUT is a requisite for probing the mechanisms underlying the dvUT functions at the atomic level. Based on the crystal structure of dvUT,12 an equilibrated dvUT channel−membrane−solution system was achieved by running molecular dynamics simulations for 16 ns. The equilibrium structure of dvUT is shown in Figure 1. To search and view the urea-permeation pore, the HOLE software was used.37,38 The identified urea pore is plotted in Figure 1A as a green tube. 5126
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 3. Snapshots describing the transport process of urea in the pore. Zu represents the Z-coordinate of the carbon atom in the urea molecule. The hydrogen bonds are represented by dashed lines. (A) A urea molecule enters the entrance of the pore from the extracellular side. (B) At the external resident site (Sext), the urea forms hydrogen bonds with O3, water, and T294. (C) The hydroxyl arm of T294, Ha1, drags the urea to the middle region of the pore. (D) The urea is located at the unstable equilibrium site near the center of the pore. (E) In the transitional configuration, the hydroxyl arm of T130, Ha2, is dragging the urea to O4. (F) The process of dragging the urea from the middle region to O4 has been completed. (G) At the middle resident site (Smid), urea forms a hydrogen bond with O4 and can be bridged with O5 and/or O6 via water. (H) At the internal resident site (Sint), the urea forms a hydrogen bond with O6 and can interact with the amide hydrogen of L129 via a bridging water molecule.
regions and allow urea molecules to pass smoothly through the urea pore.12,39 Variation of the Free Energy, Occupation Probability Profile, and Urea-Resident Sites. With the equilibrated structure of dvUT obtained from our MD simulations, we were able to characterize the transport properties of dvUT in the present study. We started with investigating the urea transport in dvUT. To find the most stable resident positions of urea in dvUT and then to investigate the dynamics of urea transport, we calculated the variation in the free energy of urea, ΔG(Z), along the pore axis Z using the adaptive biasing force (ABF) approach.34 Figure 2A shows that the curve of the variation in the free energy vs Z exhibits three minima (energy valleys) and multiple maxima (energy barriers). The first three maxima appear in the region where the hydrophobic residues are located. The minima occur mostly in proximity to the positions of the oxygen atoms in the ladders. In the middle region of the pore, there exists a high, wide, and steep barrier. On the basis of the calculated variation in the free energy, we determined the urea-resident sites in dvUT by utilizing the Monte Carlo (MC) method. In the MC simulations, the number of MC-time (i.e., the residence time) at each site of the pore was determined. The details are provided in the Supporting Information SI1. Then, the occupation and cumulative probability of urea vs the Z-coordinate along the whole pore were calculated and plotted in Figure 2B. The sites with the maximum occupation probability values were defined as the urea-resident sites. On the urea-resident sites, urea molecules bound strongly to the pore, and therefore the urea-
from one oxygen ladder to another. Moreover, the plane that contains O1, O2, and O3 forms an angle with the plane that contains O4, O5, and O6, which makes the movement of urea from one oxygen ladder to another more difficult. This difficulty is overcome by using the hydroxyl arms of T294 and T130 as a shuttle. Specifically, the hydroxyl arms (Han) of threonine residues T294 and T130 can form hydrogen bonds with the carbonyl oxygen of urea. The arms are flexible enough to follow the movement of urea as shown in our MD simulations (see Transport Process of Urea in the Pore of dvUT). Thus, the hydroxyl arms of T294 and T130 can facilitate the transition of the urea molecule from one oxygen ladder to another. Different from Figure 1B, Figure 1C shows the equilibrium structure of dvUT from the direction that is perpendicular to the plane containing the oxygen ladders. It can be seen from Figure 1C that there are six hydrophobic side chains. Among them, the phenylalanine residues, F190 and F243 from the periplasmic side and F27 and F80 from the cytoplasmic side, compress the pore into a slot-like shape. The other two hydrophobic residues, L247 and L84 positioning in the middle of the pore, also destabilize polar molecules such as water and urea in this region. In summary, our equilibrium structure of dvUT and MD simulations show that the urea pore can be considered as a composite configuration imposed by the hydrogen bonding pattern with the oxygen ladder in an otherwise hydrophobic pore. These hydrophobic regions coexist with the hydrophilic 5127
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
and the water. The hydroxyl arm, Ha2, of T130 drags the urea molecule to O4 (Figure 3E, a transitional configuration). The urea molecule is then pushed forward by the hydroxyl arm of T294, until the arm Ha2 swings to as far as it can reach (Figure 3F). As the urea molecule moves forward, the hydrogen bond between the urea and T130 breaks. Besides retaining the hydrogen bond between the urea and O4, the urea molecule may connect to O5 and/or O6 by using water molecules as a bridge (Figure 3G, Smid). In this region, the urea molecule appears to have a long dwelling time. The location with the maximum urea-occupation probability is the middle resident site. As the urea molecule crosses the energy barrier that is formed by the hydrophobic residue pair F27/F80 and reaches O6, it may be bridged to the amide hydrogen of L129 via water molecules (Figure 3H, Sint). These hydrogen bonding interactions significantly increase the dwelling time of the urea molecule, and this position becomes the internal resident site. Because O6 is close to the orifice of the pore, once the hydrogen bond between the urea and O6 dissociates, the urea molecule enters the intracellular solution. In summary, our simulations of the transport process of urea along the pore uncover the structural characteristics of the ureabinding sites. Specifically, at the binding sites, multiple hydrogen bonds are formed, particularly, between the urea molecule and dvUT. In addition, all three binding sites of urea are distant from the energy barriers resulting from highly hydrophobic residues (see Figure 2A). Our simulations of the urea transport can also explain why the energy barrier in the middle of the pore is exceedingly high and wide as seen from the calculated profile of the free-energy variation in dvUT (Figure 2A). As shown in this figure, the peak of this energy barrier locates at the position of the hydrophobic residue pair, L247 and L84. The barrier spans a wide range from −3.5 Å to 2.0 Å, which is the range that the hydroxyl arms of T294 and T130 can swing (as indicated by the two vertical dashed lines in Figure 2A). Within this range, the inversion of the hydroxyl arms of T294 and T130 drags the urea molecule across this energy barrier from one oxygen ladder to the other, as shown in Figures 3C−E. Moreover, our simulations also show that in the middle region, the pore is narrow and the urea molecule is completely dehydrated, whereas in the other regions the urea molecule is partially dehydrated. In summary, the middle energy barrier can be attributed to the repulsive effect of the hydrophobic residues L247 and L84, the energetic cost for inverting the hydroxyl arms of T294 and T130, and the desolvation cost. The width of the middle energy barrier coincides with the distance by which the urea molecule is dragged by the hydroxyl arms of T294 and T130. Due to the resulting high energy barrier, it is difficult for the urea molecule to enter the middle region and easy for the urea to leave this region, and therefore it is not surprising that the ureaoccupation probability in the middle region is tiny. Water Permeability of dvUT. Next, we studied the water permeability in dvUT by MD simulations. We started with resolving the puzzle whether dvUT is permeable to water molecules. In the crystal structure of dvUT no water molecule was found.12 However, water molecules were observed in Actinobacillus pleuropneumoniae urea transporter ApUT.14 Recent studies also suggested that urea, H2O, CO2, and NH3 pass through urea transporter B (UT-B).40 In our MD simulations, we observed that a number of water molecules were able to move into or out of the urea pore of dvUT. Specifically, within 20 ns, six water molecules passed through
resident sites are expected to be the actual urea-binding sites. Figure 2B shows that in dvUT the urea-resident sites locate approximately at Z = −5.0 Å, 4.0 Å, and 9.0 Å, which are defined as the external (Sext), middle (Smid), and internal (Sint) resident sites, respectively. The characteristics and underlying mechanisms of these resident sites will be elucidated in the next subsection. Transport Process of Urea in the Pore of dvUT. With the identified urea-resident sites, we characterized the details of the transport process of urea in the pore of dvUT using constant-velocity Steered MD (cv-SMD).29 Specifically, in the presence of the concentration gradient of urea across the membrane, urea molecules diffused along the pore driven by both thermal fluctuations and the concentration gradient. The diffusion of urea through dvUT was simulated using cv-SMD, and the results were analyzed in detail. Our MD simulations showed that the urea molecules that move through the pore are always accompanied by water molecules. Figure 3 shows representative snapshots of the diffusion of urea along the pore to elucidate the mechanisms of urea transport. Specifically, the coordinate of the carbon atom of a urea molecule, CO(NH2)2, in the Z-direction is denoted by Zu. When a urea molecule enters the pore from the extracellular side, an amide hydrogen atom, Ahn, of urea forms a hydrogen bond with the ladder oxygen atom O1 (Figure 3A) or O2. Along the oxygen ladder the urea molecule moves from O1 to O2, crosses the energy barrier that is created by the hydrophobic residue pair, F243 and F190, and then reaches O3. The urea molecule bound to O3 is located at the valley between the energy barrier of F243/F190 and the middle energy barrier. In this region, the amide hydrogen atoms, Ahn, of the urea molecule participate in hydrogen bonding with the ladder oxygen O3 and the hydroxyl group of water, respectively. When the urea molecule reaches the site at Zu = −4.9 Å, the hydroxyl arm of T294, Ha1, forms a hydrogen bond with the carbonyl oxygen, Cox, of the urea molecule (Figure 3B). This threehydrogen-bond configuration stabilizes the urea molecule in the pore, and the position becomes the external resident site (Sext). If the energy provided by thermal fluctuations is large enough to switch the hydroxyl arm of T294, the hydrogen bond between the urea and the water breaks, and the hydroxyl arm Ha1 of T294 drags the urea molecule to the middle region of the pore (Figure 3C). In the middle region, the hydroxyl arm of T130, Ha2, forms a hydrogen bond with the carbonyl oxygen, Cox, of the urea molecule. When the urea molecule reaches the center of the middle region (Figure 3D), the hydrogen bonds, formed by the urea with T294 and T130, are symmetrically arranged on both sides of the urea molecule. In addition, on the opposite side of residues T294 and T130, the hydrophobic side chains of V25 and V188 dislike the urea molecule. The competition between the hydrogen bonding from T294 and T130 and the hydrophobic effect from V25 and V188 destabilizes this equilibrium position, resulting in the peak of the middle energy barrier in Figure 2A. The urea molecule then leaves the peak position of the middle energy barrier and forms a hydrogen bond with O4. Hydrogen bonds are formed between the urea and O4, urea and T294, and urea and T130. This three-hydrogen-bond configuration results in the metastable equilibrium site at Zu = 0.1 Å (the configuration is not shown), which corresponds to the minimum of the middle energy barrier (see Figure 2A). As the urea molecule leaves the middle region, the hydrogen bond between the urea and T294 breaks, and a new hydrogen bond is formed between the urea 5128
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 4. Water molecules entering the urea pore of dvUT. (A) Water molecules are dispersed and partially bonded together via hydrogen bonds. (B) Water molecules form a continuous water chain. The water molecules in the water chain on the two sides of the central water molecule exhibit opposite dipole orientations. The hydrogen bonds are represented by dashed lines. The residues of the protein are depicted in licorice representation, and the water molecules are displayed in CPK representation. (C) The panoramic view of a water chain connecting the extracellular and intracellular solution. Residues of the dvUT monomer are represented as cyan cylinders. The lipid bilayer is removed for a better view. The oxygen and hydrogen atoms of water molecules are represented by red and white van der Waals balls, respectively.
the pore from the extracellular solution to the intracellular solution, and seven water molecules moved in the opposite direction by diffusion. The time taken by these water molecules to cross the pore ranged between 0.74 ns and 6.5 ns, and the average rate was about 0.4 H2O/ns/channel. Then what is the arrangement of these water molecules in the dvUT pore? Our MD simulations showed that in most instances, water molecules were dispersed and partially bound together via hydrogen bonds in the urea pore (Figure 4A), and in some cases, a continuous water chain was formed in the pore (Figure 4B). It can be seen from Figure 4B that in the water chain the oxygen atom of the water molecule in the middle of the pore forms two hydrogen bonds with the hydroxyls of T294 and T130. This central water molecule also forms hydrogen bonds with adjacent water molecules. For the remaining water molecules in the water chain, most of them are able to bind to the oxygen ladders. Interestingly, the water molecules on the two sides of the central water molecule exhibit opposite dipole orientations. For the water molecules arranged in a variety of configurations, we can describe their distribution using the probability density. To calculate the probability density, we took 3706 configurations for statistical analysis. First, we divided the pore into small segments of 0.5 Å in length and counted the number of the water molecules per segment for each configuration. Then we computed the average number of the water molecules in all configurations for each segment and obtained the distribution of the average number of water molecules along the pore (Z axis). Normalization of the distribution of the average number of the water molecule yielded the probability distribution of the water molecules in the dvUT pore. Dividing the probability distribution by the length of segment, we achieved the distribution of probability density of the water molecules in the pore (see Figure 5). To inspect how the structure of the dvUT pore affects the probability density distribution, variation in the cross-sectional area of pore along the Z axis is also plotted in Figure 5. The cross-sectional area was obtained by averaging the measurements from ten configurations of the dvUT pore. It can be seen from Figure 5 that the probability density of the water molecules and the cross-sectional area of the dvUT pore essentially vary in the same tendency, indicating that the distribution of the water molecules is mostly affected by the
Figure 5. Probability density of the water molecules in the dvUT pore and cross-sectional area of the dvUT pore along the Z axis. The black arrows indicate the positions of the ladder oxygen atoms. The blue arrows represent the positions of the centers of mass of the hydrophobic residue pairs. The green dotted lines indicate the positions of Z = −9.0, −1.0, and 7.0 Å.
geometry of the pore. It is well-known that water has a weaker dipole moment (1.85 D) and smaller length (1.5 Å) than urea (4.56 D and 3.9 Å);12 therefore, the influence of the hydrophobic residues and ladder oxygen atoms on water is less than their influence on urea. In contrast, the geometrical characteristic of the pore may have an important effect on the behavior of the water molecules. For example, within the confined range (−9.0 Å < Z < 7.0 Å), the width of the pore ranges between 3.0 and 4.5 Å, which is comparable to the size of a water molecule. As shown in Figure 5, in the confined range, the probability to find water molecules is much lower than that outside this range. In addition to geometry confinement, the chemical composition of the protein also plays a role. Under the influence of the hydrophobic residues, the probability density curve exhibits minimums (around L84/ L247 and F27/F80) or shows an apparent reduction (around F243/F190). On the other hand, the attraction of the ladder oxygen atoms causes a slight increase of the probability density. To provide a microscopic description of the distribution of water molecules in dvUT, we examined the behaviors of the water molecules during the MD simulation. It was found that in the middle of the pore (−2.0 Å < Z < 0 Å), where the pore is 5129
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
midpoint, there are more surrounding water molecules, and the interactions between the adjacent water dipole moments become stronger. Consequently, beyond the range of −4.0 Å < Z < 2.0 Å, the behavior of the orientations of the water dipole moments become more complex. Beyond the confined range, −9.0 Å < Z < 7.0 Å, the number of water molecules increases significantly, which greatly decreases the influence of the dvUT protein on the orientations of the dipole moments of the water molecules. Outside the pore the orientations of the water dipole moments are arbitrary, and ⟨cos θ⟩ approaches 0. Because of computational difficulty, we averaged only with 3706 configurations, and therefore outside the pore the ⟨cos θ⟩ curve did not become horizontal.
very narrow and the influence of the hydrophobic residues L84/L247 is strong, the average number of the water molecules per channel was about 0.4 and the average spacing between adjacent water molecules was larger than 5 Å. In other words, for most of the time there was no water molecule. The probability of finding one water molecule was about one-third, and the probability of finding two water molecules was about 5% in this range. It was a rare event for a water chain to be formed in the dvUT pore. Finally, the polarization behavior of the water molecules in the dvUT pore was also investigated. Specifically, we calculated the projection of the dipole moment of each water molecule in the water chain along the axis of the pore and plotted the curve in Figure 6. It can be seen from this figure that at the midpoint
4. CONCLUSIONS Using MD and MC simulation methods, we studied the structural characteristics and transport behaviors of dvUT. From binding free energy calculations we determined the positions of three urea binding sites, which are located at Z = −5.0 Å, 4.0 Å, and 9.0 Å for the external (Sext), middle (Smid), and internal (Sint) binding sites, respectively. These binding sites are separated by the energy barriers, and the sites participate in multiple hydrogen bonding interactions with the urea molecule. Our MD simulations suggested that dvUT is water permeable, and the average rate is about 0.4 H2O/ns/ channel for the water molecules to cross the dvUT pore. The movement of urea is always accompanied by water. We calculated the probability density of the water molecules in the dvUT pore and the projection of the dipole moments of the water molecules along the axis of the pore. We found that the pseudo-2-fold symmetry of the wall of the pore was the reason for the water molecules on the two sides of the midpoint of the pore to have opposite dipole moment orientations.
Figure 6. Projection of the dipole moments of the water molecules along the axis of the dvUT pore. θ is the angle between the dipole moment and the pore axis. cos θ is the projection of the dipole moment along the pore axis for a water molecule in the water chain when the dipole moment of the water molecule is set as a unit. ⟨cos θ⟩ is the statistical average of cos θ for all water molecules occurring at the same section of the pore when the water molecules are in a dispersed and partially bonded state. The magenta dashed line represents the cross-sectional area of the pore along the pore axis.
■
ASSOCIATED CONTENT
S Supporting Information *
Additional figures on MD and SMD simulations. This material is available free of charge via the Internet at http://pubs.acs.org.
■
of the pore, Z = −1.0 Å, both cos θ and ⟨cos θ⟩ are equal to zero, which indicates that the dipole moment is perpendicular to the pore axis. On the two sides of Z = −1.0 Å, the water molecules have opposite dipole moment orientations. The behavior of the orientation of the water dipole moments is determined by the structure and composition of the dvUT protein. It can be seen from the structure of dvUT (see Figure 1) that the nature and positioning of the side chains forming the wall of the pore are remarkably symmetrical: for every residue along the wall on one side of the midpoint, there exists an identical residue on the other side of the midpoint at the symmetrical location (except the Q24/E187 pair).12 Namely, the wall has a pseudo-2-fold symmetry, and the axis of the symmetry is perpendicular to the axis of the pore, passing through the midpoint of the pore and lying in the plane containing the oxygen ladders. This pseudo-2-fold symmetry makes the water dipole moment perpendicular to the axis of the pore at Z = −1.0 Å, and makes the water molecules have the opposite dipole moment orientations on the two sides of Z = −1.0 Å. By using the hydrogen bond as a bridge a water molecule have an influence not only on the nearest neighboring water molecule but also on the next nearest neighboring water molecules. When the water molecule is moving away from the
AUTHOR INFORMATION
Corresponding Authors
*Tel: (573) 882-6045. E-mail: [email protected]. *Tel: (86) 27-68766665. E-mail: [email protected]. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS We are very grateful to Dr. M. Zhou for helpful discussions. We also thank Benjamin Merideth very much for editing the manuscript. A portion of the computations were performed on the HPC resources at the University of Missouri Bioinformatics Consortium (UMBC). This work was supported by the US National Science Foundation Career Award (Grant DBI0953839) and National Natural Science Foundation of China (grants 11304123 and 11374234).
■
REFERENCES
(1) Nielsen, S.; Terris, J.; Smith, C. P.; Hediger, M. A.; Ecelbarger, C. A.; Knepper, M. A. Cellular and Subcellular Localization of the Vasopressin- Regulated Urea Transporter in Rat Kidney. Proc. Natl. Acad. Sci. U. S. A. 1996, 93, 5495−5500.
5130
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781− 1802. (25) Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant Pressure Molecular Dynamcis Algorithms. J. Chem. Phys. 1994, 101, 4177− 4189. (26) Feller, S. E.; Zhang, Y. H.; Pastor, R. W.; Brooks, B. R. Constant-pressure Molecular Dynamics Simulation- the Langevin Piston Method. J. Chem. Phys. 1995, 103, 4613−4621. (27) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−8593. (28) Miyamoto, S.; Kollman, P. A. SETTLE: An Analytical Version of the SHAKE and RATTLE Algorithm for Rigid Water Molecules. J. Comput. Chem. 1992, 13, 952−962. (29) Izrailev, S.; Stepaniants, S.; Isralewitz, B.; Kosztin, D.; Lu, H.; Molnar, F.; Wriggers, W.; Schulten, K. Steered Molecular Dynamic. In Computational Molecular Dynamics: Challenges, Methods, Ideas; Deuflhard, P.; Hermans, J.; Leimkuhler, B.; Mark, A. E.; Reich, S.; Skeel, R. D, Eds.; Springer-Verlag: Berlin, 1998; Vol. 4 of Lecture Notes in Computational Science and Engineering, pp 39−65. (30) Dehez, F.; Pebay-Peyroula, E.; Chipot, C. Binding of ADP in the Mitochondrial ADP/ATP Carrier Is Driven by an Electrostatic Funnel. J. Am. Chem. Soc. 2008, 130, 12725−12733. (31) Hénin, J.; Tajkhorshid, E.; Schulten, K.; Chipot, C. Diffusion of Glycerol Through Escherichia Coli Aquaglyceroporin GlpF. Biophys. J. 2008, 94, 832−839. (32) Gamini, R.; Sotomayor, M.; Chipot, C.; Schulten, K. Cytoplasmic Domain Filter Function in the Mechanosensitive Channel of Small Conductance. Biophys. J. 2011, 101, 80−89. (33) Modi, N.; Bárcena-Uribarri, I.; Bains, M.; Benz, R.; Hancock, R. E. W.; Kleinekathöfer, U. Role of the Central Arginine R133 Toward the Ion Selectivity of the Phosphate Specific Channel OprP: Effects of Charge and Solvation. Biochemistry 2013, 52, 5522−5532. (34) Darve, E.; Pohorille, A. Calculating Free Energies Using Average Force. J. Chem. Phys. 2001, 115, 9169−9183. (35) Hénin, J.; Chipot, C. Overcoming Free Energy Barriers Using Unconstrained Molecular Dynamics Simulations. J. Chem. Phys. 2004, 121, 2904−2914. (36) den Otter, W. K. Thermodynamic Integration of the Free Energy Along a Reaction Coordinate in Cartesian Coordinates. J. Chem. Phys. 2000, 112, 7283−7292. (37) Smart, O. S.; Neduvelil, J. G.; Wang, X.; Wallace, B. A.; Sansom, M. S. HOLE: a Program for the Analysis of the Pore Dimensions of Ion Channel Structural Models. J. Mol. Graphics 1996, 14, 354−360. (38) Smart, O. S.; Goodfellow, J. M.; Wallace, B. A. The Pore Dimensions of Gramicidin A. Biophys. J. 1993, 65, 2455−2460. (39) Levin, E. J.; Cao, Y.; Enkavi, G.; Quick, M.; Pan, Y.; Tajkhorshid, E.; Zhou, M. Structure and Permeation Mechanism of a Mammalian Urea Transporter. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 11194−11199. (40) Geyer, R. R.; Musa-Aziz, R.; Enkavi, G.; Mahinthichaichan, P.; Tajkhorshid, E.; Boron, W. F. Movement of NH3 Through the Human Urea Transporter B: A New Gas Channel. Am. J. Physiol. Renal Physiol. 2013, 304, F1447−F1457.
(2) Abreu, C.; Sanguinetti, M.; Amillis, S.; Ramon, A. UreA, the Major urea/H+ Symporter in Aspergillus Nidulans. Fungal Genet. Biol. 2010, 47, 1023−1033. (3) Strugatsky, D.; McNulty, R.; Munson, K.; Chen, C.-K.; Soltis, S. M.; Sachs, G.; Luecke, H. Structure of the Proton-gated Urea Channel from the Gastric Pathogen Helicobacter Pylori. Nature 2013, 493, 255−258. (4) Knepper, M. A.; Mindell, J. A. Molecular Coin Slots for Urea. Nature 2009, 462, 733−734. (5) You, G.; Smith, C. P.; Kanal, Y.; Lee, W.-S.; Stelzner, M.; Hedlger, M. A. Cloning and Characterization of the Vasopressinregulated Urea Transproter. Nature 1993, 365, 844−847. (6) Smith, C. P. Mammalian Urea Transporters. Exp. Physiol. 2009, 94, 180−185. (7) Mistry, A. C.; Chen, G.; Kato, A.; Nag, K.; Sands, J. M.; Hirose, S.; Chandra, A. A Novel Type of Urea Transporter, UT-C, Is Highly Expressed in Proximal Tubule of Seawater Eel Kidney. Am. J. Physiol. Renal Physiol. 2005, 288, F455−F465. (8) Weeks, D. L.; Eskandari, S.; Scott, D. R.; Sachs, G. A H+-Gated Urea Channel: The Link Between Helicobacter Pylori Urease and Gastric Colonization. Science 2000, 287, 482−485. (9) Valladares, A.; Montesinos, M. L.; Herrero, A.; Flores, E. An ABC-type, High-affinity Urea Permease Identified in Cyanobacteria. Mol. Microbiol. 2002, 43, 703−715. (10) Kojima, S.; Bohner, a.; Von Wirén, N. Molecular Mechanisms of Urea Transport in Plants. J. Membr. Biol. 2006, 212, 83−91. (11) Godara, G.; Smith, C.; Bosse, J.; Zeidel, M.; Mathai, J. Functional Characterization of Actinobacillus Pleuropneumoniae Urea Transport Protein, ApUT. Am. J. Physiol. 2009, 296, R1268−73. (12) Levin, E. J.; Quick, M.; Zhou, M. Crystal Structure of a Bacterial Homologue of the Kidney Urea Transporter. Nature 2009, 462, 757− 762. (13) Maciver, B.; Smith, C. P.; Hill, W. G.; Zeidel, M. L. Functional Characterization of Mouse Urea Transporters UT-A2 and UT-A3 Expressed in Purified Xenopus Laevis Oocyte Plasma Membranes. Am. J. Physiol. Renal Physiol. 2008, 294, F956−F964. (14) Raunser, S.; Mathai, J. C.; Abeyrathne, P. D.; Rice, A. J.; Zeidel, M. L.; Walz, T. Oligomeric Structure and Functional Characterization of the Urea Transporter from Actinobacillus Pleuropneumoniae. J. Mol. Biol. 2009, 387, 619−627. (15) Zhao, D.; Sonawane, N. D.; Levin, M. H.; Yang, B. Comparative Transport Efficiencies of Urea Analogues Through Urea Transporter UT-B. Biochim. Biophys. Acta 2007, 1768, 1815−1821. (16) Mannuzzu, L. M.; Moronne, M. M.; Macey, R. I. Estimate of the Number of Urea Transport Sites in Erythrocyte Ghosts Using a Hydrophobic Mercurial. J. Membr. Biol. 1993, 133, 85−97. (17) Yang, B.; Verkman, A. S. Analysis of Double Knockout Mice Lacking Aquaporin-1 and Urea Transporter UT-B. Evidence for UT-Bfacilitated Water Transport in Erythrocytes. J. Biol. Chem. 2002, 277, 36782−36786. (18) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (19) Fiser, A.; Sali, A. Modeller: Generation and Refinement of Homology-based Protein Structure Models. Methods Enzymol. 2003, 374, 461−491. (20) Martí-Renom, M. A.; Stuart, A. C.; Fiser, A.; Sánchez, R.; Melo, F.; Sali, A. Comparative Protein Structure Modeling of Genes and Genomes. Annu. Rev. Biophys. Biomol. Struct. 2000, 29, 291−325. (21) Sali, A.; Blundell, T. L. Comparative Protein Modelling by Satisfaction of Spatial Restraints. J. Mol. Biol. 1993, 234, 779−815. (22) MacKerell, A. D.; Bashford, D.; Bellott, M.; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; et al. All-atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586−3616. (23) Feller, S. E.; MacKerell, A. D. An Improved Empirical Potential Energy Function for Molecular Simulations of Phospholipids. J. Phys. Chem. B 2000, 104, 7510−7515. (24) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kalé, L.; Schulten, K. Scalable 5131
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Computation and Simulation of the Structural Characteristics of the Kidney Urea Transporter and Behaviors of Urea Transport Zhe Wang,† Tao Yu,†,‡ Jian-Ping Sang,*,†,‡ Xian-Wu Zou,† Chengfei Yan,§ and Xiaoqin Zou*,§ †
Department of Physics, Wuhan University, Wuhan 430072, China Department of Physics, Jianghan University, Wuhan 430056, China § Dalton Cardiovascular Research Center, Department of Physics and Astronomy, Department of Biochemistry, and Informatics Institute, University of Missouri, Columbia, Missouri 65211, United States ‡
S Supporting Information *
ABSTRACT: Urea transporters are a family of membrane proteins that transport urea molecules across cell membranes and play important roles in a variety of physiological processes. Although the crystal structure of bacterial urea channel dvUT has been solved, there lacks an understanding of the dynamics of urea transport in dvUT. In this study, by using molecular dynamics simulations, Monte Carlo methods, and the adaptive biasing force approach, we built the equilibrium structure of dvUT, calculated the variation in the free energy of urea, determined the urea-binding sites of dvUT, gained insight into the microscopic process of urea transport, and studied the water permeability in dvUT including the analysis of a water chain in the pore. The strategy used in this work can be applied to studying transport behaviors of other membrane proteins.
monomer contains 12 α-helices. Six of them, Pa and T1a− T5a, form a hemicylindrical domain. The other six helices, Pb and T1b−T5b, form the other hemicylindrical domain in an opposite orientation. Highly conserved residues from six different helices of the protein, Pa and Pb, T3a and T3b, and T5a and T5b, are brought together to form the ureapermeation pore.12 Water permeability in dvUT has never been studied experimentally, though water transport was measured in ApUT.14 Mammalian UT17 was also found to be water permeable, although the conclusion remained controversial.13 Despite big progress achieved in characterizing the structure and function of dvUT, several important questions remain unanswered: Is dvUT water permeable? How many ureabinding sites are there and where are they in dvUT? How is a urea molecule transported in the pore of dvUT? These fundamental questions will be addressed in the present study.
1. INTRODUCTION Urea, a common organic molecule, is the end-product of mammalian metabolism of nitrogenous compounds and is concentrated in the mammalian kidney.1 Urea is also a common source of nitrogen for bacteria, fungi, and plants.2 Specific species of bacteria catabolize urea into ammonia and CO2, which is useful for buffering acid loads.3 Although urea is uncharged, it has a strong dipole moment, which prohibits it from penetrating the nonpolar lipid membrane. Therefore, to transport urea molecules rapidly through the membrane, a special kind of integral membrane protein is required.4 Several specific urea transporter proteins have been found in prokaryotes and eukaryotes.5 UT-A (six isoforms, encoded by alternative splicing of the SLC14A2 gene) and UT-B (two isoforms, encoded by the SLC14A1 gene) are most common in mammalian, avian, amphibian, teleost, elasmobrach, and bacterial species.6 UT-C is observed only in two species of fish.7 Their homologues in bacteria include the H+-gated channel,3,8 the ABC-type urea active transporter,9 the iondependent urea transporter,10 the urea transporter-like ApUT (Actinobacillus pleuropneumoniae urea transporter),11 and dvUT (Desulfovibrio vulgaris urea transporter).12 Several experimental studies5,13−15 have characterized the structures and functions of a few mammalian and bacterial UTs to understand the permeation mechanisms of why and how urea transporters facilitate the diffusion of urea. These channellike proteins were found to facilitate the diffusion of urea and urea analogues along their concentration gradients at rates ranging from 104 to 106 s−1 without costing extra energy.13,16 The high-resolution crystal structure of dvUT revealed that the transmembrane protein is a homotrimer.12 Each dvUT © 2015 American Chemical Society
2. MODELS AND METHODS Construction of the dvUT System. The crystal structure of dvUT (PDB entry code: 3K3F)12 was used as the starting conformation. All the detergent molecules and heavy metal atoms were removed using VMD (visual molecular dynamics).18 The four missing residues 164−167 were modeled using Modeller.19−21 Hydrogen atoms were added using the PSFGEN module of the VMD package. The N-terminal end was capped with -NH3+, and the C-terminal end was capped Received: November 11, 2014 Revised: February 14, 2015 Published: March 17, 2015 5124
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 1. Equilibrium structure of the urea transporter dvUT. (A) Three-dimensional equilibrium structure of the dvUT monomer obtained by running MD simulations. The urea pore is colored green. The neighboring solution is colored white. (B) A view from the direction parallel to the plane containing the oxygen ladders. The oxygen atoms (On) can form hydrogen bonds with the amide hydrogen atoms (Ahn) of urea. The hydroxyl arms (Han) of T294 and T130 can form hydrogen bonds with the carbonyl oxygen of urea. (C) A view from the direction perpendicular to the plane containing the oxygen ladders. There are six hydrophobic side chains, which destabilize water and urea. The inset shows the carbonyl oxygen (Cox) and the amide hydrogen (Ahn) atoms of the urea molecule. Oxygen, red; carbon, cyan; nitrogen, blue; hydrogen, white.
with an -COO−, leaving a net charge of +e for the whole protein. To mimic the in vivo membrane environment, the protein was embedded into a 1-palmitoyl-2-oleoyl-phosphatidylethanolamine (POPE) lipid bilayer (90 Å × 90 Å), with the membrane aligned along the Z-axis. The lipid molecules that overlapped with the protein were removed. The system was then solvated using TIP3P water molecules and ionized with 200 mM NaCl by random substitution of water molecules with Na+ or Cl− ions. The dimension of the system was 90 Å × 90 Å × 80 Å. The whole system consisted of 59 176 atoms, in which 5060 atoms belonged to dvUT, 19 375 atoms to the membrane, 34 698 atoms to TIP3P water molecules, and 43 atoms to the ions. The force-field parameters for lipid, protein, ions, and water were taken from the CHARMM27 force field.22,23 Molecular Dynamic (MD) Simulation. MD simulations were performed for the whole system using the NAMD (NAnoscale Molecular Dynamics) software.24 Simulations were carried out in the NPT ensemble. Pressure was kept at 1 atm using the Langevin piston method.25,26 Temperature was kept at 298 K by coupling to a Langevin thermostat.26 The periodic boundary condition was applied to all the dimensions. Electrostatic interactions were treated using the Particle Mesh Ewald algorithm.27 van der Waals interactions and electrostatic interactions were truncated at 10 Å and smoothened at 8 Å.28 To equilibrate the initial dvUT channel−membrane− solution system, a series of simulations were performed. First, the initial structure was compacted by running energy minimizations for 5000 steps. Then the energy-minimized structure was simulated for 1.35 ns with fixed protein atoms to enable good contact between the protein and lipid molecules. Finally, the refined system was equilibrated for about 16 ns by MD simulations. The time step was set to 1 fs. The root-meansquare deviation (RMSD) of all the heavy atoms in the protein backbone as a function of MD simulation time is shown in Figure S1 of Supporting Information. Steered Molecular Dynamics (SMD) Simulation. Under normal conditions, the time scale for a urea molecule to cross the pore is 10−6 s to 10−4 s, but most MD runs simulate no more than tens of nanoseconds. To investigate urea transport by dvUT, we used the technique of constant-velocity steered molecular dynamics.29 In the SMD simulation, an external force was applied to the carbon atom of the urea molecule to make the urea travel quickly through the pore. The initial
configuration was generated by placing a urea molecule at the entrance of the equilibrium structure of dvUT (i.e., at around Z = −12.0 Å, see the Equilibrium Structure of dvUT and UreaPermeation Pore. To achieve equilibrium, the initial configuration was minimized using 5000 steps followed by 2 ns of MD simulations in the NPT ensemble. A harmonic constraint (with a spring constant of 5.0 kcal/mol·Å2) was attached to the carbon atom of the urea molecule; this carbon atom was pulled along the Z-axis at a constant velocity (10−2 Å/ps) so that the urea molecule was able to pass through the channel in ∼2.5 ns. The corresponding force fluctuation is shown in Figure S2 of Supporting Information. Adaptive Biasing Force (ABF) Calculations. The ABF method has been widely used in conjunction with molecular dynamics simulations to evaluate the free energy of binding for many biochemically interesting systems.30−33 Variation in the free energy of binding, ΔG(Z), along the Z-axis was calculated using the ABF method,34 which relies upon the integration of the average force acting on the urea molecule along the Zdirection. In the NAMD implementation of ABF,35 this force was evaluated within the classical thermodynamic integration formalism.36 The derivative of the free energy, dG(Z)/dZ, was estimated locally throughout the simulation, thereby providing a continuous update of the biasing force. The computation in the present study covered a distance ranging from −13 Å to 11 Å along the Z-axis, for which the origin (Z = 0) was located at the center of mass of the Cα atoms in the protein. In our calculations, the distance between the projected position of the carbon atom in the urea molecule along the Z-axis (i.e., the pore axis) and the center of mass of the Cα atoms in the protein was chosen as the reaction coordinate. To enhance the sampling efficiency, the full reaction coordinate along the Zaxis was divided into 24 nonoverlapping windows and the length of each window was set to 1 Å. For each window, at least 5 ns of MD trajectory were generated to obtain relatively uniform sampling and to reduce statistical errors. This protocol led to a total simulation time of at least 130 ns for the system. Instantaneous values of the force were accrued in bins of 0.1 Å width, which had been used to model transport phenomena.35 The resulting data were then integrated to produce the profile of the variation in the free energy between −13 Å to 11 Å. The starting conformation of the system for each window was taken from a SMD conformation in the corresponding window. 5125
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 2. (A) Variation in the free energy of urea along the axis of dvUT pore. The black arrows indicate the urea-resident sites, which can be considered as the actual urea-binding sites. The red arrows indicate the positions of the ladder oxygen atoms. The gray dashed lines show the swinging range of the hydroxyl arms (Ha1 of T294 and Ha2 of T130). The blue arrows represent the positions of the centers of mass of the hydrophobic residue pairs, including F243/F190, L84/L247, and F27/F80, respectively. (B) The occupation probability of urea along the pore axis (black solid line).
Monte Carlo (MC) Simulation. The MC simulation is a powerful tool for characterizing random walks of particles. To be able to apply MC simulations to urea transport, we need to validate that the urea movement in dvUT can be regarded as diffusion of particles. As shown in our MD simulations, as a urea molecule travels along the pore (Z-axis), the carbonyl group of the urea is parallel to the X-axis most of the time. Therefore, the movement of the urea molecule through the pore can be regarded as a translational motion. Namely, the urea molecule walks along the pore axis (Z-axis) with its carbonyl group parallel to the X-axis as the diffusive motion of a particle in the urea pore. The specific protocol for our MC simulations is described as follows. The urea pore, i.e., the walking space, was divided into a grid of size 0.57 Å × 0.57 Å × 0.57 Å. The orifices of the pore at Z = −11 and 11 Å, respectively, corresponded to level 0 on the extracellular side and level 38 on the intracellular side of the grid, respectively. The walker moved on the grid from site i to its neighboring site j with the jumping probability Pij ∝ exp( −ΔGij/NAkBT )
To quantitatively characterize the pore and to perform permeation simulations, we need to define the locations of the two orifices of the urea pore. For convenience, the mass center of Cα atoms in the protein was set as the origin of the coordinate system (Z = 0). It can be seen from Figure 1A that the cross section of the pore increases dramatically at Z = −11 and 11 Å, where the pore becomes two wide vestibules. Therefore, the cross section at Z = −11 Å (11 Å) was defined as the orifice of the pore to extracellular (intracellular) solution. The more distant from the orifices, the narrower the pore. In the middle region, the pore is so narrow that only one urea molecule can pass through. Consequently, outside the pore the urea molecules are fully hydrated. They are partially dehydrated when they enter the pore and become fully dehydrated when they are in the middle of the pore. Similar dehydration behavior has been reported in detail for the transport of urea molecules in UT-B.39 We further examined our equilibrated structure of dvUT for the detailed composition of the urea pore. Figure 1B and 1C shows the stereoview of the residues lining the pore in two perpendicular directions. As shown in Figure 1B, along one side of the pore, E187 and V188 from the helix Pb and V25 and Q24 from the helix Pa form two linear arrays of three oxygen atoms, which are referred to as oxygen ladders. The ladder atoms O1 and O6 are the side-chain oxygen atoms of E187 and Q24, respectively. The ladder atoms O2, O3, O4, and O5 are the backbone oxygen atoms of E187, V188, V25, and Q24, respectively. These oxygen atoms can form hydrogen bonds with the amide hydrogen atoms (Ahn) of the urea molecule. O1, O2, and O3 form one oxygen ladder, and O4, O5, and O6 form the other ladder. Within each oxygen ladder, the distances between the adjacent oxygen atoms are 3.0−4.0 Å. Therefore, urea can move successively from one oxygen site to the next along the same oxygen ladder by forming and breaking hydrogen bonds. On the opposite side of the oxygen ladders on the pore, L293 and T294 from the helix T5b and T130 and L129 from the helix T5a provide hydroxyl or amide hydrogens that can potentially form hydrogen bonds with the carbonyl oxygen atoms of urea (denoted as Cox). Intriguingly, the distance between the two oxygen ladders is large, i.e., Z(O4) − Z(O3) = 7 Å. It is therefore difficult to move a urea molecule
(1)
where ΔGij = Gj − Gi, and Gi and Gj are the free energies of the walker positioning at the jumping-out-site i and jumping-ontosite j, respectively. NA is the Avogadro constant, and kB is the Boltzmann constant. In the present case, the absolute value of the free energy of the urea molecule in dvUT is unknown; ΔGij is therefore estimated as the difference between the variations in the free energy at site i and site j.
3. RESULTS AND DISCUSSION Equilibrium Structure of dvUT and Urea-Permeation Pore. Constructing an equilibrium conformation of dvUT is a requisite for probing the mechanisms underlying the dvUT functions at the atomic level. Based on the crystal structure of dvUT,12 an equilibrated dvUT channel−membrane−solution system was achieved by running molecular dynamics simulations for 16 ns. The equilibrium structure of dvUT is shown in Figure 1. To search and view the urea-permeation pore, the HOLE software was used.37,38 The identified urea pore is plotted in Figure 1A as a green tube. 5126
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 3. Snapshots describing the transport process of urea in the pore. Zu represents the Z-coordinate of the carbon atom in the urea molecule. The hydrogen bonds are represented by dashed lines. (A) A urea molecule enters the entrance of the pore from the extracellular side. (B) At the external resident site (Sext), the urea forms hydrogen bonds with O3, water, and T294. (C) The hydroxyl arm of T294, Ha1, drags the urea to the middle region of the pore. (D) The urea is located at the unstable equilibrium site near the center of the pore. (E) In the transitional configuration, the hydroxyl arm of T130, Ha2, is dragging the urea to O4. (F) The process of dragging the urea from the middle region to O4 has been completed. (G) At the middle resident site (Smid), urea forms a hydrogen bond with O4 and can be bridged with O5 and/or O6 via water. (H) At the internal resident site (Sint), the urea forms a hydrogen bond with O6 and can interact with the amide hydrogen of L129 via a bridging water molecule.
regions and allow urea molecules to pass smoothly through the urea pore.12,39 Variation of the Free Energy, Occupation Probability Profile, and Urea-Resident Sites. With the equilibrated structure of dvUT obtained from our MD simulations, we were able to characterize the transport properties of dvUT in the present study. We started with investigating the urea transport in dvUT. To find the most stable resident positions of urea in dvUT and then to investigate the dynamics of urea transport, we calculated the variation in the free energy of urea, ΔG(Z), along the pore axis Z using the adaptive biasing force (ABF) approach.34 Figure 2A shows that the curve of the variation in the free energy vs Z exhibits three minima (energy valleys) and multiple maxima (energy barriers). The first three maxima appear in the region where the hydrophobic residues are located. The minima occur mostly in proximity to the positions of the oxygen atoms in the ladders. In the middle region of the pore, there exists a high, wide, and steep barrier. On the basis of the calculated variation in the free energy, we determined the urea-resident sites in dvUT by utilizing the Monte Carlo (MC) method. In the MC simulations, the number of MC-time (i.e., the residence time) at each site of the pore was determined. The details are provided in the Supporting Information SI1. Then, the occupation and cumulative probability of urea vs the Z-coordinate along the whole pore were calculated and plotted in Figure 2B. The sites with the maximum occupation probability values were defined as the urea-resident sites. On the urea-resident sites, urea molecules bound strongly to the pore, and therefore the urea-
from one oxygen ladder to another. Moreover, the plane that contains O1, O2, and O3 forms an angle with the plane that contains O4, O5, and O6, which makes the movement of urea from one oxygen ladder to another more difficult. This difficulty is overcome by using the hydroxyl arms of T294 and T130 as a shuttle. Specifically, the hydroxyl arms (Han) of threonine residues T294 and T130 can form hydrogen bonds with the carbonyl oxygen of urea. The arms are flexible enough to follow the movement of urea as shown in our MD simulations (see Transport Process of Urea in the Pore of dvUT). Thus, the hydroxyl arms of T294 and T130 can facilitate the transition of the urea molecule from one oxygen ladder to another. Different from Figure 1B, Figure 1C shows the equilibrium structure of dvUT from the direction that is perpendicular to the plane containing the oxygen ladders. It can be seen from Figure 1C that there are six hydrophobic side chains. Among them, the phenylalanine residues, F190 and F243 from the periplasmic side and F27 and F80 from the cytoplasmic side, compress the pore into a slot-like shape. The other two hydrophobic residues, L247 and L84 positioning in the middle of the pore, also destabilize polar molecules such as water and urea in this region. In summary, our equilibrium structure of dvUT and MD simulations show that the urea pore can be considered as a composite configuration imposed by the hydrogen bonding pattern with the oxygen ladder in an otherwise hydrophobic pore. These hydrophobic regions coexist with the hydrophilic 5127
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
and the water. The hydroxyl arm, Ha2, of T130 drags the urea molecule to O4 (Figure 3E, a transitional configuration). The urea molecule is then pushed forward by the hydroxyl arm of T294, until the arm Ha2 swings to as far as it can reach (Figure 3F). As the urea molecule moves forward, the hydrogen bond between the urea and T130 breaks. Besides retaining the hydrogen bond between the urea and O4, the urea molecule may connect to O5 and/or O6 by using water molecules as a bridge (Figure 3G, Smid). In this region, the urea molecule appears to have a long dwelling time. The location with the maximum urea-occupation probability is the middle resident site. As the urea molecule crosses the energy barrier that is formed by the hydrophobic residue pair F27/F80 and reaches O6, it may be bridged to the amide hydrogen of L129 via water molecules (Figure 3H, Sint). These hydrogen bonding interactions significantly increase the dwelling time of the urea molecule, and this position becomes the internal resident site. Because O6 is close to the orifice of the pore, once the hydrogen bond between the urea and O6 dissociates, the urea molecule enters the intracellular solution. In summary, our simulations of the transport process of urea along the pore uncover the structural characteristics of the ureabinding sites. Specifically, at the binding sites, multiple hydrogen bonds are formed, particularly, between the urea molecule and dvUT. In addition, all three binding sites of urea are distant from the energy barriers resulting from highly hydrophobic residues (see Figure 2A). Our simulations of the urea transport can also explain why the energy barrier in the middle of the pore is exceedingly high and wide as seen from the calculated profile of the free-energy variation in dvUT (Figure 2A). As shown in this figure, the peak of this energy barrier locates at the position of the hydrophobic residue pair, L247 and L84. The barrier spans a wide range from −3.5 Å to 2.0 Å, which is the range that the hydroxyl arms of T294 and T130 can swing (as indicated by the two vertical dashed lines in Figure 2A). Within this range, the inversion of the hydroxyl arms of T294 and T130 drags the urea molecule across this energy barrier from one oxygen ladder to the other, as shown in Figures 3C−E. Moreover, our simulations also show that in the middle region, the pore is narrow and the urea molecule is completely dehydrated, whereas in the other regions the urea molecule is partially dehydrated. In summary, the middle energy barrier can be attributed to the repulsive effect of the hydrophobic residues L247 and L84, the energetic cost for inverting the hydroxyl arms of T294 and T130, and the desolvation cost. The width of the middle energy barrier coincides with the distance by which the urea molecule is dragged by the hydroxyl arms of T294 and T130. Due to the resulting high energy barrier, it is difficult for the urea molecule to enter the middle region and easy for the urea to leave this region, and therefore it is not surprising that the ureaoccupation probability in the middle region is tiny. Water Permeability of dvUT. Next, we studied the water permeability in dvUT by MD simulations. We started with resolving the puzzle whether dvUT is permeable to water molecules. In the crystal structure of dvUT no water molecule was found.12 However, water molecules were observed in Actinobacillus pleuropneumoniae urea transporter ApUT.14 Recent studies also suggested that urea, H2O, CO2, and NH3 pass through urea transporter B (UT-B).40 In our MD simulations, we observed that a number of water molecules were able to move into or out of the urea pore of dvUT. Specifically, within 20 ns, six water molecules passed through
resident sites are expected to be the actual urea-binding sites. Figure 2B shows that in dvUT the urea-resident sites locate approximately at Z = −5.0 Å, 4.0 Å, and 9.0 Å, which are defined as the external (Sext), middle (Smid), and internal (Sint) resident sites, respectively. The characteristics and underlying mechanisms of these resident sites will be elucidated in the next subsection. Transport Process of Urea in the Pore of dvUT. With the identified urea-resident sites, we characterized the details of the transport process of urea in the pore of dvUT using constant-velocity Steered MD (cv-SMD).29 Specifically, in the presence of the concentration gradient of urea across the membrane, urea molecules diffused along the pore driven by both thermal fluctuations and the concentration gradient. The diffusion of urea through dvUT was simulated using cv-SMD, and the results were analyzed in detail. Our MD simulations showed that the urea molecules that move through the pore are always accompanied by water molecules. Figure 3 shows representative snapshots of the diffusion of urea along the pore to elucidate the mechanisms of urea transport. Specifically, the coordinate of the carbon atom of a urea molecule, CO(NH2)2, in the Z-direction is denoted by Zu. When a urea molecule enters the pore from the extracellular side, an amide hydrogen atom, Ahn, of urea forms a hydrogen bond with the ladder oxygen atom O1 (Figure 3A) or O2. Along the oxygen ladder the urea molecule moves from O1 to O2, crosses the energy barrier that is created by the hydrophobic residue pair, F243 and F190, and then reaches O3. The urea molecule bound to O3 is located at the valley between the energy barrier of F243/F190 and the middle energy barrier. In this region, the amide hydrogen atoms, Ahn, of the urea molecule participate in hydrogen bonding with the ladder oxygen O3 and the hydroxyl group of water, respectively. When the urea molecule reaches the site at Zu = −4.9 Å, the hydroxyl arm of T294, Ha1, forms a hydrogen bond with the carbonyl oxygen, Cox, of the urea molecule (Figure 3B). This threehydrogen-bond configuration stabilizes the urea molecule in the pore, and the position becomes the external resident site (Sext). If the energy provided by thermal fluctuations is large enough to switch the hydroxyl arm of T294, the hydrogen bond between the urea and the water breaks, and the hydroxyl arm Ha1 of T294 drags the urea molecule to the middle region of the pore (Figure 3C). In the middle region, the hydroxyl arm of T130, Ha2, forms a hydrogen bond with the carbonyl oxygen, Cox, of the urea molecule. When the urea molecule reaches the center of the middle region (Figure 3D), the hydrogen bonds, formed by the urea with T294 and T130, are symmetrically arranged on both sides of the urea molecule. In addition, on the opposite side of residues T294 and T130, the hydrophobic side chains of V25 and V188 dislike the urea molecule. The competition between the hydrogen bonding from T294 and T130 and the hydrophobic effect from V25 and V188 destabilizes this equilibrium position, resulting in the peak of the middle energy barrier in Figure 2A. The urea molecule then leaves the peak position of the middle energy barrier and forms a hydrogen bond with O4. Hydrogen bonds are formed between the urea and O4, urea and T294, and urea and T130. This three-hydrogen-bond configuration results in the metastable equilibrium site at Zu = 0.1 Å (the configuration is not shown), which corresponds to the minimum of the middle energy barrier (see Figure 2A). As the urea molecule leaves the middle region, the hydrogen bond between the urea and T294 breaks, and a new hydrogen bond is formed between the urea 5128
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Figure 4. Water molecules entering the urea pore of dvUT. (A) Water molecules are dispersed and partially bonded together via hydrogen bonds. (B) Water molecules form a continuous water chain. The water molecules in the water chain on the two sides of the central water molecule exhibit opposite dipole orientations. The hydrogen bonds are represented by dashed lines. The residues of the protein are depicted in licorice representation, and the water molecules are displayed in CPK representation. (C) The panoramic view of a water chain connecting the extracellular and intracellular solution. Residues of the dvUT monomer are represented as cyan cylinders. The lipid bilayer is removed for a better view. The oxygen and hydrogen atoms of water molecules are represented by red and white van der Waals balls, respectively.
the pore from the extracellular solution to the intracellular solution, and seven water molecules moved in the opposite direction by diffusion. The time taken by these water molecules to cross the pore ranged between 0.74 ns and 6.5 ns, and the average rate was about 0.4 H2O/ns/channel. Then what is the arrangement of these water molecules in the dvUT pore? Our MD simulations showed that in most instances, water molecules were dispersed and partially bound together via hydrogen bonds in the urea pore (Figure 4A), and in some cases, a continuous water chain was formed in the pore (Figure 4B). It can be seen from Figure 4B that in the water chain the oxygen atom of the water molecule in the middle of the pore forms two hydrogen bonds with the hydroxyls of T294 and T130. This central water molecule also forms hydrogen bonds with adjacent water molecules. For the remaining water molecules in the water chain, most of them are able to bind to the oxygen ladders. Interestingly, the water molecules on the two sides of the central water molecule exhibit opposite dipole orientations. For the water molecules arranged in a variety of configurations, we can describe their distribution using the probability density. To calculate the probability density, we took 3706 configurations for statistical analysis. First, we divided the pore into small segments of 0.5 Å in length and counted the number of the water molecules per segment for each configuration. Then we computed the average number of the water molecules in all configurations for each segment and obtained the distribution of the average number of water molecules along the pore (Z axis). Normalization of the distribution of the average number of the water molecule yielded the probability distribution of the water molecules in the dvUT pore. Dividing the probability distribution by the length of segment, we achieved the distribution of probability density of the water molecules in the pore (see Figure 5). To inspect how the structure of the dvUT pore affects the probability density distribution, variation in the cross-sectional area of pore along the Z axis is also plotted in Figure 5. The cross-sectional area was obtained by averaging the measurements from ten configurations of the dvUT pore. It can be seen from Figure 5 that the probability density of the water molecules and the cross-sectional area of the dvUT pore essentially vary in the same tendency, indicating that the distribution of the water molecules is mostly affected by the
Figure 5. Probability density of the water molecules in the dvUT pore and cross-sectional area of the dvUT pore along the Z axis. The black arrows indicate the positions of the ladder oxygen atoms. The blue arrows represent the positions of the centers of mass of the hydrophobic residue pairs. The green dotted lines indicate the positions of Z = −9.0, −1.0, and 7.0 Å.
geometry of the pore. It is well-known that water has a weaker dipole moment (1.85 D) and smaller length (1.5 Å) than urea (4.56 D and 3.9 Å);12 therefore, the influence of the hydrophobic residues and ladder oxygen atoms on water is less than their influence on urea. In contrast, the geometrical characteristic of the pore may have an important effect on the behavior of the water molecules. For example, within the confined range (−9.0 Å < Z < 7.0 Å), the width of the pore ranges between 3.0 and 4.5 Å, which is comparable to the size of a water molecule. As shown in Figure 5, in the confined range, the probability to find water molecules is much lower than that outside this range. In addition to geometry confinement, the chemical composition of the protein also plays a role. Under the influence of the hydrophobic residues, the probability density curve exhibits minimums (around L84/ L247 and F27/F80) or shows an apparent reduction (around F243/F190). On the other hand, the attraction of the ladder oxygen atoms causes a slight increase of the probability density. To provide a microscopic description of the distribution of water molecules in dvUT, we examined the behaviors of the water molecules during the MD simulation. It was found that in the middle of the pore (−2.0 Å < Z < 0 Å), where the pore is 5129
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
midpoint, there are more surrounding water molecules, and the interactions between the adjacent water dipole moments become stronger. Consequently, beyond the range of −4.0 Å < Z < 2.0 Å, the behavior of the orientations of the water dipole moments become more complex. Beyond the confined range, −9.0 Å < Z < 7.0 Å, the number of water molecules increases significantly, which greatly decreases the influence of the dvUT protein on the orientations of the dipole moments of the water molecules. Outside the pore the orientations of the water dipole moments are arbitrary, and ⟨cos θ⟩ approaches 0. Because of computational difficulty, we averaged only with 3706 configurations, and therefore outside the pore the ⟨cos θ⟩ curve did not become horizontal.
very narrow and the influence of the hydrophobic residues L84/L247 is strong, the average number of the water molecules per channel was about 0.4 and the average spacing between adjacent water molecules was larger than 5 Å. In other words, for most of the time there was no water molecule. The probability of finding one water molecule was about one-third, and the probability of finding two water molecules was about 5% in this range. It was a rare event for a water chain to be formed in the dvUT pore. Finally, the polarization behavior of the water molecules in the dvUT pore was also investigated. Specifically, we calculated the projection of the dipole moment of each water molecule in the water chain along the axis of the pore and plotted the curve in Figure 6. It can be seen from this figure that at the midpoint
4. CONCLUSIONS Using MD and MC simulation methods, we studied the structural characteristics and transport behaviors of dvUT. From binding free energy calculations we determined the positions of three urea binding sites, which are located at Z = −5.0 Å, 4.0 Å, and 9.0 Å for the external (Sext), middle (Smid), and internal (Sint) binding sites, respectively. These binding sites are separated by the energy barriers, and the sites participate in multiple hydrogen bonding interactions with the urea molecule. Our MD simulations suggested that dvUT is water permeable, and the average rate is about 0.4 H2O/ns/ channel for the water molecules to cross the dvUT pore. The movement of urea is always accompanied by water. We calculated the probability density of the water molecules in the dvUT pore and the projection of the dipole moments of the water molecules along the axis of the pore. We found that the pseudo-2-fold symmetry of the wall of the pore was the reason for the water molecules on the two sides of the midpoint of the pore to have opposite dipole moment orientations.
Figure 6. Projection of the dipole moments of the water molecules along the axis of the dvUT pore. θ is the angle between the dipole moment and the pore axis. cos θ is the projection of the dipole moment along the pore axis for a water molecule in the water chain when the dipole moment of the water molecule is set as a unit. ⟨cos θ⟩ is the statistical average of cos θ for all water molecules occurring at the same section of the pore when the water molecules are in a dispersed and partially bonded state. The magenta dashed line represents the cross-sectional area of the pore along the pore axis.
■
ASSOCIATED CONTENT
S Supporting Information *
Additional figures on MD and SMD simulations. This material is available free of charge via the Internet at http://pubs.acs.org.
■
of the pore, Z = −1.0 Å, both cos θ and ⟨cos θ⟩ are equal to zero, which indicates that the dipole moment is perpendicular to the pore axis. On the two sides of Z = −1.0 Å, the water molecules have opposite dipole moment orientations. The behavior of the orientation of the water dipole moments is determined by the structure and composition of the dvUT protein. It can be seen from the structure of dvUT (see Figure 1) that the nature and positioning of the side chains forming the wall of the pore are remarkably symmetrical: for every residue along the wall on one side of the midpoint, there exists an identical residue on the other side of the midpoint at the symmetrical location (except the Q24/E187 pair).12 Namely, the wall has a pseudo-2-fold symmetry, and the axis of the symmetry is perpendicular to the axis of the pore, passing through the midpoint of the pore and lying in the plane containing the oxygen ladders. This pseudo-2-fold symmetry makes the water dipole moment perpendicular to the axis of the pore at Z = −1.0 Å, and makes the water molecules have the opposite dipole moment orientations on the two sides of Z = −1.0 Å. By using the hydrogen bond as a bridge a water molecule have an influence not only on the nearest neighboring water molecule but also on the next nearest neighboring water molecules. When the water molecule is moving away from the
AUTHOR INFORMATION
Corresponding Authors
*Tel: (573) 882-6045. E-mail: [email protected]. *Tel: (86) 27-68766665. E-mail: [email protected]. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS We are very grateful to Dr. M. Zhou for helpful discussions. We also thank Benjamin Merideth very much for editing the manuscript. A portion of the computations were performed on the HPC resources at the University of Missouri Bioinformatics Consortium (UMBC). This work was supported by the US National Science Foundation Career Award (Grant DBI0953839) and National Natural Science Foundation of China (grants 11304123 and 11374234).
■
REFERENCES
(1) Nielsen, S.; Terris, J.; Smith, C. P.; Hediger, M. A.; Ecelbarger, C. A.; Knepper, M. A. Cellular and Subcellular Localization of the Vasopressin- Regulated Urea Transporter in Rat Kidney. Proc. Natl. Acad. Sci. U. S. A. 1996, 93, 5495−5500.
5130
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
Article
The Journal of Physical Chemistry B
Molecular Dynamics with NAMD. J. Comput. Chem. 2005, 26, 1781− 1802. (25) Martyna, G. J.; Tobias, D. J.; Klein, M. L. Constant Pressure Molecular Dynamcis Algorithms. J. Chem. Phys. 1994, 101, 4177− 4189. (26) Feller, S. E.; Zhang, Y. H.; Pastor, R. W.; Brooks, B. R. Constant-pressure Molecular Dynamics Simulation- the Langevin Piston Method. J. Chem. Phys. 1995, 103, 4613−4621. (27) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−8593. (28) Miyamoto, S.; Kollman, P. A. SETTLE: An Analytical Version of the SHAKE and RATTLE Algorithm for Rigid Water Molecules. J. Comput. Chem. 1992, 13, 952−962. (29) Izrailev, S.; Stepaniants, S.; Isralewitz, B.; Kosztin, D.; Lu, H.; Molnar, F.; Wriggers, W.; Schulten, K. Steered Molecular Dynamic. In Computational Molecular Dynamics: Challenges, Methods, Ideas; Deuflhard, P.; Hermans, J.; Leimkuhler, B.; Mark, A. E.; Reich, S.; Skeel, R. D, Eds.; Springer-Verlag: Berlin, 1998; Vol. 4 of Lecture Notes in Computational Science and Engineering, pp 39−65. (30) Dehez, F.; Pebay-Peyroula, E.; Chipot, C. Binding of ADP in the Mitochondrial ADP/ATP Carrier Is Driven by an Electrostatic Funnel. J. Am. Chem. Soc. 2008, 130, 12725−12733. (31) Hénin, J.; Tajkhorshid, E.; Schulten, K.; Chipot, C. Diffusion of Glycerol Through Escherichia Coli Aquaglyceroporin GlpF. Biophys. J. 2008, 94, 832−839. (32) Gamini, R.; Sotomayor, M.; Chipot, C.; Schulten, K. Cytoplasmic Domain Filter Function in the Mechanosensitive Channel of Small Conductance. Biophys. J. 2011, 101, 80−89. (33) Modi, N.; Bárcena-Uribarri, I.; Bains, M.; Benz, R.; Hancock, R. E. W.; Kleinekathöfer, U. Role of the Central Arginine R133 Toward the Ion Selectivity of the Phosphate Specific Channel OprP: Effects of Charge and Solvation. Biochemistry 2013, 52, 5522−5532. (34) Darve, E.; Pohorille, A. Calculating Free Energies Using Average Force. J. Chem. Phys. 2001, 115, 9169−9183. (35) Hénin, J.; Chipot, C. Overcoming Free Energy Barriers Using Unconstrained Molecular Dynamics Simulations. J. Chem. Phys. 2004, 121, 2904−2914. (36) den Otter, W. K. Thermodynamic Integration of the Free Energy Along a Reaction Coordinate in Cartesian Coordinates. J. Chem. Phys. 2000, 112, 7283−7292. (37) Smart, O. S.; Neduvelil, J. G.; Wang, X.; Wallace, B. A.; Sansom, M. S. HOLE: a Program for the Analysis of the Pore Dimensions of Ion Channel Structural Models. J. Mol. Graphics 1996, 14, 354−360. (38) Smart, O. S.; Goodfellow, J. M.; Wallace, B. A. The Pore Dimensions of Gramicidin A. Biophys. J. 1993, 65, 2455−2460. (39) Levin, E. J.; Cao, Y.; Enkavi, G.; Quick, M.; Pan, Y.; Tajkhorshid, E.; Zhou, M. Structure and Permeation Mechanism of a Mammalian Urea Transporter. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 11194−11199. (40) Geyer, R. R.; Musa-Aziz, R.; Enkavi, G.; Mahinthichaichan, P.; Tajkhorshid, E.; Boron, W. F. Movement of NH3 Through the Human Urea Transporter B: A New Gas Channel. Am. J. Physiol. Renal Physiol. 2013, 304, F1447−F1457.
(2) Abreu, C.; Sanguinetti, M.; Amillis, S.; Ramon, A. UreA, the Major urea/H+ Symporter in Aspergillus Nidulans. Fungal Genet. Biol. 2010, 47, 1023−1033. (3) Strugatsky, D.; McNulty, R.; Munson, K.; Chen, C.-K.; Soltis, S. M.; Sachs, G.; Luecke, H. Structure of the Proton-gated Urea Channel from the Gastric Pathogen Helicobacter Pylori. Nature 2013, 493, 255−258. (4) Knepper, M. A.; Mindell, J. A. Molecular Coin Slots for Urea. Nature 2009, 462, 733−734. (5) You, G.; Smith, C. P.; Kanal, Y.; Lee, W.-S.; Stelzner, M.; Hedlger, M. A. Cloning and Characterization of the Vasopressinregulated Urea Transproter. Nature 1993, 365, 844−847. (6) Smith, C. P. Mammalian Urea Transporters. Exp. Physiol. 2009, 94, 180−185. (7) Mistry, A. C.; Chen, G.; Kato, A.; Nag, K.; Sands, J. M.; Hirose, S.; Chandra, A. A Novel Type of Urea Transporter, UT-C, Is Highly Expressed in Proximal Tubule of Seawater Eel Kidney. Am. J. Physiol. Renal Physiol. 2005, 288, F455−F465. (8) Weeks, D. L.; Eskandari, S.; Scott, D. R.; Sachs, G. A H+-Gated Urea Channel: The Link Between Helicobacter Pylori Urease and Gastric Colonization. Science 2000, 287, 482−485. (9) Valladares, A.; Montesinos, M. L.; Herrero, A.; Flores, E. An ABC-type, High-affinity Urea Permease Identified in Cyanobacteria. Mol. Microbiol. 2002, 43, 703−715. (10) Kojima, S.; Bohner, a.; Von Wirén, N. Molecular Mechanisms of Urea Transport in Plants. J. Membr. Biol. 2006, 212, 83−91. (11) Godara, G.; Smith, C.; Bosse, J.; Zeidel, M.; Mathai, J. Functional Characterization of Actinobacillus Pleuropneumoniae Urea Transport Protein, ApUT. Am. J. Physiol. 2009, 296, R1268−73. (12) Levin, E. J.; Quick, M.; Zhou, M. Crystal Structure of a Bacterial Homologue of the Kidney Urea Transporter. Nature 2009, 462, 757− 762. (13) Maciver, B.; Smith, C. P.; Hill, W. G.; Zeidel, M. L. Functional Characterization of Mouse Urea Transporters UT-A2 and UT-A3 Expressed in Purified Xenopus Laevis Oocyte Plasma Membranes. Am. J. Physiol. Renal Physiol. 2008, 294, F956−F964. (14) Raunser, S.; Mathai, J. C.; Abeyrathne, P. D.; Rice, A. J.; Zeidel, M. L.; Walz, T. Oligomeric Structure and Functional Characterization of the Urea Transporter from Actinobacillus Pleuropneumoniae. J. Mol. Biol. 2009, 387, 619−627. (15) Zhao, D.; Sonawane, N. D.; Levin, M. H.; Yang, B. Comparative Transport Efficiencies of Urea Analogues Through Urea Transporter UT-B. Biochim. Biophys. Acta 2007, 1768, 1815−1821. (16) Mannuzzu, L. M.; Moronne, M. M.; Macey, R. I. Estimate of the Number of Urea Transport Sites in Erythrocyte Ghosts Using a Hydrophobic Mercurial. J. Membr. Biol. 1993, 133, 85−97. (17) Yang, B.; Verkman, A. S. Analysis of Double Knockout Mice Lacking Aquaporin-1 and Urea Transporter UT-B. Evidence for UT-Bfacilitated Water Transport in Erythrocytes. J. Biol. Chem. 2002, 277, 36782−36786. (18) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual Molecular Dynamics. J. Mol. Graphics 1996, 14, 33−38. (19) Fiser, A.; Sali, A. Modeller: Generation and Refinement of Homology-based Protein Structure Models. Methods Enzymol. 2003, 374, 461−491. (20) Martí-Renom, M. A.; Stuart, A. C.; Fiser, A.; Sánchez, R.; Melo, F.; Sali, A. Comparative Protein Structure Modeling of Genes and Genomes. Annu. Rev. Biophys. Biomol. Struct. 2000, 29, 291−325. (21) Sali, A.; Blundell, T. L. Comparative Protein Modelling by Satisfaction of Spatial Restraints. J. Mol. Biol. 1993, 234, 779−815. (22) MacKerell, A. D.; Bashford, D.; Bellott, M.; Dunbrack, R. L.; Evanseck, J. D.; Field, M. J.; Fischer, S.; Gao, J.; Guo, H.; Ha, S.; et al. All-atom Empirical Potential for Molecular Modeling and Dynamics Studies of Proteins. J. Phys. Chem. B 1998, 102, 3586−3616. (23) Feller, S. E.; MacKerell, A. D. An Improved Empirical Potential Energy Function for Molecular Simulations of Phospholipids. J. Phys. Chem. B 2000, 104, 7510−7515. (24) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.; Skeel, R. D.; Kalé, L.; Schulten, K. Scalable 5131
DOI: 10.1021/jp511300u J. Phys. Chem. B 2015, 119, 5124−5131
