On Conduction in a Bacterial Sodium Channel

Voltage-gated Na+-channels are transmembrane proteins that are responsible for the fast depolarizing phase of the action potential in nerve and muscular cells. Selective permeability of Na+ over Ca2+ or K+ ions is essential for the biological function of Na+-channels. After the emergence of the first high-resolution structure of a Na+-channel, an anionic coordination site was proposed to confer Na+ selectivity through partial dehydration of Na+ via its direct interaction with conserved glutamate side chains. By combining molecular dynamics simulations and free-energy calculations, a low-energy permeation pathway for Na+ ion translocation through the selectivity filter of the recently determined crystal structure of a prokaryotic sodium channel from Arcobacter butzleri is characterised. The picture that emerges is that of a pore preferentially occupied by two ions, which can switch between different configurations by crossing low free-energy barriers. In contrast to K+-channels, the movements of the ions appear to be weakly coupled in Na+-channels. When the free-energy maps for Na+ and K+ ions are compared, a selective site is characterised in the narrowest region of the filter, where a hydrated Na+ ion, and not a hydrated K+ ion, is energetically stable.


Introduction
Sodium channels allow the passive diffusion of Na + ions down their electrochemical gradient. They were first discovered, together with potassium selective channels, in nerve fibres where they mediate the fast depolarizing phase of action potentials [1]. Na + -channels are also involved in the initiation of action potentials in cardiac myocytes and in general, in the propagation of electrical impulses in cardiac, muscle and nerve cells. The cell membrane is exposed to a high-Na + /low-K + concentration on the extracellular side, and to a low-Na + /high-K + concentration on the intracellular side. It is the selectivity of Na + and K + channels, in the presence of these concentration gradients what makes possible the control of the membrane potential. Several atomic structures of K + -selective channels have been solved since 1998 [2,3,4,5,6]. Recently, the atomic structure of a Na + -selective channel, NavAb from the bacterium Arcobacter butzleri, has also been determined [7], and provides the opportunity to investigate how selectivity and conduction are realized at the atomic level in Na + -channels.
Na + and K + -channels share the same general architecture, with four transmembrane domains contributing to a central pore. This central pore, where ion permeation occurs, is delimited by two transmembrane helices, S5 and S6, linked by an intervening loop. The Na + -channel NavAb was crystallized in the closed state, where the four S6 helices are arranged in a conical shape, defining a bundle-crossing at the intracellular side and a water-filled cavity above it. In the open state, the S6 helices are thought to bend, and the water-filled cavity becomes a continuous with the intracellular solution, as observed in K + -channel [8,9]. The region responsible for selective permeation is located in the loop between S5 and S6. In K + -channels, each chain contributes to the selectivity filter with the signature sequence T/S-x-G-Y/F-G. The oxygen atoms of the residues from this signature sequence define a series of four binding sites where dehydrated K + ions bind [2]. The selectivity filter of the bacterial Na + -channel NavAb is much wider than that of known K + -channels ( Figure 1). Although Na + ions were not observed in the crystallographic structure, the presence of three Na + binding sites was hypothesized [7]. Residues Thr175 and Leu176 define two rings of carbonyl oxygen atoms at the intracellular entrance of the filter. Water molecules from the hydration shell of a Na + ion may interact through hydrogen bonds with these two layers of carbonyl oxygen atoms, thus defining two binding sites for hydrated Na + ions. Following the convention introduced by Payandeh et al [7], these binding sites will be referred as S IN and S CEN . A third Na + binding site, S HFS , composed by side chains of four glutamate residues, Glu177, from each of the four protein chains was also proposed [7]. The side chain of Glu177 points to the lumen of the filter, and they could contribute to the hydration shell of cations.
Molecular Dynamics (MD) simulations were crucial for understanding the mechanisms of conduction and selectivity in K + -channels. The knock-on mechanism, first hypothesized in the fifties by Hodgkin and Keynes [10], was described at atomic detail by computational studies [11,12], and free-energy calculations revealed energy barriers of the order of 2-4 kcal/mol for the concerted motion of three K + ions through the filter of K +channels [13,14,15]. These free-energy barriers along a multi-ion permeation pathway were shown to increase if one of the permeating ions was Na + instead of K + , offering some insight about selectivity in K + -channels [16,17,18]. The recent crystal structure of a prokaryotic voltage gated sodium channel NavAb from Arcobacter butzleri in combination with all-atom MD simulations can also provide important information concerning conduction and selectivity in Na + -channels. Here, simulation strategies similar to the ones used for K + -channels were adopted to analyze conduction and selectivity in the NavAb Na + -channel, with the purpose to get insights into: (i) the number of Na + binding sites, (ii) the nature of the conduction mechanism, and (iii) selectivity.

Results
Root mean square deviation of the protein backbone atoms below 1.4 Å , and 1.0 Å in the case of the backbone atoms of the filter, residue Thr175 to Ser178, were observed in a 40-ns MD trajectory, which suggests that the system is stable. The selectivity filter was initially occupied by a single Na + ion in S CEN , which remained at this position for the entire trajectory. A second Na + ion approached the filter from the extracellular side during the course of the simulation, and it stayed close to S HFS for more than 10 ns. The intracellular cavity below the selectivity filter was initially occupied by 15 water molecules. The number of water molecules inside the cavity reached equilibrium value of ,45 in the course of the simulation, due to the incorporation of water molecules from the extracellular solution across the selectivity filter. This contrasts with the situation in K + -channels, where the movement of water molecules through the selectivity filter was highly unlikely [19]. Similar observations where described by Klein et al [20].
MD simulations in the nanosecond time scale cannot reveal the mechanisms of conduction and selectivity in Na + -channels alone. In order to get further insights into the permeation process, the permeation free-energy profiles for Na + and K + ions in the Na +channel were calculated by the umbrella sampling technique (see the Methods section for the computational details). The freeenergy profile of a single Na + ion moving from the intracellular cavity to the extracellular solution displays two minima with similar energies (Figure 2). In the innermost minimum, the ion lies in the plane defined by the carbonyl oxygen atoms of residues Leu176, binding site S CEN . Figure S1 shows the positions of the permeating ion in a plane perpendicular to the permeation axis. At S CEN , a Na + ion preserves its entire hydration shell, i.e. it is coordinated only by oxygen atoms from the surrounding water molecules. The nature of the atoms of the coordination shell of a Na + ion at different positions along the permeation pore is shown in Figure S2. In the region above S CEN , the permeating ion may move ,2 Å away from the pore axis, where it interacts with a carbonyl oxygen atom of Leu176 from one of the chains. This position does not correspond to a minimum in the free-energy profile. Proceeding toward the extracellular side, a local freeenergy minimum is encountered with the ion close to the pore axis, in a region 1-2 Å below the oxygen atoms of the side chain of Glu177. At this position, the Na + ion is fully solvated by water molecules. This local minimum is followed by the second global free-energy minimum, where the ion occupies S HFS . The ion moves away from the centre of the pore axis ( Figure S1) and two of the water molecules from the hydration shell are lost, and substituted by oxygen atoms from the side chain of Glu177 and Ser178 ( Figure S2).
The free-energy profile of a single K + ion permeating the selectivity filter of the Na + -channel is markedly different from the free-energy profile of single Na + ion ( Figure 2). In the most stable minimum, a K + ion occupies a position similar to the one occupied by Na + in its innermost minimum. As observed with Na + , a K + ion at this position preserves intact its hydration shell. The second minimum is analogous to the outermost minimum encountered by Na + , with oxygen atoms from the side chains of Glu177 and Ser178 contributing to the coordination shell of the ion. In contrast to Na + , no local minimum is present in the region between the carbonyl oxygen atoms of Leu176 and the side chain oxygen atoms of Glu177 in the permeation pathway of K + . Overall, the composition of the coordination shell and the position of the ion in a plane perpendicular to the pore axis are analogous for K + and Na + ions ( Figure S1 and S2). Likewise, a free-energy barrier higher than 8 kcal/mol hampers the movement of both Na + and K + ions from S HFS to the extracellular solution. These high free-energy barriers associated with the exit of Na + and K + from the selectivity filter may originate from the binding affinity for cations of this region of the pore, where the carboxyl groups of the side chain of Glu177 line the pore on the extracellular side. Therefore, under these circumstances, to enhance conduction, it seems possible that the selectivity filter is occupied on average by one or more positive ions, and that conduction occurs when an incoming ion displaces the ion(s) already inside the pore, like described in K + -channels [14]. To test this hypothesis, free-energy maps for conduction events with two Na + ions were calculated ( Figure 3). It was found that in the most stable configuration, both S CEN and S HFS are occupied by Na + ions. Starting from this position and moving the innermost ion toward the intracellular side, a local minimum is observed with the bottom ion inside the cavity and the top ion in S HFS . The free energy of this local minimum is 4.360.4 kcal/mol higher than the free-energy of the minimum with ions in S CEN and S HFS . The differences in free energy and the energy barriers were evaluated along the minimum

Author Summary
Ion channels are integral membrane proteins that control the passive diffusion of ions down their electrochemical gradient. According to the most permeating ion species, ion channels are classified in three categories: K + -channels, Na + -channels, and Ca 2+ -channels. The atomic structure of a K + -channel was the first to be solved experimentally more than 10 years ago. This structure inspired numerous computational studies, which revealed the mechanisms of conduction and selectivity in K + -channels. Recently, the first atomic structure of a Na + selective channel has been solved. Here, molecular dynamics simulations and freeenergy calculations are described and a possible mechanism for Na + conduction is identified. In contrast to what it is observed in K + -channels, ion movements through Na +channels appeared highly uncorrelated. energy path that connects the local energy minima, and the errors were estimated dividing the trajectories into three data sets (see the Methods section for further details). A second free-energy minimum is observed with both ions close to Glu177 and Ser178. The free energy of this configuration is 1.260.6 kcal/ mol higher than the free energy of the minimum with ions in S CEN and S HFS . The free-energy barriers between the two basins are 3.560.5 kcal/mol and 2.460.3 kcal/mol respectively in the outward and inward directions. Regardless of the position of the bottom ion (S CEN or S HFS ), the movement between the extracellular solution and the pore of the outermost ion does not encounter any free-energy barriers higher than 3 kcal/mol.
In order to understand if the low free-energy barriers for ion translocation across the pore were exclusive of Na + ions, freeenergy maps for conduction events with two K + ions were also calculated ( Figure 4). Like in the case of Na + ions, the lowest freeenergy configuration has two ions occupying binding sites S CEN and S HFS , but in contrast to a pore occupied by Na + ions, this freeenergy minimum is much broader. A local minimum with both ions close to Glu117 and Ser178 can be also described. The free energy of this configuration is 2.760.5 kcal/mol higher than the free energy of the lowest minimum, and the barriers between the two basins are 5.060.4 kcal/mol and 2.360.5 kcal/mol respectively in the outward and inward directions. The inward movement of the ion at S CEN toward the intracellular cavity while the second ion is at S HFS causes an energy increase of 1.860.2 kcal/mol. The change in free energy is lower than 1 kcal/mol if the ion in S HFS leaves the filter toward the extracellular solution and the barriers for the entry/exit of K + in S IN and S HFS are below 2 kcal/mol. Under physiological conditions, the cellular membrane is exposed to high-Na + /low-K + concentration on the extracellular side, and to a low-Na + /high-K + concentration on the intracellular side. As a consequence, the selectivity filter of Na + -channels is more likely to be approached by Na + ions on the extracellular side and by K + ions on the intracellular side. To characterize this situation, the free-energy map for a mixture of K + and Na + ions in the selectivity filter was calculated, with the K + ion in the innermost position and the Na + ion in the outermost position ( Figure 5). The free-energy profile characterising the outward movement of K + through the selectivity filter is similar to that when both ions are K + . In the lowest free-energy configuration, the Na + ion is in site S HFS , while the K + ion can move in a wide region around S CEN , experiencing low-energy barriers (,2 kcal/ mol). A second local minimum exists with both ions close to Glu177 and Ser178 (2.660.8 kcal/mol higher in energy). The barriers between the two minima are similar to those observed in the case of two K + ions (4.661.1 in the outward direction, 2.060.5 in the inward direction).

Discussion
Since the number and position of Na + binding sites in the NavAb channel were only speculated using the crystallographic structure, insight from atomistic simulations is essential. Two Na + binding sites have been described in the selectivity filter in Figure 2. Potential of Mean force for conduction events with one Na + ion and one K + ion. Values on the x-axis correspond to the distance between the permeating ion and the centre of the carbonyl oxygen atoms of residues Thr175 along the pore axis. The free-energy for a permeating Na + ion is shown in a blue continuous line and K + in a pink dashed line. Snapshots from the umbrella sampling trajectories are shown for some significant configurations. Residues 175 to 178 of two opposite subunits are depicted in licorice representation, together with the permeating ion in VDW spheres, yellow and green for Na + and K + . Water molecules closer than 5 Å (licorice representation) and 15 Å (grey lines) from the ion are also shown. doi:10.1371/journal.pcbi.1002476.g002 simulations with either one or two permeating ions. These binding positions correspond to sites S CEN and S HFS , which were hypothesized experimentally. Computations confirm that a Na + ion remains close to the centre of the pore axis and it is fully hydrated in the innermost binding site. In contrast, it occupies a position slightly off the centre of the pore axis, while it interacts directly with the side chain of Glu177 in the outermost site. The side chain oxygen atom of Ser178 also contributes to the coordination shell of a Na + ion in S HFS , which was not an obvious experimental observation. The experimentally speculated binding site S IN does not appear as a well defined minimum in the freeenergy maps from atomistic simulations. A configuration with two ions in S CEN and S HFS is found to be the most stable during conduction events with two permeating ions. Surprisingly, the energetic cost of having two Na + ions in close proximity of binding site S HFS is extremely low (,1 kcal/mol). Therefore, the picture that emerges from the free-energy profiles is that of a pore preferentially occupied by two ions, which can switch between different configurations by crossing low energy-barriers. In contrast to K + -channels, the movements of the ions appear to be weakly coupled.
The nature of the conduction mechanism for K + ions in this Na + -channel was also tested, in order to understand selectivity. In the case of Na + , the strongest impediment to conduction is the raise in free energy associated with the inward movement of an ion from S CEN to the intracellular cavity. Similarly, the free energy also increases when a K + ion moves from S CEN to the intracellular cavity, but to a lower extent (,2 kcal/mol for K + ; ,4 kcal/mol for Na + ). The free energy of an ion in the cavity is likely to be affected by the conformation of the intracellular gate in NavAb, which is closed. And this effect may be different for K + and Na + , as already observed in K + -channels [21].
The biological function of Na + -channels is to allow Na + ions inside the cell at nearly the rate of free diffusion, while at the same time blocking the flux of K + ions in the outward direction. During permeation in the selectivity filter, the highest free-energy barrier that a Na + ion experiences in the inward direction is ,2 kcal/mol (independently of the nature of the innermost ion, Na + or K + ), while the highest energy barrier for the outward movement of K + ions is ,5 kcal/mol (independently of the nature of the outermost ion, Na + or K + ). Considering that discrimination of ions is performed with an efficiency of one K + ion every 10-100 Na + ions [7,22], the differences in free-energy barriers of 1-2 kcal/mol described before are perfectly in line with the experimental observations. The highest barrier for the passage of ions through the filter is localized between S CEN and S HFS , both for K + and Na + . In the case of Na + , a local free-energy minimum exists in this region, with the ion fully hydrated and aligned at the centre of the channel axis. In contrast, the same configuration is not a freeenergy minimum for K + . This difference is already evident in the free-energy profiles with a single permeating ion. The instability of a hydrated K + ion in the middle of the selectivity filter is a possible cause of the preference of the channel for Na + permeation. The shortcomings of the adopted computational methods should be clearly kept in mind when interpreting the present results. The convergence of the PMF calculations in complex molecular systems is already a well reported issue [23,24]. Here, the magnitude of the errors of the PMF was calculated by comparing different parts of the simulated trajectories. This approach is reliable if the simulated trajectories are extensive enough to sample the relevant ion configurations in the selectivity filter. However, the possibility that structural changes of the protein on a longer time scale may affect the calculated energy profiles cannot be ruled out. In contrast to K + -channels, the movements of the ions appear to be weakly coupled in Na + -channels, and the ions maintain their hydration shell. Therefore, the issue of the absolute binding energy values of the ions [11], that is, the energy of the ions in the protein relative to their energy in water should be less crucial than in the case of K + -channels. In the later, the water molecules of the ion's hydration shell are replaced with the carbonyl oxygen atoms of the protein. A particular critical point is the presence of ionisable residues in the close proximity of the permeation pathway (Glu177). The ionization state of these residues may have an important effect on the permeation properties, and need to be closely analyzed. In this respect, it is also important to remember that a non-polarisable force field has been used and polarization effects may be important in the present context [25,26,27]. Such calculations will be reported in the future.
Another crucial point in this type of studies that was clearly pointed out by Warshel et al. [28] is that 'the issue of ion selectivity cannot be resolved quantitatively without calculating the corresponding ion current, which is the actual direct observable of the system'. As they explained in their paper 'such a calculation should be able to convert the free energy profile of the system to corresponding time dependence of the ion permeation process'. In Figure 3. Potential of Mean force for conduction events with two Na + ions. Values on the x/y-axis correspond to distances along the pore axis between the permeating ions and the centre of the carbonyl oxygen atoms of residues Thr175 (d1) and Leu176 (d2) respectively. Counter lines are drawn every 1 kcal/mol. Snapshots from the umbrella sampling trajectories are shown for some significant configurations. Residues 175 to 178 of two opposite subunits are depicted in licorice representation, together with the permeating Na + ion (yellow). Water molecules closer than 5 Å (licorice representation) and 15 Å (grey lines) to the ion are shown. The minimum energy path between a minimum with one ion in the cavity and a minimum with one ion in extracellular solution is shown as a black line. The free energy along the minimum energy path is shown in Figure S3. doi:10.1371/journal.pcbi.1002476.g003 order to do this the present study should be extended, and the free energy obtained here should be coupled to Brownian dynamics simulations as those reported by Warshel et al. [28]. In such a way, under certain conditions that allow replicating experimental data and using the calculated free energy profiles as starting point, a detailed understanding of ion selectivity may then be possible.

Model Definition
Atomic coordinates of the NavAb channel were taken from the Protein Data Bank entry 3RVY [7]. Only residues 116 to 221 that correspond to the pore region of the channel were included in the model. Default protonation states were used for all the ionisable residues. N-and C-terminals were amidated and acetylated respectively. The channel was centred in the x-y plane and the permeation pathway was aligned with the z-axis. The aromatic belt defined by the amphipathic residues Trp195 was aligned with the upper layer of a pre-equilibrated bilayer of DOPC molecules. All lipid molecules closer than 2.0 Å to the protein atoms were removed. The system was solvated by ,20.000 water molecules, and 32 Na + ions and 24 Cl 2 ions were added. A Na + ion was positioned in site S CEN . Harmonic restraints were initially applied to the backbone atoms and to the Na + ion in the selectivity filter, and gradually removed during a one nanosecond period. The total production run was 40 ns.
MD trajectories were simulated with the version 2.7 of NAMD [29], using the CHARMM27 force field with CMAP corrections [30], and the TIP3P model for water molecules [31]. Standard parameters for K + and Na + in CHARMM27 force field were adopted. Simulations were performed in the NpT ensamble. Pressure was kept at 1atm by the Nose-Hoover Langevin piston method [32,33], with a damping time constant of 100 ps and a period of 200 ps. Temperature was kept at 300 K by coupling to a Langevin thermostat, with a damping coefficient of 5 ps 21 [33]. Electrostatic interactions were treated by the Particle Mesh Ewald algorithm, with grid spacing below 1 Å [34]. Van der Waals interactions were truncated at 12 Å , and smoothed at 10 Å . Hydrogen atoms were restrained by the SETTLE algorithm [35], which allowed a 2 fs time-step.

Free-Energy Calculation
Free-energy profiles for one and two ion conduction events were calculated by umbrella sampling [36]. The reaction coordinate for one ion conduction events was the distance along the z-axis between the permeating ion and the centre of the carbonyl oxygen atoms of Thr175. The same reaction coordinate was used for the bottom ion in the analyses of conduction events considering two ions. The reaction coordinate for the upper ion was the distance along the z-axis between the ion and the centre of the carbonyl oxygen atoms of Leu176. Harmonic potentials (force constant 10 kcal*mol 21 *Å 22 ) were applied to the reaction coordinates. Each umbrella sampling simulation consisted of 0.5 ns, and the first 50 ps were considered as equilibration period and discarded. 20 and 270 umbrella sampling simulations were performed respectively for conduction events considering one or two ions, moving the centres of the harmonic potentials in 1.0 Å steps. The starting configurations for the umbrella sampling simulations were defined using the final snapshot of the normal MD trajectory. The water oxygen atoms or the Na + ions closer to the centres of the reaction coordinates were selected. Then, the positions of the ions subjected to the harmonic potentials were switched with the positions of these molecules/ions closer to the centre of the harmonic potentials. The ions restrained in the umbrella sampling simulations were the ion already inside the selectivity filter in the normal MD simulation, plus an ion randomly taken from the extracellular solution in the case of two-dimensional umbrella sampling simulations. For simulations with K + ions, all the Na + ions were alchemically transformed to K + ions in the first snapshot, while for the calculation of the energy map for Na + /K + mixtures, only one ion (the innermost of the two restrained ions) was transformed to K + . Free-energy profiles were calculated by the Weighted Histogram Analysis Method [37]. The string method [38] was used to calculate the minimum energy path in free-energy profiles of conduction events with two ions ( Figure S3). The initial guess for the minimum energy path was the segment that connects the free-energy minimum with the innermost ion in the cavity, and the free-energy minimum with the outermost ion in the extracellular solution (200 equally spaced points were used for the discretization). The path evolved in the direction opposite to the free-energy gradient until the root mean square distance between two successive iterations was below 10 23 Å . The relative energies of the local minima and the energy barriers between different free-energy basins provided in the text are evaluated along this minimum energy path. In order to estimate the errors affecting the free-energy values, the umbrella sampling trajectories were divided into 3 separate sets of 150 ps each. Free-energy profiles and minimum energy paths were computed separately for each data set. The values reported in the text are the average and the standard deviation among these 3 data sets. The free-energy profiles and the minimum energy paths shown in Figures 2, 3, 4, 5 and S3 were calculated using the whole umbrella sampling trajectories after the equilibration period.
Note Added In Proof. While this paper was in revision similar results were reported by B. Corry and M. Thomas in J. Am. Chem. Soc., 2012, 134 (3), pp 1840-1846.

Supporting Information
Figure S1 Na + /K + displacement in the x-y plane. The displacement of Na + /K + ions with respect to the channel axis is shown for an ion in different positions along the pore axis: (A) from the carbonyl oxygen atoms of Thr175 to those of Leu176, (B) from the carbonyl oxygen atoms of Leu176 to 3 Å above them, (C) from 3 Å above the carbonyl oxygen atoms of Leu176 to 1 Å below the side chain oxygen atoms of Glu177, and (D) from 1 Å below to 1 Å above the side chain oxygen atoms of Glu177. All the umbrella sampling trajectories from the simulations with two ions were considered for the analysis, taking snapshots every 20 ps after the equilibration period. (TIF) Figure S2 Coordination number of Na + and K + ions. Oxygen atoms were considered part of the coordination shell of the ion if closer than 2.8 Å from a Na + ion or 3.2 Å from a K + ion. All the umbrella-sampling trajectories from the simulations with two ions were used for this analysis, taking snapshots every 20 ps after the equilibration period. The black line shows the total number of coordinating oxygen atoms, which is the sum of the oxygen atoms coming from: water molecules (blue line), protein backbone (pink line), Glu177 side chain (red line), and Ser178 side chain (green line). The average positions along the filter axis of the carbonyl oxygen atoms of Thr175, Leu176, and Ser178 are shown along the x-axis. (TIF) Figure S3 Free energy along the minimum energy path for Na + ions (A), K + ions (B), and K + /Na + mixture with K + in the innermost position (C). The minimum energy path was calculated between a local free-energy minimum with the bottom ion in the intracellular cavity and the top ion inside the selectivity filter, and a local free-energy minimum with the bottom ion inside the selectivity filter and the top ion in the extracellular solution. The path was discretized using 200 equally spaced points, and the minimum energy path was calculated using the string method. The values of the reaction coordinates, d1 and d2, along the minimum energy path are shown along the bottom and top xaxis for a subset of the local free-energy minima. d1/d2 are defined as the distance along the pore axis between the permeating ion in the inward/outward position and the centre of the carbonyl oxygen atoms of residues Thr175/Leu176. (TIF)