Molecular Mechanism of Acetate Transport Through the Acetate

5 days ago - Acetate is a central metabolite that plays a key role in almost all organisms and acetate channels are often essential for their survival...
0 downloads 0 Views 1MB Size
Subscriber access provided by UNIV OF TEXAS DALLAS

Computational Biochemistry

Molecular Mechanism of Acetate Transport Through the Acetate Channel SatP Meng Wu, Liping Sun, Qingtong Zhou, Yao Peng, Zhijie Liu, and Suwen Zhao J. Chem. Inf. Model., Just Accepted Manuscript • DOI: 10.1021/acs.jcim.8b00975 • Publication Date (Web): 07 Mar 2019 Downloaded from http://pubs.acs.org on March 12, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

Molecular Mechanism of Acetate Transport Through the Acetate Channel SatP Meng Wu,†,‡,§ Liping Sun,† Qingtong Zhou,† Yao Peng,† Zhijie Liu,†,‡ Suwen Zhao*,†,‡

†iHuman

Institute, ShanghaiTech University, 393 Middle Huaxia Road, Shanghai, 201210,

China

‡School

of Life Science and Technology, ShanghaiTech University, 393 Middle Huaxia Road,

Shanghai, 201210, China

§University

of Chinese Academy of Sciences, No. 19A, Yuquan Road, Beijing, 100049, China

*Corresponding author.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

ABSTRACT Acetate is a central metabolite that plays a key role in almost all organisms and acetate channels are often essential for their survival. Recently solved structures of the acetate channel Succinate-Acetate Permease (SatP) provide an atomic view of its closed state. However, the open state of the channel, the key residue conformational changes that trigger the channel to open, and the free energy barrier of acetate transportation remain elusive. To address these questions, we performed microsecond time scale molecular dynamics (MD) simulations and umbrella sampling. Several acetate passing events were observed in the MD trajectories with the application of an external electric field. Further analyses reveal the molecular mechanism of the channel opening, which results from the repacking of key residues, such as Gln50 and Phe17, as well as the subsequent outward movement of all transmembrane helices. Our simulations show that acetate is always surrounded by several water molecules when passing through the channel. Furthermore, a high energy barrier of 15 kcal/mol was observed from the free energy profile generated by umbrella sampling on the closed state of the channel. Our study deepens the understanding of the molecular mechanism of acetate transport through the channel SatP and is expected to facilitate the drug discovery on this target.

ACS Paragon Plus Environment

Page 2 of 32

Page 3 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

INTRODUCTION Acetate is found in almost all organisms as a key metabolite ranging from bacteria, fungi to mammalian cells. Since it can be used as an important carbon source, which is indispensable to the metabolism of carbohydrates and fatty acids,1,2 acetate is critical to the global carbon cycle,3 whereby the complex biomacromolecule is first converted by many bacteria and archaea primarily into acetate and then into carbon dioxide and methane. It is involved in both the synthesis and degradation of macromolecules such as proteins, carbohydrates and lipids.4,5 With a pKa of 4.76, acetate prevails in an anionic form at physiologically neutral pH and requires a transporter or channel for passage across the plasma membrane. There are several kinds of transporters, or channels, which have been found to be involved in the transport of acetate and other monocarboxylates. In mammals, several families of transporters are responsible for the acetate transport, including proton-linked monocarboxylate transporters (MCTs) and sodiumcoupled monocarboxylate transporters (SMCTs).6,7 In addition, the Succinate-Acetate Permease (SatP) protein superfamily, widespread in bacteria, archaea and fungi, is essential for transporting acetate and thus, for the survival of these organisms.8-10 Recently, the crystal structures of SatP from human pathogenic bacteria Citrobacter koseri (ckSatP) and E. coli (ecSatP) were reported by two research groups separately, and both were identified as acetate channels.11,12 The crystal structure of ckSatP (PDB code: 5YS3 and 5YS8) is a homo-hexamer formed by two asymmetric units, each includes three monomers.11 Four bound acetate molecules per monomer were resolved,

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

occupying the binding sites S1, S2, S3, and S4 along the pore, from the cytoplasmic end to the periplasmic end. So far, both the SatP structures solved were in their closed states. While their open states, the key residue conformational changes that trigger the channel to open, and the free energy barrier of acetate transportation, all remains unknown. Recently, molecular dynamics (MD) simulation has been used successfully to study the mechanisms of many channel proteins, including the K+ channel,13-15 Na+ channel,16-18 Cl− channel,19-21 TRPV1,22,23 aquaporin,24,25 as well as channels for organic metabolites such as the Gramicidin channel,26,27 FocA (formate channel),28-30 and HpUreI (urea channel).31 This approach is well established and has provided molecular-level insight into the structure-dynamics-function relationships of the channel proteins. In this study, we investigated the dynamics of the ckSatP monomer embedded in a phospholipid bilayer as well as the characteristics of acetate transport by means of atomistic MD simulations. In total, 6.6 μs unbiased simulations without an electric field and 19.5 μs simulations under an external electric field were performed (see Table 1). Further analysis of the trajectories revealed the repacking of some key residues and their decisive role in opening the channel. In addition, the free energy profile of acetate transport was constructed using umbrella sampling. Energy decomposition for acetate transport through this channel enabled the pinpointing of key residues along the pore and the rationalization of their role in the acetate transport process.

METHODS

ACS Paragon Plus Environment

Page 4 of 32

Page 5 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

Structure Preparation. The crystal structure of ckSatP (PDB code: 5YS3, resolution: 1.8 Å) was used here.11 The structure was a homotrimer with a pairwise RMSD of Cα atoms in the range of 0.22-0.27 Å. Chain A was used since it is the most complete among the three monomers in terms of length and the protomer of the channel is a functional unit.12 We processed the structure of chain A by removing the C terminus tag and adding missing residues in the N terminus using Modeller.32 The protonation state of all ionizable residues was assigned at a neutral pH (7.0) by PROPKA,33 and residue His27 and His151 in the periplasmic side were protonated. System Setup and Molecular Dynamics Simulations. The processed crystal structure of chain A from ckSatP in complex with four acetates was embedded in a bilayer composed of 120 1palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) lipids using the CHARMM-GUI Membrane Builder.34,35 The protein orientation in the membrane is referenced to the HpUreI urea channel (PDB code: 3UX4)36 in the Orientations of Proteins in Membranes (OPM) database.37 The protein-membrane system was solvated in a periodic 0.15 M NaCl TIP3P38 water box with a minimum water height of 15.0 Å on top and bottom of the system. Parameters for acetate were generated using antechamber39 with the application of MOPAC40 semi-empirical quantum mechanical methods. Coordinate files of the simulation system were processed by AmberTools1641 in order to perform MD simulations. All simulations were performed on a GPU cluster using the CUDA version of PMEMD (Particle Mesh Ewald Molecular Dynamics) 42-44 in Amber16.41 The protein was modeled with the ff14SB protein force field,45 acetate with the GAFF2 force field46 and lipids with the AMBER Lipid14

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

force field.47−49 The constructed system was firstly energy minimized for 10,000 steps, of which the first 5,000 steps were performed using the steepest descent method and the remaining 5,000 steps used the conjugate gradient method.50 Then the simulation system was heated from 0 to 100 K using Langevin dynamics51 with a constant box volume. Restraints were applied to protein, acetates, and lipids with a force constant of 10 kcal/mol/Å2. Subsequently, the temperature was increased to 303 K, where the periodic box was coupled accordingly using anisotropic Berendsen control52 in order to maintain the pressure at around 1 atm. Note that, the same restraint of 10 kcal/mol/Å2 was again applied to the protein, acetates and lipids. These restraints were then removed from lipids and the system was equilibrated for another 5 ns at the constant pressure and temperature ensemble (NPT). The Particle mesh Ewald (PME)53 method was used to treat all electrostatic interactions beyond a cutoff of 10 Å. The SHAKE algorithm54 was used for recording the length of bonds involving hydrogen during the simulation with an integration time step of 2 fs. Lastly, the system was equilibrated for 10 ns with restraints applied on the heavy atoms of acetates and the backbone atoms of the protein, and a further 10 ns equilibration with restraints on the heavy atoms of acetates and only Cα atoms of the protein. MD Production Runs With and Without Electric Field. After the step-by-step equilibration, all restraints were removed and multiple unbiased production runs of 500 ns for each were performed, whose trajectories were saved every 10 ps. Total simulation time for these unbiased MD simulations reached 6.6 μs (Table 1). As expected, no permeation events were observed since there was no transmembrane acetate gradient.

ACS Paragon Plus Environment

Page 6 of 32

Page 7 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

In order to mimic the effect of a transmembrane acetate gradient, starting from the equilibrated structure of the unbiased simulations above, an external electric field was applied to the z-axis of the MD system box to accelerate the acetate transport from the bottom of the box to the top, corresponding to acetate transport from the cytoplasmic side of the channel to the periplasmic side in physiological conditions. Hundreds of 50ns MD simulations were performed using three field strengths, 0.05 V/nm, 0.1 V/nm and 0.2 V/nm. All simulations carried out in our study are summarized in Table 1. Table 1. Summary of MD simulations performed in this study Condition

Number of runs

Total run time (μs)

unbiased MD, 150 mM acetate

10

4.3

unbiased MD, 400 mM acetate

6

2.3

0.05 V/nm electric field, 200 mM acetate

144

7.2

0.1 V/nm electric field, 200 mM acetate

130

6.5

0.2 V/nm electric field, 200 mM acetate

115

5.8

Potential of Mean Force. The Potential of Mean Force (PMF) was computed using umbrella sampling in order to gain insight into the energetic determinants of acetate transport through the ckSatP. Using the above unbiased equilibrated structure as the starting structure, Steered Molecular Dynamics (SMD) simulations were initially carried out by pulling an acetate molecule from -25 Å to 20 Å along the z-axis, where the acetate reached in the bulk water in both ends. The pulling rate and the force constant were chosen as 1 Å/ns and 1 kcal/mol/Å2, respectively. The SMD

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

simulation, therefore, lasted 45 ns and was performed in the NPT ensemble with the application of semi-isotropic pressure scaling. A total of 87 windows spaced 0.5 Å apart were used for umbrella sampling, where the acetate was restrained in the z-dimension using only a harmonic force constant of 2.5 kcal/mol/Å2. Each window was simulated for 10 ns, so the cumulative simulation time was 870 ns. The last 5 ns from each window was used for computing the PMF using the weighted histogram analysis method (WHAM) program.55-57 Trajectory analysis was performed using Visual Molecular Dynamics (VMD)58 and CPPTRAJ,59 and visualization was performed using VMD and UCSF Chimera.60

RESULTS No Acetate Permeation in Long Time Scale Unbiased MD Simulations. During the 4.3 μs unbiased MD simulation of ckSatP with acetate concentration at 150 mM, no acetate permeation was captured and elevating the acetate concentration to 400 mM made no difference (Figure S1). Further analysis of the trajectories showed that the conformation of the ckSatP do not change significantly during the simulation as compared to the crystal structure which was in its closed state. Firstly, the root mean square deviations (RMSD) were calculated for the Cα atoms of the channel protein against the crystal structure, showing a comparably small difference of approximately 1.5 Å (Figure S2). In addition, the pore size of the ckSatP channel was analyzed along the trajectories using HOLE.61 The resulting pore radius profile only varied slightly as compared to that of the crystal structure, in particular, the constriction site, which is located at the

ACS Paragon Plus Environment

Page 8 of 32

Page 9 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

middle of the channel and formed by the side chains of three residues Phe17, Tyr72 and Leu131 (Figure 1A). Based on the above evidence, we believe that the channel is still in the closed state after several microseconds of unbiased MD simulation. Interestingly, all four acetate molecules in the binding sites of the crystal structure left the channel from their nearest exits within the first 10 ns of the simulation. Occasionally, several acetate molecules do diffuse into the channel and occupy the binding site S1 and S4. However, binding sites S2 and S3 which near the middle of the channel have never seen any visits of acetates during the whole unbiased MD simulation (Figure S1). Figure 1B shows a three-dimensional representation of the pore of the protein, which reveals the shape and size of the pore. Since the shape of the pore gives an approximate representation of the path that a translocating acetate molecule is likely to follow, the twisted nature of the pore interior suggests the acetate molecule is required to follow a nonlinear path through the pore of the channel.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 1. Pore size and shape of the channel in different conditions. (A) Pore radius profile along the channel axis in the crystal structure (black line), in the unbiased simulations without electric field (orange line), and in the biased simulations with 0.2 V/nm electric field (blue line). The constriction site in the middle of the channel pore is labeled by the vertical dashed red lines. The error bars on orange and blue lines reflect the flexibility of the residues lining the pore during

ACS Paragon Plus Environment

Page 10 of 32

Page 11 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

simulations. The inverted triangle indicates the four acetate binding sites in the crystal structure. (B) The channel structure is represented in cartoon, using the same color code as (A). Four acetate molecules in the crystal structure are represented as sticks. Dotted surfaces describe the lumen of the pore. Red dots correspond to the constricted area where the pore radius is too tight for a water molecule, green and blue dots indicate the pore radii that allow room for a single water molecule and two water molecules, respectively.61 Acetate Transport under External Electric Field. Unbiased microseconds MD simulation of the bilayer-embedded ckSatP monomer in the presence of both 150 mM and 400 mM acetate concentration failed to capture spontaneous acetate permeating events. Limited by the simulation time scale as well as computing resources, other measures were required to accelerate the probability of acetate passage. We applied the external electric field perpendicular to the bilayer membrane and performed a series of short 50 ns MD simulations. Fortunately, among hundreds of 50ns MD trajectories, we captured three acetate-passing events in one trajectory (Figure S3). Given the possibility of changes in the protein structure caused by the applied electric field, we analyzed the RMSD of Cα atoms over the simulation and the root mean square fluctuations (RMSF) of each residue, averaged over the simulation time and compared them with the unbiased simulation (Figure S4). The results show that the structure of the channel protein remains stable throughout the simulation under the external electric field (RMSD < 3 Å, RMSF < 4 Å). We performed several additional MD simulations using the last frame of this acetate-passing trajectory as the starting point, and again several acetate permeation events were captured, which supports that the channel

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

is indeed in the open state. Here, we selected the first three permeation events for analyzing the molecular mechanism of the acetate transport. Figure 2A shows the position of acetate carboxyl carbon atoms on the z-axis as a function of the simulation time. The difference in the time scale of the three passes may indicate that the residues surrounding the passing acetate had differing effects on the acetate as it passed through the channel. Indeed, this effect is related to the orientation of the residues as well as the acetate itself. In addition to the differences in passage time, each of these three transport events had a significant delay in the lower half of the channel, i.e. at the position of about -10 Å. This indicates that the acetate molecule needed to overcome a certain energy barrier, or to take time to adjust its orientation in order to better pass through the channel. Here, to describe the changes in orientation of the passing acetate, we measured the angle formed by the vector along the symmetric axis of the acetate and the z-axis (Figure 3). It can be clearly seen that the orientation of acetate changed significantly during passage through the channel, in particular from S1 to S2. The acetate needed to shift from its original downward orientation of the carboxyl oxygen to an upward orientation. Over the course of the simulations under external electric field as well as the unbiased long-time simulation, water hydrated all parts of the pore including the constriction site (Figure S5), where it appeared that water was passing through the channel along with the acetate (Movie S1, S2 S3). Therefore, we calculated the number of water molecules in the first solvation shell around the acetate (Figure 2B). These results show that acetate is accompanied by water molecules during transport through the channel, suggesting that water has a significant role in the passage of acetate.

ACS Paragon Plus Environment

Page 12 of 32

Page 13 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

Additionally, we computed the number of hydrogen bonds formed between acetate and water based on a definition of a distance between two heavy atoms < 3.5 Å and an angle > 150° (Figure 2C). As can be seen, the acetate forms a higher number of hydrogen bonds with water as compared to the protein during the whole transport process. This suggests that water could provide positive interaction as well as energy to facilitate the transport of acetate.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2. Three permeation events accompanied by water under an external electric field. (A) The loci of the permeating acetates along the pore axis (z-axis) during their transportation in a 50 ns trajectory. (B) The number of water molecules within the first solvation shell of each permeating acetate. The first solvation shell around the acetate based on a distance criterion of 3.4 Å. (C) The number of hydrogen bonds between the acetate and protein/water during the acetate permeation. The hydrogen bond is defined by a distance between the two heavy atoms < 3.5 Å and an angle > 150°.

Figure 3. Re-orientation of acetates during the unbiased simulation (orange circles in A) and during the three permeation events under external electric field (blue circles in B, C, D,

ACS Paragon Plus Environment

Page 14 of 32

Page 15 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

respectively). The re-orientation of the acetate is characterized by angles between the symmetric axis of acetate (from the methyl carbon to the midpoint of two oxygen atoms) and the z-axis. Black triangles represent the orientations of the four bound acetates in each of the three monomers in the ckSatP structure. Potential of Mean Force for the Passage of Acetate. Umbrella sampling was performed to study the energetics of acetate transport through the ckSatP channel. The initial starting structure for each umbrella window was retrieved from the SMD simulations, wherein the studied acetate was pulled through the channel along the pore axis with a constant velocity. The potential of mean force (PMF) was then calculated as a function of the acetate position along the pore axis. We checked the convergence of the PMF with three different sampling time (6 ns, 8 ns and 10 ns) per window (Figure S6). The PMF converged to within 0.6 kcal/mol from 8 ns/window to 10 ns/window, suggesting that the sampling time of 10 ns/window was adequate. The PMF with 10 ns/window is shown in Figure 4A. Since acetate transport in ckSatP proceeds from the cytoplasmic side to the periplasmic side of the channel, the figure shows the cytoplasmic side (starting from z = -25 Å) on the left and the periplasmic side on the right (ending at z = 20 Å). The whole PMF curve is above 0 kcal/mol along the channel and the highest energy point is up to nearly 15 kcal/mol. This indicates that the closed state of the channel represents a relatively high energy barrier, which explains why an acetate transport event were not captured in the long-term unbiased MD process. Of note is that the PMF profile has multiple free energy minima, suggesting the presence of multiple acetate-binding sites along the pore. As expected, the highest energy

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

barrier corresponds to the constriction site of the channel. Here, we compared the major residues that make up the constriction site both at the crystal structure and the representative structure for the highest energy barrier (Figure 4B). In order to accommodate the passage of acetate, residues F17 and Q50 need to change the orientation of their side chains to open the constriction site. It is expected that such reorientation would require a significant amount of free energy due to the limited space near the constriction area with a pore radius of only 1 or 2 Å. An energy decomposition was applied to determine the influence of the residues lining the pore through the channel on the passage of acetate (Figure S7). For each umbrella window simulation, we calculated the linear interaction energy that includes the electrostatic and van der Waals components between acetate and the surrounding residues. As the acetate enters the channel pore from the cytoplasmic side, positively charged residue K60 could have a relatively strong electrostatic interaction with the negatively charged carboxyl end of acetate (Figure S7A). Then the acetate gradually approaches the central core region of the channel, the energy barrier continues to increase until it reaches the constriction site formed by the side chains of F17, Y72, and L131, at which point the energy barrier reaches its maximum. As acetate reaches the constriction site, several important residues inside the channel such as W76, T21, N25 and N28, contribute to its passage in terms of electrostatic energy (Figure S7A). Considering that acetate itself is a small π system and there are many aromatic residues along the pore, stacking interactions are an positive energy contribution that cannot be ignored to help acetate transport. This is partly supported by the van der Waals component of the interaction energy between the acetate and main aromatic

ACS Paragon Plus Environment

Page 16 of 32

Page 17 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

residues lining the pore (Figure S7B). The negative charge of the acetate may even delocalize to the aromatic rings lining up the pore, which may lead to a smaller radius of the acetate and facilitate its permeation.21 Once free of this constriction site, the energy barrier can be gradually reduced, and the passage of the acetate becomes much smoother.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 4. (A) Potential of Mean Force (PMF) profile of the acetate transport through the channel (z-axis). The mean and standard deviations were calculated by three umbrella sampling simulations with different initial configurations. (B) Side view of the superimposed crystal structure (light gray) and the channel structure at the peak of the PMF profile (red). Top view of the orientations of four residues in the constriction site are shown in the enlarged boxes. Acetate Transport Mechanism. Since no spontaneous permeation of acetate was captured during the long-time unbiased MD simulations, we applied an extra external electric field in MD simulations to accelerate acetate transportation on the condition of not destroying the structure of the channel protein. Luckily, three permeation events from hundreds of 50 ns trajectories were observed with the electric field. From the PMF results, we learned that if the acetate is to pass through the channel, the largest energy barrier at the constriction site must be overcome. In the crystal structure, the constriction site was found to be too narrow to allow the passage of acetate or even water through the channel. Therefore, the crystal structure is stabilized in the closed state under the crystallization condition. More importantly, the constrictive residues F17, Y72, and L131 also called “FLY” are almost identical to those of other transporters in the AceTr transporter family, suggesting that members have highly conserved constrictive gates.12 Therefore, we focused on the key residues in the constriction area in order to investigate their effects on the passage of acetate. Among these, F17 has a rotatable aromatic ring which blocks the channel to a certain extent in the crystal structure. Thus, we calculated the dihedral angle of the F17 side chain in the acetate-passing trajectory under external electric field and in the unbiased MD simulations

ACS Paragon Plus Environment

Page 18 of 32

Page 19 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

(Figures 5 and S8). From these data, it can be seen that in the trajectory through which the acetate passes, the χ1 angle of the F17 side chain is more often deflected than in the unbiased simulations. This deflection further drives a corresponding change in the dihedral angle of the adjacent Q50 side chain, which in turn, pushes out the TM1 helix where residue F17 is located. There are two adjacent glycine (G46, G47) near residue Q50, and it is well known that the double glycine motif may be a helical blocker that forms a kink.62 The dihedral angle changes of F17 and Q50 dramatically lead to a local conformational change in the TM2 helix. This local motion of the helix further induces the adjacent TM3 helix to undergo an outward flip change (Figure 6), which directly drives the outward movement of another constriction residue Y72. Ultimately, the outward movement of residues F17 and Q50 (Figure S9) open the constriction site (Figure S10) and thus providing the acetate with an opportunity to pass through the gate.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 5. Dihedral angle of F17 and Q50 side chains in the acetate-passing trajectory under an electric field (blue), and unbiased MD simulation (orange). Superposition of the helical structure with residue F17 and Q50 shown as sticks and colored in blue for the open-state structure under 0.2 V/nm electric field (A), and orange for the close-state structure without electric field (B). Their side chains’ dihedral angles are calculated along the simulation time (C). Changes in the dihedral angles are labeled by numbers in accordance with the numbers in the structure representation at the left panel, with F17-χ1 in the open state labeled as 1 and 2, F17-χ2 in the closed state labeled as 3 and 4 and Q50-χ2 in the closed state as 5 and 6.

ACS Paragon Plus Environment

Page 20 of 32

Page 21 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

Figure 6. The superposition of the crystal structure (gray) and the open-state structure extracted from simulations with external electric field (blue). The open-state structures are high similar for the first two acetate-permeating events and a representative is shown in the left panel, whereas the open-state structure for the third acetate-permeating event is shown in the right panel. Six transmembrane helices are labeled from TM1 to TM6 accordingly. (A) Top view of the two superposition structures and the key residues related to acetate transport, F17 and Q50, are shown. (B) Bottom view of the corresponding superimposed structure in Figure A.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

CONCLUSIONS The characteristics of acetate passage through the ckSatP were studied using atomistic molecular dynamics simulations and as a result a detailed molecular mechanism of acetate transport through the channel has been elucidated. With a time-scale of microseconds unbiased simulations with 150 and 400 mM acetate concentrations, we were able to catch a glimpse of the dynamics of the channel protein but were not able to capture any acetate transport events. In order to increase the likelihood of acetate passage, we applied an external electric field perpendicular to the bilayer membrane. Ensuring that this electric field did not alter the protein structure, we ran hundreds of simulations under different field strengths. A total of three acetate permeation events were captured during these simulations, allowing the investigation of the pathway and its molecular mechanism. Firstly, we found that acetate is transported along with water, which hydrates all parts of the channel pore including the constriction site (Figures 2B and S5). This phenomenon has also been reported in other channel proteins such as the HpUreI urea channel31 and the voltage-gated sodium channels.18 The ckSatP crystal structure revealed one constriction site in the center of the pore formed by residues F17, Y72 and L131, which open transiently in order to allow the passage of acetate. To understand the mechanism behind the opening of this constriction site (i.e., what is the residue level trigger of the transition from closed state to open state), we analyzed side chain rotations of several residues near of constriction site. We found the deflection of the χ1 angle of the F17 side chain drives a corresponding change in the dihedral angle of the Q50 side chain (Figures 5 and

ACS Paragon Plus Environment

Page 22 of 32

Page 23 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

S8). This led to a dramatic local conformational change in the TM2 helix that subsequently pushed the adjacent TM3 helix to undergo an outward flip change (Figure 6). The outward movement of F17 and Q50 (Figure S9) eventually opened the constriction site (Figure S10), thus allowing the passage of acetate through the constriction site. The potential of mean force for acetate passage through the channel was obtained from a series of umbrella sampling simulations. The PMF profile exhibits a significantly high energy barrier of ~15 kcal/mol, which is consistent with the difficulty we experienced in capturing spontaneous passage events in unbiased MD simulations. However, the open state of the channel may have a significantly lower energy barrier considering that the pore radius is much larger around the constriction site in this state. Key residues lining the pore interior were recognized such as K60, W76, F75, T21, N25 and N28 (Figure S7). These residues positively contribute to the passage of acetate through electrostatic interactions, hydrogen bonding or stacking interactions. Once the acetate passes through the constriction site, the remainder of the channel becomes more open and energy-friendly. In summary, we proposed the molecular mechanism of acetate transport through the ckSatP channel and it deepens our understanding of organic anion transport. In addition, key residues for the passage of acetate were identified, which could guide the experimental mechanism study of this channel in the future.

ASSOCIATED CONTENT

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Supporting Information Positions of acetate on z-axis during the unbiased MD simulations (Figure S1); RMSD of C-alpha atoms in the channel protein during unbiased MD simulations (Figure S2); acetate distribution in unbiased MD simulations and simulations under external electric field (Figure S3); RMSD and RMSF of C-alpha atoms in the protein from the unbiased MD simulation and simulation under external electric field (Figure S4); water distribution probability along the channel pore (Figure S5); PMF profiles derived by umbrella sampling with different sampling time per window (Figure S6); nonbonded interactions between the channel protein and acetate (Figure S7); representative conformations of the dihedral angle of F17 and Q50 side chains (Figure S8); distances between L131 to F17, Q50 and Y72 showing the dihedral angle change and the helix movement (Figure S9); cross sectional view of the structure surface showing the opening of the constriction site (Figure S10). Movie S1, S2, S3. This information is available free of charge via the Internet at http://pubs.acs.org

AUTHOR INFORMATION Corresponding Author *Email: [email protected]. Author Contributions

ACS Paragon Plus Environment

Page 24 of 32

Page 25 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

M.W. performed the MD simulations. All authors analyzed the data. M.W., L.S. and S.Z. wrote the manuscript. Z.L. provided valuable suggestions. S.Z. conceived and supervised the project. Notes The authors declare no competing financial interest.

ACKNOWLEDGMENT We thank Dr. Andrej Sali, Dr. Jun Liao, Kang Ding and Yu Guo for their helpful suggestions and discussions. S.Z. was supported by ShanghaiTech University, the NSF of Shanghai grant 16ZR1448500, and National Key Research and Development Program of China grants 2016YFC0905900 and 2018YFA0507000. Y.P. was supported by the NSFC grant 31800633.

REFERENCES (1) Starai, V. J.; Escalante-Semerena, J. C. Acetyl-coenzyme A synthetase (AMP forming). Cell. Mol. Life Sci. 2004, 61, 2020-2030. (2) Wolfe, A. J. The acetate switch. Microbiol. Mol. Biol. Rev. 2005, 69, 12-50. (3) Ferry, J. G. Acetate metabolism in anaerobes from the Domain Archaea. Life.2015, 5, 14541471. (4) Baldwin, J. E.; Krebs, H. The evolution of metabolic cycles. Nature 1981, 291, 381-382. (5) Clark D. P.; Cronan J. E. Two-carbon compounds and fatty acids as carbon sources. EcoSal Plus 2005, 1, 1-34. (6) Halestrap, A. P.; Price, N. T. The proton-linked monocarboxylate transporter (MCT) family: structure, function and regulation. Biochem. J. 1999, 343, 281-299.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(7) Ganapathy, V.; Thangaraju, M.; Gopal, E.; Martin, P. M.; Itagaki, S.; Miyauchi, S.; Prasad, P. D. Sodium-coupled monocarboxylate transporters in normal tissues and in cancer. AAPS. J. 2008, 10, 193-199. (8) Sa-Pessoa, J.; Paiva, S.; Ribas, D.; Silva, I. J.; Viegas, S. C.; Arraiano, C. M.; Casal, M. SATP (YaaH), a succinate-acetate transporter protein in Escherichia coli. Biochem. J. 2013, 454, 585595. (9) Rabitsch, K. P.; Toth, A.; Galova, M.; Schleiffer, A.; Schaffner, G.; Aigner, E.; Rupp, C.; Penkner, A. M.; Moreno-Borchart, A. C.; Primig, M.; Esposito, R. E.; Klein, F.; Knop, M.; Nasmyth, K. A screen for genes required for meiosis and spore formation based on whole-genome expression. Curr. Biol. 2001, 11, 1001-1009. (10) Casal, M.; Queiros, O.; Talaia, G.; Ribas, D.; Paiva, S. Carboxylic acids plasma membrane transporters in Saccharomyces cerevisiae. Adv. Exp. Med. Biol. 2016, 892, 229-251. (11) Qiu, B.; Xia, B.; Zhou, Q.; Lu, Y.; He, M.; Hasegawa, K.; Ma, Z.; Zhang, F.; Gu, L.; Mao, Q.; Wang, F.; Zhao, S.; Gao, Z.; Liao, J. Succinate-acetate permease from Citrobacter koseri is an anion channel that unidirectionally translocates acetate. Cell Res. 2018, 28, 644-654. (12) Sun, P.; Li, J.; Zhang, X.; Guan, Z.; Xiao, Q.; Zhao, C.; Song, M.; Zhou, Y.; Mou, L.; Ke, M.; Guo, L.; Geng, J.; Deng, D. Crystal structure of the bacterial acetate transporter SatP reveals that it forms a hexameric channel. J. Biol. Chem. 2018, 118, 1-17. (13) Ostmeyer, J.; Chakrapani, S.; Pan, A.C.; Perozo, E.; Roux, B. Recovery from slow inactivation in K+ channels is controlled by water molecules. Nature 2013, 501, 121-124. (14) Fowler, P. W.; Abad, E.; Beckstein, O.; Sansom, M. S. Energetics of multi-ion conduction pathways in potassium ion channels. J. Chem. Theory Comput. 2013, 9, 5176-5189. (15) Fowler, P. W.; Beckstein, O.; Abad, E.; Sansom, M. S. Detailed examination of a single conduction event in a potassium channel. J. Phys. Chem. Lett. 2013, 4, 3104-3109. (16) Corry, B.; Thomas, M. Mechanism of ion permeation and selectivity in a voltage-gated sodium channel. J. Am. Chem. Soc. 2012, 134, 1840-1846.

ACS Paragon Plus Environment

Page 26 of 32

Page 27 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

(17) Oakes V.; Furini S.; Domene C. Voltage-gated sodium channels: mechanistic insights from atomistic molecular dynamics simulations. Curr Top Membr. 2016, 78, 183-214. (18) Ulmschneider, M. B.; Bagneris, C.; McCusker, E. C.; Decaen, P. G.; Delling, M.; Clapham, D. E.; Ulmschneider, J. P.; Wallace, B. A. Molecular dynamics of ion transport through the open conformation of a bacterial voltage-gated sodium channel. Proc. Nat. Acad. Sci. USA. 2013, 110, 6364-6369. (19) Cheng, M. H.; Coalson, R. D. Molecular dynamics investigation of Cl− and water transport through a eukaryotic CLC transporter. Biophys. J. 2012, 102, 1363-1371. (20) McKiernan, K. A.; Koster, A. K.; Maduke, M.; Pande, V. S. Dynamical model of the CLC2 ion channel exhibits a two-step gating mechanism. bioRxiv 2018, 228163, 1-25. (21) Wang, C. H.; Duster, A. W.; Aydintug, B. O.; Zarecki, M. G.; Lin, H. Chloride ion transport by the E. coli CLC Cl−/H+ antiporter: a combined quantum-mechanical and molecular-mechanical study. Front. Chem. 2018, 6, 1-17. (22) Chugunov, A. O.; Volynsky, P. E.; Krylov, N. A.; Nolde, D. E.; Efremov, R. G. Temperature-sensitive gating of TRPV1 channel as probed by atomistic simulations of its transand juxtamembrane domains. Sci. Rep. 2016, 6, 1-16. (23) Hanson S. M.; Newstead S.; Swartz K. J.; Sansom M. S. P. Capsaicin interaction with TRPV1 channels in a lipid bilayer: molecular dynamics simulation. Biophys. J. 2015, 108, 14251434. (24) Tajkhorshid, E.; Nollert, P.; Jensen, M. O.; Miercke, L. J.; O'Connell, J.; Stroud, R. M.; Schulten, K. Control of the selectivity of the aquaporin water channel family by global orientational tuning. Science 2002, 296, 525-530. (25) Zhu, F.; Tajkhorshid, E.; Schulten, K. Theory and simulation of water permeation in aquaporin-1. Biophys. J. 2004, 86, 50-57. (26) B. Roux. Computational studies of the gramicidin channel. Acc. Chem. Res. 2002, 35, 366375.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(27) Mustafa, M.; Henderson, D. J.; Busath, D. D. Computational studies of gramicidin permeation: an entry way sulfonate enhances cation occupancy at entry sites. Biochim. Biophys. Acta 2009, 1788, 1404-1412. (28) Feng, Z.; Hou, T.; Li, Y.. Concerted movement in pH-dependent gating of FocA from molecular dynamics simulations. J. Chem. Inf. Model. 2012, 52, 2119-2131. (29) Lv, X.; Liu, H.; Ke, M.; Gong, H. Exploring the pH-dependent substrate transport mechanism of FocA using molecular dynamics simulation. Biophys. J. 2013, 105, 2714-2723. (30) Padhi, S.; Reddy, L. K.; Priyakumar, U. D. pH-mediated gating and formate transport mechanism in the Escherichia coli formate channel. Mol. Simul. 2017, 43, 1300-1306. (31) McNulty, R.; Ulmschneider, J. P.; Luecke, H.; Ulmschneider, M. B. Mechanisms of molecular transport through the urea channel of Helicobacter pylori. Nat. Commun. 2013, 4, 1-10. (32) Sali, A.; Blundell, T. L. Comparative protein modelling by satisfaction of spatial restraints. J. Mol. Biol. 1993, 234, 779-815. (33) Olsson, M. H.; Sondergaard, C. R.; Rostkowski, M.; Jensen, J. H. PROPKA3: consistent treatment of internal and surface residues in empirical pKa predictions. J. Chem. Theory Comput. 2011, 7, 525-537. (34) Jo, S.; Kim, T.; Iyer, V. G.; Im, W. CHARMM-GUI: a web-based graphical user interface for CHARMM. J. Comput. Chem. 2008, 29, 1859-1865. (35) Wu, E. L.; Cheng, X.; Jo, S.; Rui, H.; Song, K. C.; Davila- ́Contreras, E. M.; Qi, Y.; Lee, J.; Monje-Galvan, V.; Venable, R. M.;Klauda, J. B.; Im, W. CHARMM-GUI Membrane Builder toward realistic biological membrane simulations. J. Comput. Chem. 2014, 35,1997-2004. (36) Strugatsky, D.; McNulty, R.; Munson, K.; Chen, C. K.; Soltis, S. M.; Sachs, G., and Luecke, H. Structure of the proton-gated urea channel from the gastric pathogen Helicobacter pylori. Nature 2013, 493, 255-258.

ACS Paragon Plus Environment

Page 28 of 32

Page 29 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

(37) Lomize, M. A.; Pogozheva, I. D.; Joo, H.; Mosberg, H. I.; Lomize, A. L. OPM database and PPM web server: resources for positioning of proteins in membranes. Nucleic Acids Res. 2012, 40, 370-376. (38) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 1983, 79, 926-935. (39) Wang, J.; Wang, W.; Kollman, P. A.; Case, D. A. Automatic atom type and bond type perception in molecular mechanical calculations. J. Mol. Graphics Modell. 2006, 25, 247-260. (40) Walker, R. C.; Crowley, M. F.; Case, D. A. The implementation of a fast and accurate QM/MM potential method in Amber. J. Comput. Chem. 2008, 29, 1019-1031. (41) D.A. Case, R.M. Betz, D.S. Cerutti, T.E. Cheatham, III, T.A. Darden, R.E. Duke, T.J. Giese, H. Gohlke, A.W. Goetz, N. Homeyer, S. Izadi, P. Janowski, J. Kaus, A. Kovalenko, T.S. Lee, S. LeGrand, P. Li, C. Lin, T. Luchko, R. Luo, B. Madej, D. Mermelstein, K.M. Merz, G. Monard, H. Nguyen, H.T. Nguyen, I. Omelyan, A. Onufriev, D.R. Roe, A. Roitberg, C. Sagui, C.L. Simmerling, W.M. Botello-Smith, J. Swails, R.C. Walker, J. Wang, R.M. Wolf, X. Wu, L. Xiao and P.A. Kollman (2016), AMBER 2016, University of California, San Francisco. (42) Gotz, A. W.; Williamson, M. J.; Xu, D.; Poole, D.; Le Grand, S.; Walker, R. C. Routine microsecond molecular dynamics simulations with AMBER on GPUs. 1. Generalized Born. J. Chem. Theory Comput. 2012, 8, 1542-1555. (43) Salomon-Ferrer, R.; Gotz, A. W.; Poole, D.; Le Grand, S.; Walker, R. C. Routine microsecond molecular dynamics simulations with AMBER on GPUs. 2. explicit solvent Particle Mesh Ewald. J. Chem. Theory Comput. 2013, 9, 3878-3888. (44) Le Grand, S.; Gotz, A. W.; Walker, R. C. SPFP: Speed without compromise-A mixed precision model for GPU accelerated molecular dynamics simulations. Comput. Phys. Commun. 2013, 184, 374-380.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(45) Maier, J. A.; Martinez, C.; Kasavajhala, K.; Wickstrom, L.; Hauser, K. E.; Simmerling, C. ff14SB: Improving the accuracy of protein side chain and backbone parameters from ff99SB. J. Chem. Theory Comput. 2015, 11, 3696-3713. (46) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and testing of a general amber force field. J. Comput. Chem. 2004, 25, 1157-1174. (47) Dickson, C. J.; Rosso, L.; Betz, R. M.; Walker, R. C.; Gould, I. R. GAFFlipid: a General Amber Force Field for the accurate molecular dynamics simulation of phospholipid. Soft Matter. 2012, 8, 9617-9627. (48) Skjevik, A. A.; Madej, B. D.; Walker, R. C.; Teigen, K. LIPID11: a modular framework for lipid simulations using amber. J. Phys. Chem. B 2012, 116, 11124-11136. (49) Dickson, C. J.; Madej, B. D.; Skjevik, A. A.; Betz, R. M.; Teigen, K.; Gould, I. R.; Walker, R. C. Lipid14: the Amber lipid force field. J. Chem. Theory Comput. 2014, 10, 865-879. (50) Press, W. H.; Teukolsky, S. A.; Vetterling, W. T.; Flannery, B. P. Conjugate gradient methods in multi-dimensions in numerical recipes: The Art of Scientific Computing, 3rd ed.; Cambridge University Press: New York, 2007; Chapter 10.8, pp 515-518. (51) Pastor, R.; Brooks, B.; Szabo, A. An analysis of the accuracy of Langevin and molecular dynamics algorithms. Mol. Phys. 1988, 65, 1409-1419. (52) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. Molecular dynamics with coupling to an external bath. J. Chem. Phys. 1984, 81, 3684-3690. (53) Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An Nlog(N) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089-10092. (54) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys. 1977, 23, 327.

ACS Paragon Plus Environment

Page 30 of 32

Page 31 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

(55) Kumar, S.; Rosenberg, J. M.; Bouzida, D.; Swendsen, R. H.; Kollman, P. A. The weighted histogram analysis method for free-energy calculations on biomolecules. I. The Method. J. Comput. Chem. 1992, 13, 1011-1021. (56) Roux, B. The calculation of the potential of mean force using computer simulations. Comput. Phys. Commun. 1995, 91, 275-282. (57) Grossfield, A. WHAM: The Weighted Histogram Analysis Method, version 2.0.9, http://membrane.urmc.rochester.edu/content/wham. (58) Humphrey, W., Dalke, A., and Schulten, K. VMD: Visual molecular dynamics. J. Mol. Graphics 1996, 14, 33-38. (59) Roe, D. R.; Cheatham, T. E., 3rd. PTRAJ and CPPTRAJ: Software for processing and analysis of molecular dynamics trajectory data. J. Chem. Theory Comput. 2013, 9, 3084-3095. (60) Pettersen, E. F.; Goddard, T. D.; Huang, C. C.; Couch, G. S.; Greenblatt, D. M.; Meng, E. C.; Ferrin, T. E. UCSF Chimera--a visualization system for exploratory research and analysis. J. Comput. Chem. 2004, 25, 1605-1612. (61) Smart, O. S.; Neduvelil, J. G.; Wang, X.; Wallace, B. A.; Sansom, M. S. P. HOLE: A program for the analysis of the pore dimensions of ion channel structural models. J. Mol. Graphics 1996, 14, 354-360. (62 Aurora, R.; Srinivasan, R.; Rose, G. D. Rules for alpha helix termination by glycine. Science 1994, 264, 1126-1130.

ACS Paragon Plus Environment

Journal of Chemical Information and Modeling 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

For Table of Contents Only

ACS Paragon Plus Environment

Page 32 of 32