The selectivity of the Na+/K+-pump is controlled by binding site protonation and self-correcting occlusion

Abstract
eLife digest
Introduction
Results
Discussion
Materials and methods
References
Article and author information
Metrics

Abstract

The Na⁺/K⁺-pump maintains the physiological K⁺ and Na⁺ electrochemical gradients across the cell membrane. It operates via an 'alternating-access' mechanism, making iterative transitions between inward-facing (E₁) and outward-facing (E₂) conformations. Although the general features of the transport cycle are known, the detailed physicochemical factors governing the binding site selectivity remain mysterious. Free energy molecular dynamics simulations show that the ion binding sites switch their binding specificity in E₁ and E₂. This is accompanied by small structural arrangements and changes in protonation states of the coordinating residues. Additional computations on structural models of the intermediate states along the conformational transition pathway reveal that the free energy barrier toward the occlusion step is considerably increased when the wrong type of ion is loaded into the binding pocket, prohibiting the pump cycle from proceeding forward. This self-correcting mechanism strengthens the overall transport selectivity and protects the stoichiometry of the pump cycle.

https://doi.org/10.7554/eLife.16616.001

eLife digest

A protein called the sodium-potassium pump resides in the membrane that surrounds living cells. The role of this protein is to 'pump' sodium and potassium ions across the membrane to help restore their concentration inside and outside of the cell. About 25% of the body's energy is used to keep the pump going, rising to nearly 70% in the brain. Problems that affect the pump have been linked to several disorders, including heart, kidney and metabolic diseases, as well as severe neurological conditions.

The sodium-potassium pump must be able to effectively pick out the correct ions to transport from a mixture of many different types of ions. However, it was not clear how the pump succeeds in doing this efficiently. Rui et al. have now used a computational method called molecular dynamics simulations to model how the sodium-potassium pump transports the desired ions across the cell membrane.

The pump works via a so-called 'alternating-access' mechanism, repeatedly transitioning between inward-facing and outward-facing conformations. In each cycle, it binds three sodium ions from the cell’s interior and exports them to the outside. Then, the pump binds to two potassium ions from outside the cell and imports them inside. Although the bound sodium and potassium ions interact with similar binding sites in the pump, the pump sometimes preferentially binds sodium, and sometimes potassium. The study performed by Rui et al. shows that this preference is driven by how protons (hydrogen ions) bind to the amino acids that make up the binding site.

The simulations also suggest that the pump uses a ‘self-correcting’ mechanism to prevent the pump from transporting the wrong types of ions. When incorrect ions are present at the binding sites, the pump cycle pauses temporarily until the ions detach from the pump. Only when the correct ions are bound will the pump cycle continue again.

In the future, Rui et al. hope to use long time-scale molecular dynamics simulations to show the conformational transition in action. In addition, the 'self-correcting' mechanism will be directly tested by letting the wrong and correct ions compete for the binding sites to see whether the pump will transport only the correct ions.

https://doi.org/10.7554/eLife.16616.002

Introduction

The Na⁺/K⁺-pump is a primary active membrane transporter present in nearly all animal cells. It belongs to the P-type ATPase family, which utilizes the energy released from ATP hydrolysis to move ions against their concentration gradients across a membrane barrier. The ion species transported by the pump under physiological conditions are Na⁺ and K⁺. Like many other membrane transporters, the Na⁺/K⁺-pump works according to an 'alternating-access' ion transport mechanism, with the bound ions accessible from only one side of the membrane at a time. The consensus scheme of the pump cycle is known as the 'Albers-Post' cycle (Albers, 1967; Post et al., 1969). It involves two major conformations, E₁ and E₂, with inward- and outward-facing ion binding sites, respectively. In each cycle, the E₁ conformation binds three Na⁺ from the cytosol and exports them using the energy from ATP hydrolysis. After external release of Na⁺, binding of extracellular K⁺ promotes the structural transition to the occluded E₂ state, which finally imports two K⁺ as binding of ATP returns the pump to the E₁ conformation. E₁/E₂ conformational change and the accompanied electrogenic active transport are facilitated by the phosphorylation and dephosphorylation of an aspartate residue in the cytoplasmic domain, which is a hallmark of the P-type ATPase family (Axelsen et al., 1998).

The presence of the Na⁺/K⁺-pump is essential for excitability and secondary active transport. More than 40% of the energy produced in mammals is consumed by the Na⁺/K⁺-pump (Milligan et al., 1985). Although it is a machine designed for such a precise and important function, it has been shown that many cations, including congeners of K⁺ and some organic cations, can bind to the same sites used by the pump to bind and transport K⁺ ions (Mahmmoud et al., 2015; Ratheal et al., 2010). Competitive binding between Na⁺ and other cations at the cytoplasmic side of the membrane has also been observed (Schneeberger et al., 2001). An unsolved puzzle, therefore, is how the Na⁺/K⁺-pump is able to recognize and accept two K⁺ from the extracellular matrix, where Na⁺ concentration is much higher, and how it selectively binds Na⁺ from the cytoplasm to keep the pump cycle running forward.

Structural studies have provided tremendous insights into the function of the Na⁺/K⁺-pump, which consists of two obligatory subunits, α (catalytic) and β (auxiliary), and sometimes a tissue specific regulatory subunit from the FXYD protein family (Geering, 2006; Mercer et al., 1993). The transmembrane region of the α-subunit contains the ion binding sites within its 10 helices (called M1-M10). Recent crystal structures of the pump in its E₁ and E₂ states reveal the locations of the three Na⁺ binding sites in E₁ and the two K⁺ binding sites in E₂ (Kanai et al., 2013; Laursen et al., 2013; Morth et al., 2007; Shinoda et al., 2009). From structural alignments based on the heavy atom positions in the binding sites, it becomes clear that the binding pocket harboring sites I and II in E₁ overlaps with those in E₂ (Figure 1). Site III is only formed in E₁. It is located between the transmembrane helices M5, M6, and M8 and is thought to exclusively bind Na⁺ but appears to catalyze H⁺ transport (Poulsen et al., 2010; Ratheal et al., 2010) in a manner that presents a complex dependence on the concentrations of Na⁺, K⁺, and H⁺ (Mitchell et al., 2014). Previous studies have indicated that protonation at the ion binding pocket may play a role in the selectivity of these sites for K⁺ when the pump is in its E₂ state (Yu et al., 2011). Biochemical assays on the mutant D926N, which is often used as a surrogate for protonated D926, also show that it induces distinct effects on Na⁺ and K⁺ binding (Einholm et al., 2010; Jewell-Motz et al., 1993), suggesting a potential change in its protonation state upon the E₁/E₂ transition. Taken together, the evidence, although indirect, suggests the possibility of an E₁ specific protonation state that favors Na⁺ binding to the pump. The nature of this protonation state is, however, unclear.

Figure 1

Download asset Open asset

The ion binding sites in (A) Na₃·E₁·(ADP·Pi) (PDBID 3WGV) and (B) E₂(K₂) (PDBID 2ZXE) states.

Only the transmembrane helices M4, M5, M6, M8, and M9 from the α subunit are shown. Residues in the binding site are highlighted in stick presentation with those protonatable colored in yellow. Na⁺ (yellow) and K⁺ (magenta) ions are in spheres. The binding site number indices are presented on top of the ions. The view is from the extracellular side towards the intracellular side. The figure is produced with PyMOL (DeLano, 2002).

https://doi.org/10.7554/eLife.16616.003

In the present study, we started with the recently published crystal structure of the Na⁺/K⁺-pump in a partially occluded Na₃·E₁(ADP·Pi) state and used molecular dynamics (MD) simulations to show that there is a correlation between the binding pocket protonation and the Na⁺ selectivity in E₁. The binding sites in Na₃·E₁(ADP·Pi) were tested one by one by 'alchemically transforming' the bound Na⁺ into a K⁺ while keeping the other two sites occupied by Na⁺. A similar approach was previously used to study the selectivity in other conformational states of the pump (Yu et al., 2011) (see also Materials and methods section). The results show that the Na⁺ selectivity at all three sites is realized only when a specific set of binding pocket residues are protonated. This set of residues is, by comparison, different from that in the K⁺ selective E₂ state. The implication is that the protonation state and the selectivity of the pump are tightly coupled; when the pump undergoes a transition between E₁ and E₂, a protonation state switch occurs. The present findings also show that the effective selectivity of the pump is reinforced by a self-correcting mechanism, which prevents the occlusion step on either side of the membrane if the wrong type of ions were loaded into the binding sites. This mechanism ensures that counterproductive transport cycles do not occur.

Results

pKa of the binding pocket protonatable residues

The crystal structures of the pump in its E₁ (PDBID 3WGV [Kanai et al., 2013]) and E₂ (PDBID 2ZXE [Shinoda et al., 2009]) conformations show that there is a large structural overlap between the Na⁺ and K⁺ binding pockets. The sites I and II are in the main binding chamber formed by helices M4, M5, and M6, while site III, located in between M5, M6, and M8 (Figure 1), is Na⁺ exclusive and appears in E₁ only. The coordination of Na⁺ and K⁺ at sites I and II are similar. Many residues are found to coordinate both Na⁺ and K⁺ at these two sites in E₁ and E₂ (Table 1). Although the structural difference between the ion binding pockets in E₁ and E₂ may not be particularly large, the local physicochemical environment can display considrable variations. In particular, the latter could affect the pKa and the protonation states of the six protonatable residues in the binding pocket (Figure 1). Their pKa values calculated using PROPKA3.1 (Olsson et al., 2011) are listed in Table 2.

Table 1

Atoms coordinating the binding site ions in the crystal structures and from the MD simulations. O is the backbone carbonyl oxygen atom. OG and OG1 are the hydroxyl oxygen atoms in serine and threonine. OD1 and OD2 are the carboxyl oxygen atoms in asparate. OE1 and OE2 are the carboxyl oxygen atoms in glutamate. OH2 is the water oxygen.

https://doi.org/10.7554/eLife.16616.004

E₁	Site I				Site II				Site III
E₁	x-ray		MD		x-ray		MD		x-ray		MD
	A323	O	T772	OG1	V322	O	E779	OE1	Y771	O	Y771	O
	E779	OE1	T772	O	V325	O	D804	OD1	T774	O	T774	O
	D808	OD1	N776	OD1	E327	OE2	D808	OD1	Q923	OE1	Q923	OE1
			D808	OD1	D804	OD1	Water	OH2			D926	OD1
			D808	OD2	Water	OH2

E₂	Site I				Site II
E₂	x-ray		MD		x-ray		MD
	T772	O	S775	OG	V322	O	A323	O
	S775	OG	N776	OD1	V325	O	V325	O
	N776	OD1	D804	OD2	E779	OE2	E779	OE1
	D804	OD2			D804	OD2	D804	OD1

Table 2

pKa values of binding site titratable residues calculated from the crystal structures. The crystal structure resolution is given below the PDB ID.

https://doi.org/10.7554/eLife.16616.005

	E₁		E₂
	3WGV 2.8 Å	4HQJ 4.3 Å	2ZXE 2.4 Å	3B8E 3.5 Å
D804	5.9(6.2)*	11.1 (11.2)	3.7	0.8 (2.1)
D808	3.5(3.1)	3.7 (3.7)	5.8	6.8 (6.6)
D926	6.4(7.2)	5.6 (5.6)	8.9	7.4 (8.4)
E327	11.0(11.3)	5.7 (5.6)	8.3	10.8 (9.9)
E779	9.9(8.4)	7.4 (7.3)	10.7	9.6 (8.0)
E954	9.6(9.7)	9.2 (9.2)	10.3	10.7 (10.3)

*If two chains are present in the same asymmetric unit, the pKa of the same residue in the other chain is shown inside the parenthesis.

It is important to note that the pKa values calculated with empirical methods like PROPKA are sensitive to local structural perturbations. Even when the structures assume the identical conformational state and are taken from the same asymmetric unit from a crystal, the pKa values of the same residue can differ by more than one pH unit. For example, the pKa of E779 is 9.9 in chain A of the crystal structure 3WGV but the value is 8.4 in chain C, yet, the structural difference between the two chains is minimal (Table 2). Because of this, only the structures with the highest resolution for the E₁ and E₂ states of the pump, 3WGV and 2ZXE, were used to guide the protonation state assignments in order to avoid false assignments of protonation state originated from structural uncertainty in lower resolution crystal structures. Interestingly, the protonation states of E327 and D804 assigned based on the 3WGV and 4HQJ structures are reversed. According to the PROPKA results, E327 appears to be deprotonated and D804 protonated in 3WGV, and vice versa in 4HQJ. While this inconsistency could be due to the lower resolution in the crystal structure of 4HQJ, it is worth noting that these two residues are in close proximity from one another suggesting that a proton could be passed back and forth between them in the E₁ state.

Protonation states and E₁ binding pocket stability

Previous calculations have indicated that a specific set of residues has to be protonated for the E₂ binding sites to be K⁺ selective (Yu et al., 2011). Although the protonation states of D926 and E954 were not considered in the previous study because they lie outside of the two K⁺ binding sites, they are likely to be protonated in E₂ according to the predicted pKa (Table 2). The binding pocket protonation derived from the crystal structure of the pump in its partially occluded Na₃·E₁(ADP·Pi) state differs from that in E₂. According to their pKa values, D804, D808, and D926 should be deprotonated and the glutamates, E327, E786, and E954, should be protonated. Arguments based on pKa values predicted by PROPKA, however, have to be taken with caution, because the empirical method is highly sensitive to the local structural variation. For example, a deprotonated acidic side chain could have a PROPKA predicted pKa value at slightly higher than 7 because of the structural snapshot used to make the prediction. To seek a more robust assessment of these factors, molecular dynamics (MD) simulations were conducted to examine the structural stability of the binding pocket for different protonation configurations. One goal is to determine the protonation states of D804 and D926, both of which have a predicted pKa value within 1.5 pH unit to 7. The protonation state of D808 was also investigated since it affects the K⁺ selectivity of the binding sites in the E₂(K₂) state. A total of eight MD simulations were carried out (Table 3), starting from all possible protonation state configurations accessible by these three residues. The protocol for setting up these systems is given in details in the Materials and methods section.

Table 3

Summary of the all-atom simulation systems and the FEP/H-REMD reduced systems. The binding site residues E327, E779, and E954 were kept protonated in all the systems.

https://doi.org/10.7554/eLife.16616.006

Systems	Binding site residues			Simulation time (ns)
				MD	FEP/H-REMD
Wildtype
E1_S0	D804-	D808-	D926-	140	2 × 128
E1_S1	D804p	D808-	D926-	300	2 × 128
E1_S2	D804-	D808p	D926-	139	2 × 128
E1_S3	D804-	D808-	D926p	503	2 × 128
E1_S4	D804p	D808p	D926-	94
E1_S5	D804p	D808-	D926p	87
E1_S6	D804-	D808p	D926p	277
E1_S7	D804p	D808p	D926p	141
E2_S0	D804-	D808p	D926-	193	2 × 128
E2_S1	D804-	D808p	D926p	350	2 × 128
P-E2_S0	D804-	D808p	D926p	100	2 × 128
Mutants
E1_S1M	D804N	D808-	D926-	40	2 × 128
E1_S2M	D804-	D808N	D926-	40	2 × 128
E1_S3M	D804-	D808-	D926N	40	2 × 128
E2_S1M	D804N	D808-	D926p	40	2 × 128
E2_S2M	D804-	D808N	D926p	40	2 × 128
E2_S3M	D804-	D808-	D926N	40	2 × 128
P-E2_S1M	D804N	D808-	D926p	40	2 × 128
P-E2_S2M	D804-	D808N	D926p	40	2 × 128

Figure 2 shows the binding pocket conformation in the E₁ systems (Table 3) at the end of the all-atom MD simulations. Out of the three aspartates when two or more are protonated the binding pocket becomes unstable, showing a large deviation from the crystal structure either with one of the bound Na⁺ expelled from its binding site (systems E1_S5-7), or with a Cl^- ion attracted into the binding pocket (system E1_S4). These scenarios are not likely to happen in the normal function cycle of the pump. On the other hand, when the number of protonated aspartates is less than two the binding pocket remains structurally close to that in the crystal structure (systems E1_S0-3). The heavy atom root mean squared deviations (RMSD) between the crystal structure and the structure snapshots from the last 50 ns of the simulations are less than 1.6 Å.

Figure 2

Download asset Open asset

Comparison of snapshots at the end of the MD simulations (green) and the crystal structure of Na₃·E₁·(ADP·Pi) (PDBID 3WGV) (white).

The binding site residues are shown in stick presentation and the ions are shown as spheres. Binding site Na⁺ ions from the MD simulation snapshots are in yellow, and the crystal Na⁺ are in orange. A Cl^- ion has entered the binding site in system E1_S4 during the simulation and is shown in green. The view is from the extracellular side towards the intracellular side. The figure is produced with PyMOL (DeLano, 2002).

https://doi.org/10.7554/eLife.16616.007

Protonation state and the E₁ binding pocket selectivity

Analysis of the structural stability of the binding sites based on the MD simulations indicates that all of the four protonation states producing stable ion binding pockets may coexist when the pump is in the Na₃·E₁(ADP·Pi) state. The relative population of these states in the state Na₃·E₁(ADP·Pi), however, would differ from one another. The overall selectivity of the binding sites has contributions from all the states, and the predominant protonation state configuration should produce binding sites that are Na⁺ selective. The determination of this protonation state configuration starts from the protein-membrane systems generated by the restraint-free MD simulations. A reduced structural model of the binding pocket is derived from each of these systems (Table 3) and the binding free energy differences between Na⁺ and K⁺ ( $Δ Δ G_{N a \to K}$ ) at the binding sites are calculated (See Materials and methods). The value of $Δ Δ G_{N a \to K}$ reflects the affinity difference between K⁺ and Na⁺ binding as $Δ Δ G_{N a \to K} = Δ G_{k} - Δ G_{N a} = R T l n (K_{D, K} / K_{D, N a})$ . The results are plotted as the logarithm of the affinity ratio, $l n (K_{D, K} / K_{D, N a})$ , in Figure 3. The values of $Δ Δ G_{N a \to K}$ are shown in Table 4.

Figure 3

Download asset Open asset

Ion binding sites selectivity characterized by $l n (K_{D, K} / K_{D, N a})$ in states Na₃E₁·(ADP·Pi) (blue) and E₂(K₂) (red).

Sites I (square), II (circle), and III (triangle) are distinguished by their shapes. Values from the previous calculations with a smaller reduced region are shown as empty symbols. All the binding site glutamates (i.e., E327, E779, and E954) are kept protonated. The protonation states of the binding site aspartates are indicated below the plot.

https://doi.org/10.7554/eLife.16616.008

Table 4

The binding free energy difference ( $Δ Δ G_{N a \to K}$ ) at all the binding sites calculated from FEP/H-REMD simulations. The energy values are in kcal/mol.

https://doi.org/10.7554/eLife.16616.009

Systems	Binding sites
	I	II	III
Wildtype
E1_S0	0.1 ± 0.1	4.5 ± 0.1	1.8 ± 0.1
E1_S1	3.7 ± 0.2	1.7 ± 0.1	4.7 ± 0.1
E1_S2	-3.5 ± 1.0	3.8 ± 0.3	2.9 ± 0.6
E1_S3	-1.5 ± 0.2	4.7 ± 0.4	-0.5 ± 0.3
E2_S0	-4.4 ± 0.3	-6.1 ± 0.0
E2_S1	-3.5 ± 0.1	-6.1 ± 0.1
P-E2_S0	-1.0 ± 0.6	0.7 ± 0.4
Mutants
E1_S1M	1.4 ± 0.1	5.7 ± 0.0	2.4 ± 0.1
E1_S2M	4.3 ± 0.7	6.6 ± 0.6	5.2 ± 0.9
E1_S3M	-0.6 ± 0.1	1.7 ± 0.1	0.0 ± 0.0
E2_S1M	-3.3 ± 0.3	-5.7 ± 0.0
E2_S2M	-5.4 ± 0.1	-6.4 ± 0.1
E2_S3M	-3.8 ± 0.1	-7.2 ± 0.1
P-E2_S1M	0.8 ± 0.9	1.9 ± 0.4
P-E2_S2M	1.6 ± 0.5	1.9 ± 0.4

As shown in Figure 3, the protonation state of the aspartates greatly affects both the E₁ and E₂ binding sites selectivity. The protonation state at E₂(K₂) binding pocket proposed in the earlier study (Yu et al., 2011) is confirmed by the present calculations, in which the reduced system includes a larger region of the protein and more residues at the outskirt of the main binding pocket, including D926 and E954. It is evident that a different set of residues is protonated in the Na⁺ selective E₁ as compared to in the K⁺ selective E₂ (Figure 3). While the glutamates (i.e., E327, E779, and E954) remain protonated in both the Na₃·E₁(ADP·Pi) and E₂(K₂) states, the binding pocket aspartates take on opposite protonation states. In Na₃·E₁(ADP·Pi) D804 is protonated and both D808 and D926 are deprotonated, and their protonation states are reversed in E₂(K₂). Free energy perturbation (FEP) calculations on the protonation/deprotonation of D804 estimated a pKa value of 9.2 (See Materials and methods), further supporting its protonated form under physiological conditions in E₁.

Charge-neutralizing mutations and shift in selectivity

The impact of charge-neutralizing mutations of the key aspartate residues on the binding site selectivity was examined. D to N mutations have previously been used in experimental studies as a strategy to ascertain the possible charged state of these residues. The results of the free energy calculations, summarized in Figure 4, indicate that all the three sites remain Na⁺ selective in the partly occluded Na₃·E₁(ADP·Pi) in both D804N and D808N, whilst D926N causes the sites to lose most of the Na⁺ selectivity. Using Na⁺ dependent phosphorylation and ATP binding assays, it was shown that the cytoplasmic binding affinity for both Na⁺ and K⁺ is reduced in D804N and D808N (Jorgensen et al., 2001; Pedersen et al., 1997). The D804N mutation affects the cytoplasmic K⁺ and Na⁺ binding differently. The mutation’s impact on cytoplasmic Na⁺ binding is not as prominent as on K⁺ binding and D804N would appear more Na⁺ selective than the wildtype. This is reflected in the calculation results. The K_D,K/K_D,Na at sites I and III are reduced by a factor of ~50 from 476.6 and 2523.3 in the wildtype to 10.3 and 54.6 in the D804N mutant, but the sites are still Na⁺ selective. The selectivity of site II is increased to a much larger extent, from 17.0 in the wildtype to 13359.7 in the mutant. This is an ~1000 fold increase in selectivity and it implies that D804N is likely more Na⁺ selective than the wildtype pump. The D808N mutant also becomes more Na⁺ selective compared to the wildtype pump. The most dramatic increase in Na⁺ selectivity happens at site II as in the case of D804N. Interestingly, site I, which prefers to bind K⁺ in the wildtype with the protonated D808 (Figure 3), is now Na⁺ selective in D808N. It could be that the spatial packing is the major contributing factor to the site I ion selectivity. In the calculations, both sites II and III lose their Na⁺ selectivity in the D926N mutant. The shifting of selectivity towards K⁺ at site III in this mutant is a bit surprising, and worth commenting on. There are two possible explanations. First, substituting a negatively charged carbonyl oxygen with a bulkier but neutral –NH₂ at the D926 sidechain could prevent entrance of an ion to this site. Thus, the simulated conformation with an ion at this site could be energetically inaccessible. In other words, this is a metastable state with the absolute free energy of Na⁺ and K⁺ binding to this conformation equally prohibitive. An alternative explanation is that the D926N mutation alters the available conformational space accessible by helix M5 and this allows K⁺ to go into site III as suggested in reference (Kanai et al., 2013). The binding of this K⁺, however, prevents the further binding of Na⁺, and results in the compromised Na⁺ selectivity in this mutant. Even though the calculations show D926N with decreased Na⁺ selectivity, they do not explain why experimentally the D926N mutant has compromised Na⁺ binding ability in the absence of K⁺ (Einholm et al., 2010). The loss in selectivity seen in this D926N mutant, however, could account for the strong inhibition by high K⁺ on Na⁺/K⁺-ATPase activity (cf. Figure 2A in reference [Einholm et al., 2010]).

Figure 4

Download asset Open asset

Charge-neutralizing mutations and their impact on binding site selectivity.

The wildtype protein is colored black and the mutations in E₁ (blue) and E₂ (red) are colored differently. Sites I (square), II (circle), and III (triangle) are distinguished by their shapes. The empty symbols represent values calculated from the outward facing P-E₂ model.

https://doi.org/10.7554/eLife.16616.010

All the D to N mutations tested for the occluded E₂(K₂) state appear to have little impact on the K⁺ selectivity, in apparent contradiction with the experimental observations showing that both D804N and D808N display decreased selectivity for external K⁺ (Kuntzweiler et al., 1996; Pedersen et al., 1997). In a previous computational study, the D808N mutation was also shown to minimally affect the K⁺ selectivity in E₂(K₂). Discrepancies between the calculations and the selectivity inferred from biochemical experiments have been noted previously, though the reason was not clear. A plausible explanation could be that the experimentally measured selectivity is an outcome from the relevant states including both P-E₂ and E₂(K₂). However, any contribution from the P-E₂ state to the observed selectivity was not taken into account by the calculations because of the missing crystal structure of the outward-facing state. To test this idea, we generated a model structure of the P-E₂ state based on the homolog structures from the Ca²⁺ SERCA pump (see below). This model was then used to calculate the $Δ Δ G_{N a \to K}$ in the wildtype and the mutant pumps. The results show that both sites I and II have increased preference for external Na⁺ in the P-E₂ state of D808N (Figure 4), which explains previous discrepancies between the calculations and the experiments.

Ion selectivity along the pump cycle

Using the crystal structures of the Ca²⁺ SERCA pump as templates and a coarse grained transition pathway calculation method, ANMpathway (Das et al., 2014), we generated structural models of the Na⁺– and K⁺–loaded outward-facing P-E₂ state (see Materials and methods). The models appear to be structurally similar to the recently published P-E₂ like structures of the Na⁺/K⁺-pump (Gregersen et al., 2016; Laursen et al., 2013) and the vanadate-Inhibited, P-E₂ mimic of the Ca²⁺ SERCA pump (Clausen et al., 2016). The MD simulations of the models are also able to predict with considerable accuracy the experimentally measured gating charge upon ion binding (Castillo et al., 2015). The structural transition pathways leading to these models provide a view of the intermediate structures along the pump cycle. Using these models, we calculated the $Δ Δ G_{N a \to K}$ for the binding site of the intermediate state structures. The entire atomistic protein-membrane systems were used in the calculations. The results are presented in Figure 5 and Table 5. The calculations are valid as the $Δ Δ G_{N a \to K}$ values mirror those computed from the reduced systems with the same protonation states (Figure 3).

Figure 5

Download asset Open asset

The binding site ion selectivity along the pump cycle.

Sites I (square), II (circle), and III (triangle) are distinguished by their shapes. Different colors indicate whether there are two (red) or three (blue) sites that are included in the calculations. The conformational states are numbered and stamped along the pump cycle in the top panel. The protonation states of the aspartates are indicated below.

https://doi.org/10.7554/eLife.16616.011

Table 5

Selectivity in the form of at the binding sites along the pump cycle from state P-E₂·K₂ to P-E₂·Na₃. The energies are in kcal/mol.

https://doi.org/10.7554/eLife.16616.012

	*p/-	-/p	p/p	p/p	p/p
	s/-	-/s	s/p	p/s	s/s
(1) P-E₂.K₂	-0.7 ± 1.0	0.1 ± 0.5	-1.0 ± 0.6	0.7 ± 0.4	-0.6 ± 0.7
(2) Intermediate	0.0 ± 0.4	1.8 ± 1.3	-0.5 ± 0.3	-5.1 ± 0.3	-4.0 ± 0.9
(3) E₂(K₂)	-4.9 ± 0.8	-3.0 ± 0.5	-5.0 ± 0.2	-4.2 ± 0.1	-9.3 ± 0.5
(4) Intermediate	0.6 ± 0.7	-0.8 ± 0.3	2.0 ± 0.3	-2.7 ± 0.2	-1.6 ± 0.4
(5) K₂.E₁	1.5 ± 0.3	1.1 ± 0.2	1.1 ± 0.1	-3.8 ± 0.3	-1.2 ± 0.5
(6) K₂.E₁*	0.6 ± 0.5	0.7 ± 0.4	1.3 ± 0.2	-0.4 ± 0.8	0.1 ± 0.7

	*s/-/-	-/s/-	s/s/-	s/s/-	s/s/-	s/s/s	s/s/s	s/s/s	s/s/s	s/s/s	s/s/s	s/s/s
	p/-/-	-/p/-	p/s/-	s/p/-	p/p/-	p/s/s	s/p/s	s/s/p	p/p/s	p/s/p	s/p/p	p/p/p
(7) Na₃·E₁(ADP·Pi)	-1.0 ± 0.5	1.7 ± 0.2	0.7 ± 0.7	2.2 ± 0.2	3.0 ± 1.0	0.8 ± 0.5	2.6 ± 0.3	6.0 ± 0.3	2.1 ± 1.5	5.8 ± 1.1	7.4 ± 0.8	7.5 ± 2.0
(8) Intermediate	2.9 ± 1.0	-0.5 ± 1.2	1.5 ± 1.6	2.6 ± 1.3	2.5 ± 1.3	3.0 ± 0.7	0.6 ± 0.2	0.8 ± 0.5	4.0 ± 0.8	2.8 ± 1.1	2.7 ± 0.5	3.5 ± 1.5
(9) P-E₂·Na₃	2.3 ± 0.3	1.6 ± 0.5	2.1 ± 0.3	2.4 ± 0.3	4.6 ± 0.5	1.4 ± 0.3	1.7 ± 0.3	-2.2 ± 0.3	3.1 ± 0.8	-0.6 ± 0.6	-0.8 ± 0.6	-1.2 ± 1.2

*The top and bottom rows represent the starting and ending binding site ion configurations. A 'p' represents a K⁺ (potassium) ion and an 's' represents a Na⁺ (sodium) ion. The binding sites I, II, and III in this order are separated by '/'.

The results reveal a fascinating feature of the binding site selectivity along the pump cycle. Initially, the two K⁺ binding sites in the outward-facing model P-E₂·K₂ (state ① in Figure 5) appear to be non-selective or only weakly K⁺ selective. However, the selectivity for K⁺ over Na⁺ increases as the pump occludes to form the K⁺–bound occluded state E₂(K₂) (state ③). After the dephosphorylation of P-E₂, the binding pocket opens up to the cytoplasmic side. Accompanying this structural transition is the selectivity reversal at site I, switching the site from K⁺ to Na⁺ selective (③ to ④ in Figure 5). The protonation state change further shifts the ion selectivity at the sites in E₁ towards Na⁺. Changing the protonation state in K₂E₁ (④ to ⑤ in Figure 5) to the protonation state dominant in Na₃·E₁(ADP·Pi), i.e., state ⑥ in Figure 5 with protonated D804 and deprotonated D808 and D926, further reduces the K⁺ selectivity at both sites I and II. When the pump is in this state, site I is Na⁺ selective and site II is only weakly K⁺ selective. Together, the E₂ to E₁ structural transition and the protonation state change promote the release of K⁺ to the cytoplasmic side. Transient changes in the binding pocket protonation state upon occlusion/deocclusion are possible and could lead to variations in the magnitude of $Δ Δ G_{N a \to K}$ for the intermediate states. Nonetheless, the general trend in free energy changes upon occlusion/deocclusion should remain because the $Δ Δ G_{N a \to K}$ of the occluded state E₂(K₂) has a much larger magnitude than that of the open-access outward- and inward-facing states.

The free energy calculations corroborate the notion that the two K⁺ binding from the extracellular side is sequential and possibly cooperative. The sequential binding of extracellular K⁺ was initially demonstrated by Forbush (Forbush, 1987). Recently using crystallography and isotopic measurements Ogawa et al. presented strong evidence that the first K⁺ binds at site I and the second K⁺ at site II (Ogawa et al., 2015). Whether there is cooperativity upon extracellular K⁺ binding, however, is not clear from the crystal structures. Two electrode voltage clamp experiments have shown that the two extracellular K⁺ binding events are relatively independent in the absence of extracellular Na⁺, but there is positive cooperativity of K⁺ binding when extracellular Na⁺ ions are present (Jaisser et al., 1994). The current $Δ Δ G_{N a \to K}$ calculations of the two ion bound P-E₂ (Table 5) provide an interpretation of this phenomenon. With a K⁺ ion at site I, site II is relatively nondiscriminatory ( $Δ Δ G_{N a \to K}$ = 0.7 kcal/mol, p/p to p/s in state ① P-E₂.K₂, Table 5), but it becomes much more Na⁺ selective ( $Δ Δ G_{N a \to K}$ = 2.4 kcal/mol, s/s to s/p in state ⑨ P-E₂·Na₃, Table 5) when there is a Na⁺ ion bound at site I. In the presence of extracellular Na⁺, K⁺ ions have to first compete with Na⁺ to bind at site I in order for the subsequent K⁺ to bind, therefore resulting in the observed binding cooperativity in the presence of extracellular Na⁺.

A similar phenomenon is seen in the Na⁺ branch of the pump cycle. Initially, site I in the inward-facing Na₃·E₁(ADP·Pi) structure (state ⑦) is weakly Na⁺ selective when the other two sites are filled with Na⁺. However, the intermediate state ⑧ during the occlusion process shows that its site I is highly discriminating against K⁺ with a $Δ Δ G_{N a \to K}$ of 2.9 kcal/mol (Figure 5 and Table 5). Therefore, if a K⁺ is bound in this site, the free energy barrier toward occlusion would be much higher than when a Na⁺ ion was bound and the subsequent occlusion is not likely. In the case when a K⁺ ion replaces a Na⁺ at site II or III, although the energy barrier of occlusion is not as prominent, such a binding pocket ion configuration is energetically less favorable than the native configuration and its appearance is unlikely in the first place.

Unlike the K⁺ branch, the $Δ Δ G_{N a \to K}$ calculations along the Na⁺ branch offers limited insights on the sequence of Na⁺ binding from the cytoplasmic side. The results are indicative of the selectivity at the sites, but not the absolute binding affinity. Based on the calculations alone, it is not possible to determine which site is the first to bind cytoplasmic Na⁺. It could be site III, as suggested by Kanai et al., and this prepares the other two sites for the following Na⁺ binding (Kanai et al., 2013). Or, alternatively, the first two Na⁺ bind at sites I and II in the main binding pocket and the last Na⁺ enters, takes over site II, and pushes the two previously bound ions further, so that the ions in sites II and I now move to sites I and III, in a process reminiscent of the “knock-on” mechanism occurring in potassium channels (Hodgkin et al., 1955). A direct assessment of these proposed binding sequences will require further experiments.

Discussion

The results from the MD simulations support the notion that the protonation state of the binding pocket and its selectivity are closely related. Because the binding pocket in the Na⁺/K⁺-pump displays considerable flexibility, it is worth pausing to reflect on the possible mechanism that underlies this relation. While the selectivity of a very rigid binding site is first and foremost predetermined by its atomic geometry, the selectivity of a flexible binding site is strongly affected by local ion-ligand and ligand-ligand interactions. In such flexible systems, selectivity is controlled by the both number and the physicochemical properties of ion-coordinating ligands (Yu et al., 2010). For example, high-field ligands, such as deprotonated acidic side chains tend to favor Na⁺ binding and protonation revert those to low-field ligands, which tend to favor K⁺. Hence, in the Na⁺/K⁺-pump protonation is exploited to modulate selectivity by altering the electrostatic properties of several of residues in the binding pocket. This is also consistent with the results from our previous study on the ion selectivity in E₂(K₂), which concluded that changes in the electrostatic properties of the protonatable residues was the likely mechanism responsible for the K⁺ selectivity in the E₂ state of the pump (Yu et al., 2011). Even though the local structural rearrangements at the binding sites are small upon the change in protonation state, it is enough to change the physical properties of the coordinating ligand, thus giving rise to discernable differences in the ion selectivity. The crystal structures of the Na⁺/K⁺-pump in its Na₃.E₁·(ADP·Pi) and E₂(K₂) states show similarity in their binding pocket configurations (Figure 1), including the coordination patterns of the bound ions at sites I and II (Table 1). The heavy atom RMSD between the binding pocket residues is 2.5 Å and the structural difference remains after hundreds of ns of simulations. Empirical pKa and FEP calculations based on the MD simulation equilibrated structures indicate that the binding pocket glutamates (i.e., E327, E779, and E954) are likely protonated in both Na₃.E₁·(ADP·Pi) and E₂(K₂), although previous mutagenesis experiments showed that charge-neutralizing mutation E327Q affects pump function, possibly by altering ion-pump interactions and the kinetics of the occlusion/deocclusion reactions along the pump cycle (Jorgensen et al., 2001; Kuntzweiler et al., 1995; Li et al., 2006; Nielsen et al., 1998). The calculations also suggest that the protonation states of D804, D808, and D926 are different in Na₃·E₁·(ADP·Pi) and E₂(K₂). Among the three aspartates only D804 is protonated in order for all three sites to stay Na⁺ selective in Na₃·E₁·(ADP·Pi). The E₁ binding pocket devoid of ions carries a net charge of -2, consistent with that deduced from previous fluorescence studies (Schneeberger et al., 1999). Mutagenesis experiments are also consistent with the unlikely protonation of D808 and D926 when the pump is in E₁ trying to bind cytoplasmic Na⁺ (Jewell-Motz et al., 1993; Pedersen et al., 1997). The protonation states of these three residues are reversed in E₂(K₂) with D804 deprotonated and the other two protonated, a result that is supported by a previous computational study (Yu et al., 2011). The different protonation states for E₁ and E₂ also agree well with the 'four-site' model proposed by Skou and Esmann more than 30 years ago with the K⁺-bound E₂ state carrying two sidechain protons (H₂EK₂) and the Na⁺-bound E₁ state carrying only one (HENa₃) (Skou et al., 1980).

It is worth pointing out that the second proton (i.e., the proton on D926) in the K⁺ bound E₂ state must come from and return to the same side of the membrane during the pump cycle so that the net charge moved per cycle by the pump is one. It seems unlikely that this proton could come from the extracellular side because altering extracellular pH would then protonate or deprotonate D926, causing major changes in Na⁺ binding affinity and the maximum pumping turnover rate. This is not observed over an external pH range 9.6 to 5.6 (Mitchell et al., 2014; Vasilyev et al., 2004; Yu et al., 2011). Therefore, it is more likely that the D926 proton comes from the cytoplasmic side. This is supported by MD simulations of the P-E₂·Na₂ and Na₂·E₁·(ADP·Pi) revealing the existence of aqueous pathways connecting the cytoplasm and D926 (Figure 6). One of the water pathways is located between the helices M5, M7, and M8 (Figure 6B), similar to the C-terminal proton pathway previously proposed (Poulsen et al., 2010). A proton could enter through this pathway to protonate D926 and then leave through the same pathway during the E₂ to E₁ transition, or through an alternative path passing the main binding pocket along with the dissociating ions as seen in the Na₂·E₁·(ADP·Pi) simulation (Figure 6A).

Figure 6

Download asset Open asset

Water pathways from the cytoplasm to D926 in (A) Na₂·E₁·(ADP·Pi) and (B) P-E₂·Na₂.

The top (top) and side (bottom) views are shown. D926 are shown in sphere representation. Water path connecting the cytoplasm and the D926 are in surface representation colored in blue. Na⁺ (yellow) are shown as spheres.

https://doi.org/10.7554/eLife.16616.013

Under physiological conditions, the challenge faced by the Na⁺/K⁺-pump is to effectively pick out the correct ion species from a solution much more concentrated with other types of ions. Even when the binding affinity of the ion species being transported is a few orders of magnitude higher than that of the other ion ( $Δ Δ G_{b i n d}$ = 1 to 3 kcal/mol), the advantage in binding is undermined by the higher concentration of the competing ion. The pump overcomes this problem by raising the free energy barrier for the occlusion step in the presence of incorrect ions in the binding pocket, thus preventing the faulty transport (Figure 5 and Table 5). This self-correcting mechanism makes sure that the intracellular and the extracellular gates do not close (i.e., occlude), which is required for the pump cycle to go forward, unless the correct ion configuration is present at the binding pocket. Therefore, this is an inherent part of the ping-pong mechanism the pump uses to transport ions with high selectivity. Although other energetic barriers, like pathways the ions use to travel to the binding site, could contribute to the ion binding specificity, it is not likely the case here as both Na⁺ and K⁺ can access the binding pocket in the simulations of the E₁ and P-E₂ states without ions bound.

The concept of a self-correcting pumping mechanism sheds new light on a number of puzzling functional observations. For instance, it explains why E₁ displays a strong apparent affinity for K⁺ in competitive binding assays (Schneeberger et al., 2001), even though the pump still binds and occludes Na⁺ from the cytoplasmic side. The calculations show that all the sites in the Na₃·E₁·(ADP·Pi) state are selective for Na⁺, likely because this is a partially occluded state. A wrong combination of binding site ions would not have been occluded and reached such a state. In a fully inward-open conformation, the selectivity of the sites ought to be fairly weak. The high affinity for K⁺ observed in E₁ in these experiments is due to the backward occlusion leading to the K⁺-bound state E₂(K₂). The proposed self-correcting mechanism parallels the suggestion based on experiments that the phosphorylation and occlusion of E₁ requires 3 Na⁺ bound, and this increases the apparent affinity for Na⁺ in the normal pump cycle (Schneeberger et al., 2001). The self-correcting mechanism can also account for the different effects of Na⁺ and K⁺ congeners on the release rate of the occluded ⁸⁶Rb⁺ (a K⁺ congener) to the extracellular side upon the “backdoor” phosphorylation in the presence of Pi reported by Forbush (Forbush, 1987). This backward deocclusion is accelerated in the presence of Na⁺, and it follows a single phase, while a much slower deocclusion process with a second slow phase is observed in the presence of K⁺ or Rb⁺. According to the self-correcting mechanism, it is likely that upon the replacement of one of the two occluded ⁴²K⁺ or ⁸⁶Rb⁺ by cold K⁺ at one site, the pump succeeded to occlude again, resulting in the second slow phase. These explanations are sensible, in turn lending support to the proposed self-correcting mechanism. Given that the selectivity in the outward/inward facing states is after all not as strong as previously thought, such a mechanism must be in place.

This analysis suggests that an important component of the overall selectivity emerges from the increased free energy barrier associated with the occlusion process while the wrong type of ion is bound to a site. Although the pump would be kinetically more efficient by preventing the wrong ions from reaching the binding sites to begin with, the free energy calculations do not support the notion of highly selective binding sites for the open-access states. While such mechanism may seem inefficient because the pump must try to figure out that there is an issue with selectivity only after binding of several ions of the wrong type, it is important to realize that the Na⁺/K⁺-pump is not a particular fast molecular machine. The estimated turnover rate is less than 100 per second and decreases even further as the transmembrane voltage becomes more negative (Heyse et al., 1994). The implication is that the system has plenty of time, at the molecular level, to function near thermodynamic equilibrium. In other words, the pump has not been evolutionarily optimized to be a particularly fast molecular machine, but an energetically efficient one. Thus, even though the idea of a self-correcting ion selectivity mechanism seems counterintuitive and inefficient, it is consistent with the physical conditions under which the pump has to operate.

In summary, the present study highlights the tight coupling between the Na⁺/K⁺ selectivity of the binding sites, the protonation state of the coordinating residues, and the conformational state of the pump. The important functional consequence of such tight coupling is the necessity to have the correct type and number of ions in the binding pocket for the pumping cycle to proceed forward toward the occlusion step. Because the ion binding selectivity is strongly dependent upon the protein conformation, a self-correcting mechanism counteracting the effect of the ion concentration in the environment ensures an efficient function of the pump.

Materials and methods

Preparing the all-atom membrane-Na⁺/K⁺-pump simulation systems

Request a detailed protocol

The crystal structures, 3WGV (Kanai et al., 2013) and 2ZXE (Shinoda et al., 2009), representing the Na⁺/K⁺-pump in its Na₃·E₁(ADP·Pi) and E₂(K₂) states were used to build the simulation systems. The structure 3WGV contains two copies of the pump assembly in the asymmetric unit. Since the structural variations between the two assemblies are minimal, only the copy including chains A, B, and G was kept. Several small molecule ligands were co-crystalized with the pump in 3WGV, including an ADP, a $A l F_{4}^{-}$ ion, an oligomycin A, two Mg²⁺, three cholesterol, four Na⁺, and five 1,2-diacyl-sn-glycero-3-phosphocholine molecules for each copy of the pump structure. Among these ligands, oligomycin A was not included in the simulation system. The ion was replaced by a $P O_{4}^{3 -}$ and POPC molecules were used in place of the 1,2-diacyl-sn-glycero-3-phosphocholine. The structure of 2ZXE also contains small molecule ligands including cholesterol, K⁺, Mg²⁺, and $M g F_{4}^{2 -}$ . Similarly, a $P O_{4}^{3 -}$ was used to replace the $M g F_{4}^{2 -}$ ion. Since the Na⁺/K⁺-pump crystal structures were obtained from different organisms, the residue numbers differ slightly. For the sake of convenience, the numbering scheme in the newly resolved 3WGV from pig kidney was adopted in this study.

After removing irregularities from the PDB files, the Na⁺/K⁺-pump subunits were capped with acetylated N-termini and amidated C-termini. The ectodomain of the β-subunit was not included in simulation systems to reduce the computational cost. The orientation was chosen to be the same as in the OPM database (Lomize et al., 2006). At this stage, the protonation states of the binding site residues were assigned with the PATCH command in CHARMM (Brooks et al., 2009). Eight different protonation states were considered in Na₃·E₁(ADP·Pi) state systems and both the protonated and charged D926 were included in the E₂(K₂) state systems. This resulted in a total of ten systems before proceeding to the next step. Table 3 shows the system names and the associated binding site protonation. We used the membrane builder module in CHARMM-GUI (Jo et al., 2007, 2008, 2009; Wu et al., 2014) to generate POPC bilayers around the pump structures. Once this was completed, each protein-membrane complex was solvated by an equal molar mixture of KCl and NaCl. The final system had a combined cation concentration of 0.3 M (i.e., [K⁺] = [Na⁺] = 0.15 M). At the end the E1 system was 84 × 110 × 158 Å³ in size and contained ~138,000 atoms, while the dimension of the E2 system was 85 × 109 × 105 Å³ with ~134,000 atoms. Each completed system was subjected to a 675-ps equilibration with reducing restraints on the heavy atoms to relax the initially uncorrelated components, followed by a 10-ns unrestrained pre-production using the simulation package NAMD2.9 (Phillips et al., 2005). After the systems were well equilibrated, they were simulated longer using the special-purpose supercomputer Anton (Shaw et al., 2009), which is designed for long time scale MD simulations.

Generating the D to N mutant systems

Request a detailed protocol

Experimentally asparagine and glutamine are used as surrogates for protonated aspartate to study the effect of protonation. The effects of these charge-neutralizing mutants on the selectivity of the pump can be investigated computationally. The most interesting mutations in the context of this study are the single D to N mutations at the binding pocket, including residues D804, D808, and D926. The mutations were made on the crystal structures by replacing the proton on the OD2 atom in the aspartate with an –NH₂ amine group. D804N and D808N mutations were also made in the outward facing P-E₂ structural model (see below) to study how they affect the external K⁺ binding. The protonation states of the other titratable residues were determined with PROPKA. Each of these systems was equilibrated without any restraints for 40 ns. The system snapshot with the least structural deviation of the pump to the averaged structure during the simulation was used to generate the reduced system. The mutant systems are shown in Table 3.

Selectivity at the binding site

Request a detailed protocol

The absolute free energy of an ion i binding to a binding site inside a protein has the following form,

Δ G_{i, b i n d} = (G_{i, i n t}^{s i t e} - G_{i, i n t}^{b u l k}) + [- k_{B} T I n (F_{t} C^{\circ}) - G_{i, t r a n s}^{s i t e}] .

The difference in the first term, $G_{i, i n t}^{s i t e} - G_{i, i n t}^{b u l k}$ , represents the nonbonded interaction component of the free energy change upon moving the ion from the bulk solution to the binding site. The subtraction in the second term, $- k_{B} T I n (F_{t} C^{\circ}) - G_{i, t r a n s}^{s i t e}$ , reflects the lost of translational freedom. The translational freedom factor F_t in bulk solution can be evaluated analytically under a rigid rotor approximation (Deng et al., 2006). Its final form depends on the force constants and the equilibrium values in the distance and angle restraints applied on the ion and the surrounding protein atoms. Based on Equation (1), the selectivity of a binding site can be expressed as the binding free energy difference of two ion species. For example, the binding free energy change upon changing ion i to j at the binding site is

\begin{array}{ll} Δ Δ G_{i \to j} & = (G_{j, i n t}^{s i t e} - G_{j, i n t}^{b u l k}) - (G_{i, i n t}^{s i t e} - G_{i, i n t}^{b u l k}) + (G_{i, t r a n s}^{s i t e} - G_{j, t r a n s}^{s i t e}) \\ = Δ G_{i \to j}^{s i t e} - G_{i \to j}^{b u l k} - G_{i \to j}^{t r a n s} . \end{array}

A negative value of $Δ Δ G_{i \to j}$ indicates that ion j binds more favorably than ion i and the site is j selective. There are three terms to be evaluated in Equation (2). $Δ G_{i \to j}^{s i t e}$ and $Δ G_{i \to j}^{t r a n s}$ are calculated using the reduce binding site model, while $Δ G_{i \to j}^{b u l k}$ is computed using a water sphere with the impact from the bulk solution factored in with a boundary potential (see below).

Generating the reduced binding site systems

Request a detailed protocol

To generate the reduced binding site, the all-atom system prepared according to the procedures above was divided into an inner region and an outer region. The inner region was defined as residues and water molecules within 15 Å to the center of mass of the bound ions. Everything within this region was treated explicitly. An extended inner region was specified by a 3-Å thick shell continuing from the inner region outwards to create a smooth spherical dielectric cavity. Water molecules in this region were removed and their impact on the binding site was included implicitly. The coordinates of atoms in these extended region and those linked to them within three atomic bonds in the inner region were held fixed during the simulations. The rest of the system was considered the outer region, in which the atoms were removed and their impact on the inner region atoms was represented by the General Solvent Boundary Potential (GSBP) in the form of a solvent-shielded static field and a solvent-induced reaction field (Im et al., 2001). The reaction field arising from changes in charge distribution in the inner region was expressed in terms of a generalized multipole expansion using 11 spherical harmonic functions. Both the solvent-shielded static field and the reaction field matrix were independent of the inner region configuration, and therefore were calculated only once before further simulations using the finite-difference Poisson−Boltzmann (PB) method with the PBEQ module (Im et al., 1998) in CHARMM. In these calculations a dielectric constant of 1 was used for the inside of the protein within the inner and outer regions, whereas the rest of the system had a dielectric constant of 80. The atomic Born radii used in the PB calculations were determined by free energy calculations in explicit solvent (Nina et al., 1997). All these reduced systems were further equilibrated for 200 ps using Langevin dynamics at 303.15 K with a friction coefficient of 5 ps⁻¹ assigned to all non-hydrogen atoms. All bonds involving hydrogen atoms were fixed with the SHAKE algorithm (Ryckaert et al., 1977). Nonbonded interactions within 14.5 Å were accounted for explicitly, while everything else beyond this distance was treated with an extended electrostatics method (Stote et al., 1991). The simulation program CHARMM (Brooks et al., 2009) was used for the equilibration. Table 3 lists all the reduced systems.

Free energy perturbation simulations with boosting potential

Request a detailed protocol

It is known that the binding free energy calculated from FEP/MD simulations suffers from convergence issues because the residue sidechains at the binding site only sample a few of the several accessible rotameric states. This problem can be alleviated by introducing a replica-exchange scheme aiming at enhancing the sampling of sidechains (Jiang et al., 2009, 2010). This scheme allows exchange between the neighboring λ windows. Each λ window is coupled with a set of replica systems with a boosting potential of increasing strength used to accelerate the inter-conversion between sidechain rotameric states. We employed this hybrid FEP/H-REMD scheme implemented in the REPDSTR module in CHARMM (Jiang et al., 2010) to calculate the $Δ Δ G_{i \to j}$ .

The boosting potential for the χ₁ sidechain torsion of a binding site residue was obtained by fitting the potential of mean force as a function of the torsion χ, $𝒲 (χ)$ , with a cosine Fourier series in the form of

U_{B P} (χ) = \sum_{n = 1}^{n} K_{n} {1 + \cos [n (χ - χ_{0, n})]} .

The angle χ₁ is dihedral formed by the bonded atoms N, CA, CB, and CG. The total number of the cosine terms, N, varies from 3 to 6, depending on which one produce a better fit to the $𝒲 (χ)$ . A residue is considered in the binding site if any of its sidechain heavy atoms is within 4.5 Å to the bound Na⁺ or K⁺. Table 6 lists all the binding site residues included in the boosting potential calculations.

Table 6

Binding site residues and the fitted parameters k_i and χ_0,i.

https://doi.org/10.7554/eLife.16616.014

	k₁	χ_0,1	k₂	χ_0,2	k₃	χ_0,3	k₄	χ_0,4	k₅	χ_0,5	k₆	χ_0,6
M4
E327p	-0.932	53.11	1.055	61.08	1.844	62.05	-0.312	56.44	0.392	89.04	0.225	82.28
M5
Y771	-2.19	62.51	1.219	89.27	2.46	60.91	0.662	40.85	0.721	44.75	0.639	50.42
T774	-4.25	35.82	2.437	101.3	3.602	61.71	1.128	44.86	1.071	91.96	0.857	76.09
S775	3.068	-26.37	0.124	38.67	1.291	62.22
N776	2.664	-166.23	0.918	87.15	1.709	62.85	-0.53	71.79	-0.495	50.76	-0.406	75
E779p	-2.508	29.17	1.792	88.01	2.591	62.28	0.487	24.81	0.476	73.74	0.462	53.74
M6
D804-	3.895	-172.95	0.489	80.04	1.897	60.56	-0.197	82.91	-0.271	82.85	-0.193	79.47
D804p	2.968	-157.98	0.532	81.76	1.835	60.45	-0.247	31.05	-0.296	89.84	-0.278	83.37
N804	2.981	-160.18	0.371	84	1.451	59.42
D808-	3.215	-159.94	0.21	57.87	1.744	59.13
D808p	-2.298	55.97	1.204	66.88	2.44	62.53	0.746	37.74	0.715	55.74	0.627	63.72
N808	-2.344	29.85	1.771	63.6	2.654	64.99	0.914	39.19	0.73	58.37	0.63	61.7
M8
Q923	-2.943	33.4	2.217	87.95	2.671	61.48	0.603	74.2	0.535	30.65	0.551	44.23
D926-	4.361	-169.1	0.454	95.51	1.902	60.16	-0.283	92.02	-0.306	88.18	-0.215	85.37
D926p	-2.704	28.74	1.477	93.38	2.744	61.08	0.816	58.99	0.852	58.37	0.809	58.93
N926	3.674	-157.18	0.479	89.17	1.377	61.78	-0.625	82.59

Figure 7

Download asset Open asset

Fitting the potential of mean force $𝒲 (χ)$ (red) with the boosting potential, $U_{B P} (χ)$ (black).
https://doi.org/10.7554/eLife.16616.015

We performed umbrella sampling simulations to obtain the $𝒲 (χ)$ for each binding site residue. First, the entire transmembrane helix containing the residue of interest was taken from the crystal structure with its orientation kept the same as in the OPM database. This was to provide a proper secondary structure environment. Next, we replaced the sidechains of all the other residues on the helix to hydrogen atoms to remove their steric effects. An implicit planar membrane model, EEF1/IMM1 (Lazaridis, 2003; Schneeberger et al., 1999) was used in place of a membrane bilayer to solvate the helix. 72 umbrella windows along χ were set up with a window spacing of of 5°. A harmonic bias potential was applied in each window with a force constant of 100 kcal/mol/rad². After the windows were generated, each of them was simulated for 50 ns at 303.15 K using Langevin dynamics with a 2 fs time-step. All bond lengths involving hydrogen were fixed with the SHAKE algorithm (Ryckaert et al., 1977). The cutoff distance for nonbonded interactions was set to 11 Å. The simulation program CHARMM (Brooks et al., 2009) was used to conduct the simulations. The resulting distributions of χ were unbiased using the weighted histogram analysis method (WHAM) (Kumar et al., 1992). Figure 7 shows the $𝒲 (χ)$ and the fitted curves. The fitting parameters and are given in Table 6. An appropriate boosting potential can be easily applied using the CONS DIHE command in CHARMM, with the sign in front of $U_{B P} (χ)$ reversed. This will effectively cancel out the potential barrier between the rotamers of a given residue sidechain.

Before computing $Δ G_{i \to j}^{s i t e}$ and $Δ G_{i \to j}^{t r a n s}$ , restraints were set up to restrict the translational freedom of the bound ion of interest. First, three points within the protein region were picked out. They were used along with the position of the bound ion to set up the restraints. We employed the same protocol as in the Ligand Binder module (Jo et al., 2013) in CHARMM-GUI to select the three protein anchor points, p₁, p₂, and p₃. The relative position of the bound ions to the protein can be defined by the distance r from l₁ to p₁, the angle θ between l₁, p₁ and p₂, and the torsion ψ along l₁, p₁, p₂, and p₃. The translational restraint potential took the form of

u_{t r a n s} = \frac{1}{2} [k_{r} {(r - r_{0})}^{2} + k_{θ} {(θ - θ_{0})}^{2} + k_{ψ} {(ψ - ψ_{0})}^{2}] .

The force constant regarding distance (i.e., k_r) was set to 10 kcal/mol/Å², and the rest of the force constants were 200 kcal/mol/rad². The equilibrium values in u_trans were taken from the 200 ps equilibration of the reduced binding site.

$Δ G_{N a \to K}^{s i t e}$ was computed with a set of FEP/H-REMD simulations. The translational restraint, u_trans, was applied to restrict the translational freedom of the ion. The Na⁺ ion was changed alchemically to a K⁺ with an alchemical coupling factor, λ. 16 λ windows were used. Each λ window included 8 replicas with the strength of boosting potential scaled from 0 to 1. Exchange attempts were made every 0.2 ps and were only allowed between neighboring replicas with different boosting potential strengths in the same λ window and between the neighboring λ windows with zero boosting potential. A total of 128 replica systems were included in the calculations. Each replica system was simulated for 2 ns. To compute $Δ G_{N a \to K}^{t r a n s}$ at a given binding site, two sets of simulations were set up, one with Na⁺ at the binding site and the other with K⁺ for the calculation of $Δ G_{N a, t r a n s}^{s i t e}$ and $Δ G_{K, t r a n s}^{s i t e}$ . The translational and orientational restraints were decoupled gradually (Beglov et al., 1994) using a coupling factor λ. The λ window and the boosting potential setup were similar to those used in the $Δ G_{N a \to K}^{s i t e}$ calculations. All the FEP/H-REMD calculations were performed using the REPDSTR module in CHARMM (Jiang et al., 2010). The output energies from the zero boosting potential systems were collected and processed using WHAM (Kumar et al., 1992). To compute the standard error of $Δ Δ G$ , we divided the trajectories into 10 blocks. The standard deviations of the averaged $Δ Δ G$ from all the blocks are computed and reported in Tables 4 and 5 and Figures 3–5. Cautions should be taken when comparing the calculated and experimentally measured $Δ Δ G$ , because the latter usually contains contributions from multiple states of the pump while the calculated $Δ Δ G$ is done using a single state.

The bulk ion-water system and the calculation of $Δ G_{N a \to K}^{b u l k}$

Request a detailed protocol

The bulk system was generated by building a water sphere of 10 Å radius and centering it at the origin. The water sphere was made form pre-equilibrated water boxes with TIP3P water molecules. A Na⁺ ion was placed at the center of the sphere. Any water molecules overlapping with the ion were deleted. The influence of the remaining bulk was approximated by the spherical solvent boundary potential (Beglov et al., 1994). The system was equilibrated for 200 ps at 303.15 K. Other simulation options were kept the same as described in the reduced binding site system. During the equilibration the position of the ion was restrained using a weak harmonic bias potential with a force constant of 0.5 kcal/mol/Å². Once equilibrated, the Na⁺ was gradually changed to a K⁺ using the PERT module in CHARMM with 11 λ windows. The simulation time for each λ window was 1 ns. The resulting $Δ G_{N a \to K}^{b u l k}$ is 18.34 kcal/mol after unbiasing the energy outputs with WHAM. The calculated $Δ G_{N a \to K}$ for all the binding sites in the systems summarized in Table 3 are shown in Table 4.

D804 pKa calculations with explicit solvent

Request a detailed protocol

To further confirm the protonation state of D804, we evaluated its pKa shift with additional simulations in explicit solvent using the following formula:

Δ p K a = (Δ G_{s i t e}^{d e p r o t} - Δ G_{b u l k}^{d e p r o t}) / 2.303 k_{B} T .

$Δ G_{s i t e}^{d e p r o t}$ and $Δ G_{s i t e}^{d e p r o t}$ are the free energy change of aspartate deprotonation in the protein environment and in bulk water, respectively. The reduced system of E1_S1 (Table 3) was used to compute the deprotonation free energy of D804 at the ion binding site ( $Δ G_{s i t e}^{d e p r o t}$ ). To calculate the $Δ G_{b u l k}^{d e p r o t}$ , an aspartic acid residue with an acetylated N-terminus and an amidated C-terminus was put into a pre-equilibrated water sphere of 15 Å radius. As in the calculations of $Δ G_{N a \to K}^{b u l k}$ and $Δ G_{K \to N a}^{b u l k}$ , the impact of the bulk solution beyond the current system was incorporated by the spherical solvent boundary potential (Beglov et al., 1994). Alchemical FEP calculations were carried out in both systems. The λ windows were evenly spaced to gradually deprotonate the aspartate. The numbers of windows used in the $Δ G_{s i t e}^{d e p r o t}$ and $Δ G_{b u l k}^{d e p r o t}$ calculations were 24 and 10 respectively. Each window was equilibrated for 200 ps and contained 5 ns of sampling. The calculated $Δ G_{s i t e}^{d e p r o t}$ and $Δ G_{b u l k}^{d e p r o t}$ are −44.8 kcal/mol and −51.9 kcal/mol, respectively. Using Equation (5), the ΔpKa is 5.1, meaning the pKa of D804 is right shifted by 5.1 pK unit when it is in the Na⁺/K⁺ pump ion binding site. The final calculated pKa of D804 is 4.1 + 5.1 = 9.2 with 4.1 being the pKa value of an aspartate in bulk water (Berg et al., 2002), reinforcing the conclusion that D804 is protonated in the E₁ state of the pump.

Transition pathway calculations

Request a detailed protocol

Three conformational transition pathways were generated based on the shared homology between the Na⁺/K⁺-pump and the Ca²⁺ SERCA pump. One of such pathways connects states E₂(K₂) and P-E₂·K₂. The second connects states E₂(K₂) and K₂·E₁, and the third describes the transition between states Na₃·E₁-P and P-E₂·Na₃. First, A Cα-atom-only transition pathway for each of these was generated using the ANMPathway online server (http://anmpathway.lcrc.anl.gov/anmpathway.cgi) (Das et al., 2014) based on the SERCA pump crystal structures including 3B9B (Olesen et al., 2007), 1WPG (Toyoshima et al., 2004), and 1VFP (Toyoshima et al., 2004). Since both the Na⁺/K⁺-pump and the Ca²⁺ SERCA pump are P-type ATPases, it is reasonable to assume that they have similar pumping mechanisms and they share the same set of states along the pump cycle. Although crystal structures of the Na⁺/K⁺-pump are scarce, multiple high-resolution structures of the SERCA pump in different states are available (Karlsen et al., 2016). Among them, 3B9B represents the outward facing P-E₂ state. All the other SERCA pump structures were aligned based on their transmembrane region Cα positions to the Na⁺/K⁺-pump and the Cα RMSD were computed to find those that resemble the structural states captured by the Na⁺/K⁺-pump crystal structures 3WGV (Kanai et al., 2013) and 2ZXE (Shinoda et al., 2009). The two SERCA pump structures representing states Na₃·E₁ATP (3WGV) and E₂(K₂) (2ZXE) are 1VFP and 1WPG, respectively. Each generated coarse-grained pathway is made of a sequence of structural snapshots containing Cα atoms only. These snapshots are called images and are distributed at equal RMSD intervals. The images from the coarse-grained pathways were then used in the all-atom targeted molecular dynamics (TMD) simulations to guide each transition. The starting system in the first and second pathway TMD simulations was taken from the well-equilibrated MD simulation system E2_S1 starting from the crystal structure 2ZXE. The third pathway was realized using the MD simulation system E1_S1 built from the crystal structure 3WGV, with the catalytic D369 phosphorylated. The protocol of the TMD simulations followed those published before (Castillo et al., 2015)

Selectivity calculations along the pump cycle

Request a detailed protocol

Binding site ion selectivity was evaluated for nine systems representing different stages along the pump cycle. Among these systems are well-defined states: P-E₂·K₂, E₂(K₂), K₂·E₁, Na₃·E₁(ADP·Pi), and P-E₂·Na₃. The structures used for them are ① the generated P-E₂·K₂ model, ③ the equilibrated crystal structure 2ZXE, the equilibrated crystal structure 3WGV with ⑦ bound Na⁺ and with ⑤ bound Na⁺ replaced by K⁺ at sites I and II, and ⑨ the generated P-E₂·Na₃ model. Several intermediate states were also included. Intermediate state ② is between states ① P-E₂·K₂ and ③ E₂(K₂), taken from the TMD simulations of the first path after image 35. State ④ was taken at image 33 along the transition from state ③ E₂(K₂) to ⑤ K₂·E₁. Along the same vein, state ⑧ represented the intermediate at image 29 between states ⑦ Na₃·E₁(ADP.Pi) and ⑨ P-E₂·Na₃. The intermediates were taken near the midpoints of the transitions. An additional state ⑥ K₂·E₁^* with altered binding site residue protonation was also included to reflect the protonation state change upon the E₂ to E₁ transition. All the states included in the selectivity calculations and their relative placement along the pump cycle can be seen in Figure 5.

To evaluate the selectivity at a given binding site in a system, three simulations were performed with slightly varied Lennard-Jones (LJ) parameters for the ion of interest. In the first simulation, the parameters of the ion were unaltered. A linear interpolation of the LJ and the NBFIX parameters between Na⁺ and K⁺ was made. In the second simulation, the LJ and NBFIX parameters of the ion were replaced by those of a Na/K hybrid at the middle point of the interpolation. In the last simulation, the ion’s parameters were completely changed to represent the other ion type, either K⁺ or Na⁺, depending on the starting ion type. 10 ns trajectories were generated for each simulation using the simulation package OpenMM6.2 (Eastman et al., 2013). Energies were calculated from each trajectory using two parameter sets. The difference was fed into the WHAM equation (Kumar et al., 1992) and the free energy changes upon mutating the ion to the intermediate as well as upon mutating the intermediate to the other ion type were solved self-consistently. The sum of the two free energies gives the $Δ G^{s i t e}$ . The same strategy was used to compute the $Δ G^{b u l k}$ in water ( $Δ G^{b u l k}$ kcal/mol). The binding free energy difference is given by $Δ Δ G_{b i n d} = Δ G^{s i t e} - G^{b u l k}$ . The selectivity of binding site with different occupancies was also included in the calculations. The $Δ Δ G_{b i n d}$ are shown in Table 5.

NAMD simulation protocol

Request a detailed protocol

The initial relaxation and the restraint-free equilibration of the membrane pump systems were performed using the NAMD2.9 (Phillips et al., 2005) simulation package with the input scripts from the CHARMM-GUI membrane builder module. The NAMD simulation temperature was set to 303.15 K using Langevin Dynamics with a damping coefficient of 10 ps⁻¹ during the relaxation and 1 ps⁻¹ during the restraint-free equilibration. The van der Waals interactions were smoothly switched off at 10–12 Å by a force switching function (Steinbach et al., 1994) and the electrostatic interactions were calculated using the particle-mesh Ewald method with a mesh size of ~1 Å.

Anton simulation protocol

Request a detailed protocol

After a short equilibration with NAMD, each all-atom system listed in Table 3 were subjected to a hundred ns scale long production without any restraints using the special-purpose supercomputer Anton (Shaw et al., 2009). The volume of the periodic cell was kept constant and the temperature was set to 303.15K using the Nosé-Hoover thermostats (Martyna et al., 1992). The lengths of all bonds involving hydrogen atoms were constrained using M-SHAKE (Kräutler et al., 2001). The cutoff of the van der Waals and short-range electrostatic interactions was set to an optimal value suggested by the Anton guesser script, guess_chem. Long-range electrostatic interactions were evaluated using the k-space Gaussian split Ewald method (Shan et al., 2005) with a 64 × 64 × 64 mesh. The integration time step was 2 fs. The r-RESPA integration method (Tuckerman et al., 1992) was employed and long-range electrostatics were evaluated every 6 fs.

OpenMM simulation protocol

Request a detailed protocol

Simulations used to compute ion selectivity along the pump cycle were conducted using the simulation package OpenMM6.2 (Eastman et al., 2013). The constant pressure and temperature (NPT) ensemble was used for these simulations. The pressure was maintained at 1 atm with a Monte Carlo barostat, which attempts to adjust the system volume every 0.2 ps. The Langevin dynamics algorithm with a 1.0 ps⁻¹ friction coefficient was used to hold the simulation temperature, which is at 303.15 K. The lengths of all bonds involving hydrogen were fixed and the integration time step was set to 2 fs. A force switching function was applied from 10 to 12 Å to gradually turn off the van der Waals interactions and the particle-mesh Ewald method with an error tolerance of 0.0005 was used to evaluate the electrostatic interactions.

Simulation force field parameters

Request a detailed protocol

The same force field parameters were used in all other simulations in this study. The PARAM27 all-atom force field of CHARMM (MacKerell et al., 1998) with a modified version of dihedral cross-term correction (Mackerell et al., 2004) was used for the protein and the C36 lipid force field (Klauda et al., 2010) was used for POPC. Water molecules were modeled with the TIP3P potential (Jorgensen et al., 1983).

References

1. Albers RW
(1967) Biochemical aspects of active transport
Annual Review of Biochemistry 36:727–756.

https://doi.org/10.1146/annurev.bi.36.070167.003455
- Google Scholar
1. Axelsen KB
2. Palmgren MG
(1998) Evolution of substrate specificities in the P-type ATPase superfamily
Journal of Molecular Evolution 46:84–101.

https://doi.org/10.1007/PL00006286
- Google Scholar
1. Beglov D
2. Roux B
(1994) Finite representation of an infinite bulk system: Solvent boundary potential for computer simulations
The Journal of Chemical Physics 100:9050–9063.

https://doi.org/10.1063/1.466711
- Google Scholar
Book
(2002)
Biochemistry (5th ed)

W.H. Freeman.
- Google Scholar
1. Brooks BR
2. Brooks CL
3. Mackerell AD
4. Nilsson L
5. Petrella RJ
6. Roux B
7. Won Y
8. Archontis G
9. Bartels C
10. Boresch S
11. Caflisch A
12. Caves L
13. Cui Q
14. Dinner AR
15. Feig M
16. Fischer S
17. Gao J
18. Hodoscek M
19. Im W
20. Kuczera K
21. Lazaridis T
22. Ma J
23. Ovchinnikov V
24. Paci E
25. Pastor RW
26. Post CB
27. Pu JZ
28. Schaefer M
29. Tidor B
30. Venable RM
31. Woodcock HL
32. Wu X
33. Yang W
34. York DM
35. Karplus M
(2009) CHARMM: the biomolecular simulation program
Journal of Computational Chemistry 30:1545–1614.

https://doi.org/10.1002/jcc.21287
- Google Scholar
1. Castillo JP
2. Rui H
3. Basilio D
4. Das A
5. Roux B
6. Latorre R
7. Bezanilla F
8. Holmgren M
(2015) Mechanism of potassium ion uptake by the Na(+)/K(+)-ATPase
Nature Communications 6:7622.

https://doi.org/10.1038/ncomms8622
- Google Scholar
1. Clausen JD
2. Bublitz M
3. Arnou B
4. Olesen C
5. Andersen JP
6. Møller JV
7. Nissen P
(2016) Crystal structure of the vanadate-inhibited Ca(2+)-ATPase
Structure 24:617–623.

https://doi.org/10.1016/j.str.2016.02.018
- Google Scholar
1. Das A
2. Gur M
3. Cheng MH
4. Jo S
5. Bahar I
6. Roux B
(2014) Exploring the conformational transitions of biomolecular systems using a simple two-state anisotropic network model
PLoS Computational Biology 10:e1003521.

https://doi.org/10.1371/journal.pcbi.1003521
- Google Scholar
Software
1. DeLano WL
(2002)
The Pymol Molecular Graphics System

Schrödinger, Palo Alto, CA.
1. Deng Y
2. Roux B
(2006) Calculation of standard binding free energies: Aromatic molecules in the T4 Lysozyme L99A mutant
Journal of Chemical Theory and Computation 2:1255–1273.

https://doi.org/10.1021/ct060037v
- Google Scholar
1. Eastman P
2. Friedrichs MS
3. Chodera JD
4. Radmer RJ
5. Bruns CM
6. Ku JP
7. Beauchamp KA
8. Lane TJ
9. Wang LP
10. Shukla D
11. Tye T
12. Houston M
13. Stich T
14. Klein C
15. Shirts MR
16. Pande VS
(2013) OpenMM 4: A reusable, extensible, hardware independent library for high performance molecular simulation
Journal of Chemical Theory and Computation 9:461–469.

https://doi.org/10.1021/ct300857j
- Google Scholar
(2010) The rapid-onset dystonia parkinsonism mutation D923N of the Na+, K+-ATPase alpha3 isoform disrupts Na+ interaction at the third Na+ site
Journal of Biological Chemistry 285:26245–26254.

https://doi.org/10.1074/jbc.M110.123976
- Google Scholar
1. Forbush B
(1987)
Rapid release of 42K or 86Rb from two distinct transport sites on the Na,K-pump in the presence of Pi or vanadate

The Journal of Biological Chemistry 262:11116–11127.
- Google Scholar
1. Geering K
(2006) Fxyd proteins: New regulators of Na-k-atpase
American Journal of Physiology. Renal Physiology 290:F241–250.

https://doi.org/10.1152/ajprenal.00126.2005
- Google Scholar
(2016) Isolation, crystallizationand crystal structure determination of bovine kidney na(+),k(+)-atpase
Acta Crystallographica Section F Structural Biology Communications 72:282–287.

https://doi.org/10.1107/S2053230X1600279X
- Google Scholar
1. Heyse S
2. Wuddel I
3. Apell HJ
4. Stürmer W
(1994) Partial reactions of the Na,K-ATPase: determination of rate constants
The Journal of General Physiology 104:197–240.

https://doi.org/10.1085/jgp.104.2.197
- Google Scholar
1. Hodgkin AL
2. Keynes RD
(1955) The potassium permeability of a giant nerve fibre
The Journal of Physiology 128:61–88.

https://doi.org/10.1113/jphysiol.1955.sp005291
- Google Scholar
1. Im W
2. Beglov D
3. Roux B
(1998) Continuum solvation model: Computation of electrostatic forces from numerical solutions to the Poisson-Boltzmann equation
Computer Physics Communications 111:59–75.

https://doi.org/10.1016/S0010-4655(98)00016-2
- Google Scholar
(2001) Generalized solvent boundary potential for computer simulations
The Journal of Chemical Physics 114:2924–2937.

https://doi.org/10.1063/1.1336570
- Google Scholar
(1994) Modulation of the Na,K-pump function by beta subunit isoforms
The Journal of General Physiology 103:605–623.

https://doi.org/10.1085/jgp.103.4.605
- Google Scholar
1. Jewell-Motz EA
2. Lingrel JB
(1993) Site-directed mutagenesis of the Na,K-ATPase: consequences of substitutions of negatively-charged amino acids localized in the transmembrane domains
Biochemistry 32:13523–13530.

https://doi.org/10.1021/bi00212a018
- Google Scholar
(2009) Computation of absolute hydration and binding free energy with free energy perturbation distributed replica-exchange molecular dynamics (FEP/REMD)
Journal of Chemical Theory and Computation 5:2583–2588.

https://doi.org/10.1021/ct900223z
- Google Scholar
1. Jiang W
2. Roux B
(2010) Free energy perturbation hamiltonian replica-exchange molecular dynamics (fep/h-remd) for absolute ligand binding free energy calculations
Journal of Chemical Theory and Computation 6:2559–2565.

https://doi.org/10.1021/ct1001768
- Google Scholar
1. Jo S
2. Jiang W
3. Lee HS
4. Roux B
5. Im W
(2013) CHARMM-GUI Ligand Binder for absolute binding free energy calculations and its application
Journal of Chemical Information and Modeling 53:267–277.

https://doi.org/10.1021/ci300505n
- Google Scholar
1. Jo S
2. Kim T
3. Im W
(2007) Automated builder and database of protein/membrane complexes for molecular dynamics simulations
PLoS One 2:e880.

https://doi.org/10.1371/journal.pone.0000880
- Google Scholar
1. Jo S
2. Kim T
3. Iyer VG
4. Im W
(2008) Software news and updates - charnim-gui: A web-based grraphical user interface for charmm
Journal of Computational Chemistry 29:1859–1865.

https://doi.org/10.1002/jcc.20945
- Google Scholar
1. Jo S
2. Lim JB
3. Klauda JB
4. Im W
(2009) CHARMM-GUI Membrane Builder for mixed bilayers and its application to yeast membranes
Biophysical Journal 97:50–58.

https://doi.org/10.1016/j.bpj.2009.04.013
- Google Scholar
1. Jorgensen PL
2. Pedersen PA
(2001) Structure–function relationships of Na+, K+, ATP, or Mg2+ binding and energy transduction in Na,K-ATPase
Biochimica Et Biophysica Acta (BBA) - Bioenergetics 1505:57–74.

https://doi.org/10.1016/S0005-2728(00)00277-2
- Google Scholar
(1983) Comparison of simple potential functions for simulating liquid water
The Journal of Chemical Physics 79:926–935.

https://doi.org/10.1063/1.445869
- Google Scholar
1. Kanai R
2. Ogawa H
3. Vilsen B
4. Cornelius F
5. Toyoshima C
(2013) Crystal structure of a Na+-bound Na+,K+-ATPase preceding the E1P state
Nature 502:201–206.

https://doi.org/10.1038/nature12578
- Google Scholar
1. Karlsen JL
2. Bublitz M
(2016) How to compare, analyze, and morph between crystal structures of different conformations: The P-Type ATPase example
Methods in Molecular Biology 1377:523–539.

https://doi.org/10.1007/978-1-4939-3179-8_43
- Google Scholar
(2010) Update of the CHARMM all-atom additive force field for lipids: validation on six lipid types
The Journal of Physical Chemistry B 114:7830–7843.

https://doi.org/10.1021/jp101759q
- Google Scholar
(2001) A fast SHAKE algorithm to solve distance constraint equations for small molecules in molecular dynamics simulations
Journal of Computational Chemistry 22:501–508.

https://doi.org/10.1002/1096-987X(20010415)22:5<501::AID-JCC1021>3.0.CO;2-V
- Google Scholar
(1992) THE weighted histogram analysis method for free-energy calculations on biomolecules. I. The method
Journal of Computational Chemistry 13:1011–1021.

https://doi.org/10.1002/jcc.540130812
- Google Scholar
(1996) Asp804 and Asp808 in the transmembrane domain of the Na,K-ATPase alpha subunit are cation coordinating residues
Journal of Biological Chemistry 271:29682–29687.

https://doi.org/10.1074/jbc.271.47.29682
- Google Scholar
(1995) Glutamic acid 327 in the sheep alpha 1 isoform of Na+,K(+)-ATPase stabilizes a K(+)-induced conformational change
The Journal of Biological Chemistry 270:2993–3000.

https://doi.org/10.1074/jbc.270.7.2993
- Google Scholar
(2013) Crystal structure of the high-affinity Na+,K+-ATPase-ouabain complex with Mg2+ bound in the cation binding site
PNAS 110:10958–10963.

https://doi.org/10.1073/pnas.1222308110
- Google Scholar
1. Lazaridis T
2. Karplus M
(1999) Effective energy function for proteins in solution
Proteins: Structure, Function, and Genetics 35:133–152.

https://doi.org/10.1002/(SICI)1097-0134(19990501)35:2<133::AID-PROT1>3.0.CO;2-N
- Google Scholar
1. Lazaridis T
(2003) Effective energy function for proteins in lipid membranes
Proteins: Structure, Function, and Genetics 52:176–192.

https://doi.org/10.1002/prot.10410
- Google Scholar
(2006) The third sodium binding site of Na,K-ATPase is functionally linked to acidic pH-activated inward current
Journal of Membrane Biology 213:1–9.

https://doi.org/10.1007/s00232-006-0035-0
- Google Scholar
(2006) OPM: orientations of proteins in membranes database
Bioinformatics 22:623–625.

https://doi.org/10.1093/bioinformatics/btk023
- Google Scholar
1. MacKerell AD
2. Bashford D
3. Bellott M
4. Dunbrack RL
5. Evanseck JD
6. Field MJ
7. Fischer S
8. Gao J
9. Guo H
10. Ha S
11. Joseph-McCarthy D
12. Kuchnir L
13. Kuczera K
14. Lau FT
15. Mattos C
16. Michnick S
17. Ngo T
18. Nguyen DT
19. Prodhom B
20. Reiher WE
21. Roux B
22. Schlenkrich M
23. Smith JC
24. Stote R
25. Straub J
26. Watanabe M
27. Wiórkiewicz-Kuczera J
28. Yin D
29. Karplus M
(1998) All-atom empirical potential for molecular modeling and dynamics studies of proteins
The Journal of Physical Chemistry. B 102:3586–3616.

https://doi.org/10.1021/jp973084f
- Google Scholar
(2004) Extending the treatment of backbone energetics in protein force fields: limitations of gas-phase quantum mechanics in reproducing protein conformational distributions in molecular dynamics simulations
Journal of Computational Chemistry 25:1400–1415.

https://doi.org/10.1002/jcc.20065
- Google Scholar
(2015) K+ congeners that do not compromise Na+ activation of the Na+,K+-ATPase: hydration of the ion binding cavity likely controls ion selectivity
The Journal of Biological Chemistry 290:3720–3731.

https://doi.org/10.1074/jbc.M114.577486
- Google Scholar
(1992) Nosé–Hoover chains: The canonical ensemble via continuous dynamics
The Journal of Chemical Physics 97:2635–2643.

https://doi.org/10.1063/1.463940
- Google Scholar
(1993) Molecular cloning and immunological characterization of the gamma polypeptide, a small protein associated with the Na,K-ATPase
The Journal of Cell Biology 121:579–586.

https://doi.org/10.1083/jcb.121.3.579
- Google Scholar
1. Milligan LP
2. McBride BW
(1985)
Energy costs of ion pumping by animal tissues

The Journal of Nutrition 115:1374–1382.
- Google Scholar
(2014) Sodium and proton effects on inward proton transport through Na/K pumps
Biophysical Journal 106:2555–2565.

https://doi.org/10.1016/j.bpj.2014.04.053
- Google Scholar
(2007) Crystal structure of the sodium-potassium pump
Nature 450:1043–1049.

https://doi.org/10.1038/nature06419
- Google Scholar
(1998) Importance of intramembrane carboxylic acids for occlusion of K+ ions at equilibrium in renal Na,K-ATPase
Biochemistry 37:1961–1968.

https://doi.org/10.1021/bi972524q
- Google Scholar
1. Nina M
2. Beglov D
3. Roux B
(1997) Atomic radii for continuum electrostatics calculations based on molecular dynamics free energy simulations
The Journal of Physical Chemistry B 101:5239–5248.

https://doi.org/10.1021/jp970736r
- Google Scholar
(2015) Sequential substitution of K(+) bound to Na(+),K(+)-ATPase visualized by X-ray crystallography
Nature Communications 6:8004.

https://doi.org/10.1038/ncomms9004
- Google Scholar
1. Olesen C
2. Picard M
3. Winther AM
4. Gyrup C
5. Morth JP
6. Oxvig C
7. Møller JV
8. Nissen P
(2007) The structural basis of calcium transport by the calcium pump
Nature 450:1036–1042.

https://doi.org/10.1038/nature06418
- Google Scholar
(2011) PROPKA3: Consistent treatment of internal and surface residues in empirical pKa predictions
Journal of Chemical Theory and Computation 7:525–537.

https://doi.org/10.1021/ct100578z
- Google Scholar
(1997) Identification of Asp804 and Asp808 as Na+ and K+ coordinating residues in alpha-subunit of renal Na,K-ATPase
FEBS Letters 400:206–210.

https://doi.org/10.1016/s0014-5793(96)01381-6
- Google Scholar
1. Phillips JC
2. Braun R
3. Wang W
4. Gumbart J
5. Tajkhorshid E
6. Villa E
7. Chipot C
8. Skeel RD
9. Kalé L
10. Schulten K
(2005) Scalable molecular dynamics with NAMD
Journal of Computational Chemistry 26:1781–1802.

https://doi.org/10.1002/jcc.20289
- Google Scholar
1. Post RL
2. Kume S
3. Tobin T
4. Orcutt B
5. Sen AK
(1969) Flexibility of an active center in sodium-plus-potassium adenosine triphosphatase
The Journal of General Physiology 54:306–326.

https://doi.org/10.1085/jgp.54.1.306
- Google Scholar
1. Poulsen H
2. Khandelia H
3. Morth JP
4. Bublitz M
5. Mouritsen OG
6. Egebjerg J
7. Nissen P
(2010) Neurological disease mutations compromise a C-terminal ion pathway in the Na(+)/K(+)-ATPase
Nature 467:99–102.

https://doi.org/10.1038/nature09309
- Google Scholar
1. Ratheal IM
2. Virgin GK
3. Yu H
4. Roux B
5. Gatto C
6. Artigas P
(2010) Selectivity of externally facing ion-binding sites in the Na/K pump to alkali metals and organic cations
PNAS 107:18718–18723.

https://doi.org/10.1073/pnas.1004214107
- Google Scholar
(1977) Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes
Journal of Computational Physics 23:327–341.

https://doi.org/10.1016/0021-9991(77)90098-5
- Google Scholar
1. Schneeberger A
2. Apell HJ
(1999) Ion selectivity of the cytoplasmic binding sites of the Na,K-ATPase: I. Sodium binding is associated with a conformational rearrangement
Journal of Membrane Biology 168:221–228.

https://doi.org/10.1007/s002329900511
- Google Scholar
1. Schneeberger A
2. Apell HJ
(2001) Ion selectivity of the cytoplasmic binding sites of the Na,K-ATPase: II. Competition of various cations
Journal of Membrane Biology 179:263–273.

https://doi.org/10.1007/s002320010051
- Google Scholar
1. Shan Y
2. Klepeis JL
3. Eastwood MP
4. Dror RO
5. Shaw DE
(2005) Gaussian split Ewald: A fast Ewald mesh method for molecular simulation
The Journal of Chemical Physics 122:054101.

https://doi.org/10.1063/1.1839571
- Google Scholar
Conference
1. Shaw DE
2. Dror RO
3. Salmon JK
4. Grossman JP
5. Mackenzie KM
6. Bank JA
7. Young C
8. Deneroff MM
9. Batson B
10. Bowers KJ
11. Chow E
12. Eastwood MP
13. Ierardi DJ
14. Klepeis JL
15. Kuskin JS
16. Larson RH
17. Lindorff-Larsen K
18. Maragakis P
19. Moraes MA
20. Piana S
21. Shan Y
22. Towles B
(2009) Millisecond-scale molecular dynamics simulations on Anton
Proceedings of the Conference on High Performance Computing Networking.

https://doi.org/10.1145/1654059.1654126
- Google Scholar
(2009) Crystal structure of the sodium-potassium pump at 2.4 A resolution
Nature 459:446–450.

https://doi.org/10.1038/nature07939
- Google Scholar
1. Skou JC
2. Esmann M
(1980) Effects of ATP and protons on the Na : K selectivity of the (Na+ + K+)-ATPase studied by ligand effects on intrinsic and extrinsic fluorescence
Biochimica Et Biophysica Acta (BBA) - Biomembranes 601:386–402.

https://doi.org/10.1016/0005-2736(80)90543-X
- Google Scholar
1. Steinbach PJ
2. Brooks BR
(1994) New spherical-cutoff methods for long-range forces in macromolecular simulation
Journal of Computational Chemistry 15:667–683.

https://doi.org/10.1002/jcc.540150702
- Google Scholar
(1991) On the treatment of electrostatic interactions in biomolecular simulation
Journal De Chimie Physique Et De Physico-Chimie Biologique 88:2419–2433.

https://doi.org/10.1063/1.41345
- Google Scholar
1. Toyoshima C
2. Mizutani T
(2004) Crystal structure of the calcium pump with a bound ATP analogue
Nature 430:529–535.

https://doi.org/10.1038/nature02680
- Google Scholar
(2004) Lumenal gating mechanism revealed in calcium pump crystal structures with phosphate analogues
Nature 432:361–368.

https://doi.org/10.1038/nature02981
- Google Scholar
(1992) Reversible multiple time scale molecular dynamics
The Journal of Chemical Physics 97:1990–2001.

https://doi.org/10.1063/1.463137
- Google Scholar
(2004) Effect of extracellular pH on presteady-state and steady-state current mediated by the Na+/K+ pump
Journal of Membrane Biology 198:65–76.

https://doi.org/10.1007/s00232-004-0660-4
- Google Scholar
1. Wu EL
2. Cheng X
3. Jo S
4. Rui H
5. Song KC
6. Dávila-Contreras EM
7. Qi Y
8. Lee J
9. Monje-Galvan V
10. Venable RM
11. Klauda JB
12. Im W
(2014) CHARMM-GUI Membrane Builder toward realistic biological membrane simulations
Journal of Computational Chemistry 35:1997–2004.

https://doi.org/10.1002/jcc.23702
- Google Scholar
1. Yu H
2. Noskov SY
3. Roux B
(2010) Two mechanisms of ion selectivity in protein binding sites
PNAS 107:20329–20334.

https://doi.org/10.1073/pnas.1007150107
- Google Scholar
1. Yu H
2. Ratheal IM
3. Artigas P
4. Roux B
(2011) Protonation of key acidic residues is critical for the K⁺-selectivity of the Na/K pump
Nature Structural & Molecular Biology 18:1159–1163.

https://doi.org/10.1038/nsmb.2113
- Google Scholar

Article and author information

Author details

Huan Rui

Department of Biochemistry and Molecular Biology, The University of Chicago, Chicago, United States

Contribution
HR, Conception and design, Acquisition of data, Analysis and interpretation of data, Drafting or revising the article

Competing interests
The authors declare that no competing interests exist.

"This ORCID iD identifies the author of this article:" 0000-0002-6459-9871
Pablo Artigas

Department of Cell Physiology and Molecular Biophysics, Texas Tech University Health Sciences Center, Lubbock, United States

Contribution
PA, Analysis and interpretation of data, Drafting or revising the article

Competing interests
The authors declare that no competing interests exist.
Benoît Roux

Department of Biochemistry and Molecular Biology, The University of Chicago, Chicago, United States

Contribution
BR, Conception and design, Drafting or revising the article

For correspondence
roux@uchicago.edu

Competing interests
The authors declare that no competing interests exist.

"This ORCID iD identifies the author of this article:" 0000-0002-5254-2712

Funding

National Institutes of Health (P41GM103712-S1)

Huan Rui
Benoît Roux

National Science Foundation (MCB-1515434)

Pablo Artigas
Benoît Roux

National Institutes of Health (R15, NS081570-01A)

Pablo Artigas

National Institutes of Health (U54-GM087519)

Benoît Roux

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Acknowledgements

This work was supported by grant U54-GM087519 (to BR) and R15(NS081570-01A to PA) from NIH, and by grant MCB-1515434 (to PA and BR) from NSF. Anton computer time was provided by the National Center for Multiscale Modeling of Biological Systems (MMBioS) through Grant P41GM103712-S1 from the National Institutes of Health and the Pittsburgh Supercomputing Center (PSC). The Anton machine at PSC was generously made available by D.E Shaw Research.

Copyright

This article is distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use and redistribution provided that the original author and source are credited.

Metrics

8,121

views
576

downloads
36

citations

Views, downloads and citations are aggregated across all versions of this paper published by eLife.

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Article PDF

Open citations (links to open the citations from this article in various online reference manager services)

Mendeley

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

Huan Rui
Pablo Artigas
Benoît Roux

(2016)

The selectivity of the Na⁺/K⁺-pump is controlled by binding site protonation and self-correcting occlusion

eLife 5:e16616.

https://doi.org/10.7554/eLife.16616

Share this article

Cite this article

The ion binding sites in (A) Na3·E1·(ADP·Pi) (PDBID 3WGV) and (B) E2(K2) (PDBID 2ZXE) states.

Comparison of snapshots at the end of the MD simulations (green) and the crystal structure of Na3·E1·(ADP·Pi) (PDBID 3WGV) (white).

Ion binding sites selectivity characterized by ln(KD,K/KD,Na) in states Na3E1·(ADP·Pi) (blue) and E2(K2) (red).

Charge-neutralizing mutations and their impact on binding site selectivity.

The binding site ion selectivity along the pump cycle.

Water pathways from the cytoplasm to D926 in (A) Na2·E1·(ADP·Pi) and (B) P-E2·Na2.

Fitting the potential of mean force 𝒲(χ) (red) with the boosting potential, UBP(χ) (black).

Author details

Huan Rui

Contribution

Competing interests

Pablo Artigas

Contribution

Competing interests

Benoît Roux

Contribution

For correspondence

Competing interests

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

Categories and tags

Further reading

The ion binding sites in (A) Na₃·E₁·(ADP·Pi) (PDBID 3WGV) and (B) E₂(K₂) (PDBID 2ZXE) states.

Comparison of snapshots at the end of the MD simulations (green) and the crystal structure of Na₃·E₁·(ADP·Pi) (PDBID 3WGV) (white).

Ion binding sites selectivity characterized by $l n (K_{D, K} / K_{D, N a})$ in states Na₃E₁·(ADP·Pi) (blue) and E₂(K₂) (red).

Water pathways from the cytoplasm to D926 in (A) Na₂·E₁·(ADP·Pi) and (B) P-E₂·Na₂.

Fitting the potential of mean force $𝒲 (χ)$ (red) with the boosting potential, $U_{B P} (χ)$ (black).