Identification of Navβ1 Residues Involved in the Modulation of the Sodium Channel Nav1.4

Voltage-gated sodium channels (VGSCs) are heteromeric protein complexes that initiate action potentials in excitable cells. The voltage-gated sodium channel accessory subunit, Navβ1, allosterically modulates the α subunit pore structure upon binding. To date, the molecular determinants of the interface remain unknown. We made use of sequence, knowledge and structure-based methods to identify residues critical to the association of the α and β1 Nav1.4 subunits. The Navβ1 point mutant C43A disrupted the modulation of voltage dependence of activation and inactivation and delayed the peak current decay, the recovery from inactivation, and induced a use-dependent decay upon depolarisation at 1 Hz. The Navβ1 mutant R89A selectively delayed channel inactivation and recovery from inactivation and had no effect on voltage dependence or repetitive depolarisations. Navβ1 mutants Y32A and G33M selectively modified the half voltage of inactivation without altering the kinetics. Despite low sequence identity, highly conserved structural elements were identified. Our models were consistent with published data and may help relate pathologies associated with VGSCs to the Navβ1 subunit.

Sodium channels are multi heteromeric integral proteins composed of one pore-forming a subunit that associates with one or more auxiliary b subunit. The a subunit (260 kDa) contains 4 homologous domains (I-IV), and each domain contains 6 membrane-spanning sequences (S1-S6). These 4 domains form the voltage-dependent, cation-selective sodium channel [4]. Extracellular and intracellular loops of variable sizes connect the 6 transmembrane helices of the 4 domains [4].
The b1 subunit modulates the voltage dependence and kinetics of the pore-forming a subunit and can increase inactivation rates and promote faster recovery from inactivation in an isoform and cell-type specific manner. Previous studies have indicated that non-covalent extracellular interactions between a and b1 subunits mediate the functional modulation [6,7,8,9,10,11,12,13], however the molecular nature of the interactions remain poorly understood [2].
The Myelin protein P0 has been used to model the association of the a and b1 subunits in potassium channels [14]. Using multiple sequence alignments (MSAs) and evolutionary information, conserved amino acids can be identified. Often, these sites of conservation are important for the structure or function of proteins. In addition, family-dependent conserved positions that differ in the chemical type of amino acids can dictate binding to other proteins [15]. Studies that combine computational predictive methods, including MSA, knowledge-based evolutionary information, and secondary structure predictions yield an estimated 85% success rate for predicting amino acids involved in protein-protein associations [16].
Early studies from our laboratory, using oocytes co-injected with b1-b2 chimeras, were consistent with previous findings that indicate the extracellular domain of b1 accounts for the regulation of Nav1.4 channels [13,17]. This b1 domain belongs to the antibody variable domain family of Ig-like proteins [12]. In the present study, we combine electrophysiology, site-directed mutagenesis, and in silico methods to investigate the structure-function relationships of the VGSC b1 on the Nav1.4 a isoform. To assess the importance of a putative intramolecular disulphide bridge and a salt bridge in Navb1 subunit we generated alanine substitutions C43A and R89A that would disrupt these bonds. According to our model, the turn of a loop on b1, may also be important in the modulation of the sodium channel. To test this hypothesis we neutralize the aromatic side chain of Y32 and introduce a bulky side chain in G33 at this loop turn with the mutants Y32A and G33M, respectively and tested their effects on Nav1.4 channels.

Site-directed Mutagenesis and Heterologous Expression
The rNav1.4 cDNA cloned into the pGW1H vector was kindly provided by Eduardo Marban (John Hopkins University). Several point mutants were generated on the rNavb1 construct (Scnb1: Q00954) using the QuikChangeH II XL site-directed mutagenesis kit (Stratagene, La Jolla, CA, USA). The integrity of the PCR products were verified by electrophoresis and transformed by heat shock into ultracompetent One ShotH bacteria for amplification and nick repair. After antibiotic selection and culture, the mutated plasmids were obtained using the HiSpeedH Plasmid Purification kit (Qiagen, Mexico city, Mexico). Successful mutagenesis was confirmed by DNA sequencing with the Applied Biosystems 3730 DNA Analyzer (UNAM, Cuernavaca, Mexico).
All surgical procedures were performed in accordance with the Guide for the Care and Use of Laboratory Animals from the Mexican Council for Animal Care (Norma Oficial Mexicana NOM-062-ZOO-1999) and the National Institutes of Health Guide for the Care and Use of Laboratory Animals. This study was approved by the internal ethics committee of the Instituto de Fisología of the Universidad Autónoma de Puebla. All efforts were made to reduce the number of animals used and minimize animal suffering.

Electrophysiological Recording and Data Analysis
Two electrode voltage-clamp recordings were performed at room temperature (20-22uC) using an OC-725 C amplifier (Warner, New Haven, CT). The electrodes were pulled on a horizontal P-97 puller (Sutter Instruments, Novato, CA) and filled with 3 M KCl with a resistance of 0.6 to 1.2 MV. The sodiumcurrent signals were filtered at 1 kHz, digitized at a sampling rate of 10 kHz by a Digidata 1200 analogue-to-digital converter (Axon Instruments, Foster City, CA), and stored on a computer for analysis with pClamp software Version 8.02 (Axon Instruments, Foster City, CA). The sodium currents (I Na ) were generated by step depolarisations from a holding potential of 2100 mV at 0.1 Hz, unless otherwise indicated. To minimize voltage-clamp errors, only oocytes with minimal leak (#0.1 mA) were used in the study. Capacitive currents were manually compensated. Sodium conductances (gNa) were obtained from the peak currents generated by step depolarisations of 10 mV increments (30 ms durations) from a holding potential of 2100 mV to 0 mV and calculated with equation: where Em was the membrane potential, I Na was the current amplitude, and E eq was the equilibrium potential calculated for each cell. Steady-state inactivation data were obtained using a two-pulse protocol. First, a variable voltage conditioning pulse (from 2120 mV to 0 mV, 1000 ms duration) was given from a holding potential of 2100 mV and a gap of 2 ms. Second, a test pulse (25 ms duration) at 220 mV was given. The peak current from the second pulse was plotted as a function of the conditioning pulse potential. Activation and steady-state inactivation curves from every group were fitted using the Boltzmann distribution equation: where V was the potential of the voltage pre-pulse, V 1/2 was the half voltage of inactivation, dx was the slope, and A was a residual linear component. The time course of inactivation data from the peak current at 220 mV was fitted to a single exponential equation: where A1 was the relative fraction of current inactivation, t was the time constant, and x was the time. The recovery from inactivation was examined using a 500 ms conditioning pulse at 220 mV from a holding potential of 2100 mV followed by a variable recovery interval (Dt = 1 to 10000 ms) and a test pulse at 220 mV. Recovery data from each cell were fitted by a double exponential equation: where A 1 and A 2 were relative fractions of the recovery currents, t1 and t2 were time constants, and y0 was the amplitude of the steady-state component. The decay current after repetitive stimulation was determined by applying trains of 30 ms pulses from a holding potential of 2100 mV to a test potential of 220 mV at 1 Hz. Current amplitudes were normalized with respect to the first pulse. All the currents were analyzed using the pClamp version 10.2 software (Axon Instruments, Foster City, CA). Values were shown as mean 6 SEM. Statistical comparisons between mutant and wild type mean values were performed using unpaired Student's t-test. Comparisons of more than two mean values were performed using one-way analysis of variance (ANOVA), and pairwise multiple comparisons were performed using the Holm-Sidak method. The graphs were built and fitted using Sigmaplot 11.0 (SPSS, Inc., Chicago, Il) and Origin 8.02 (OriginLab Corp., Northampton, MA, USA).

Homology and Molecular Modeling
Various atomic models of rNavb1 were generated including fully and partially unattended homology models using three approaches: (1) I-TASSER [19], (2) Geno3D [20] and (3) CPHmodels-3.0 [21]. Model refinements were performed under the AMBER Force Field (FF) in VEGA ZZ 2.3.2 [22] and the UCSF Chimera [23]. The stereochemical qualities of the models were assessed with the VADAR server [24]. For model comparisons, sequential hydropathy analyses were performed using the MPex 3.2 software [25], and secondary structure predictions were performed using the web-based Jpred3 [26]. An inter-species MSA of b1 was performed to identify highly conserved domains involved in VGSC modulation. Sequence alignments were obtained using several algorithms and best scores were obtained using Multiple Sequence Comparisons by Log-Expectation MUSCLE [27]. Three-dimensional Smith-Waterman alignments using the BLOSUM-62 matrix with weighted secondary structure similarities (30%) were used to superpose structures [28]. Electrostatic surface potential (ESP) was calculated with Chimera under AMBER FF and coloured with the Coulombic scale for display with equation: where Q was the potential, q was the atomic partial charges, and d was the distance from the atoms i. The dielectric constant e = 4 was used to represent screening by the medium or solvent. For the purpose of this study the residue numbers included the signal peptide. However, this region is cleaved to form the quaternary structure of the membrane spanning protein, and it has been excluded from models and simulations.

Results
Several homology models of b1 were generated using the crystallized myelin protein P0 as the template, which had the highest degree of homology to b1. This protein associates with the VGSC subunits without affecting its electrophysiological behaviour [10]. The rat sequence of myelin protein P0 (PDB code: 1neu) was 25% identical to the extracellular Ig-like domain of rb1 with 63% positive substitutions of 113 paired target residues comprising 80% of the entire extracellular domain. Figure 1A shows the structural model of the extracellular domain of the sodium channel b1 subunit, and figure 1 B shows this three-dimensional representation of Navb1 (in cyan) superposed to its template (in gray) and 3 Ig-like crystal structures including a Tyrosine kinase, a T-cell receptor, and an antibody (PDB codes: 1wwb_X, 1hxm_B, 1kxq_H, respectively). Based on structure, knowledge, and interspecies MSAs ( Figure S1), 2 sites that may form intramolecular bonds critical to the modulation of the sodium channel were identified. First, cysteine 43 is conserved across species in b1 and the model of the tertiary structure of b1 indicated that it lies close to C21, a cysteine near the N-terminus. Previous data indicates the cysteine corresponding to C43 forms a disulphide bond in b3 [29]. Second, highly conserved arginine 89 was proposed to form an exposed intramolecular salt bridge with an adjacent acidic residue. To assess the importance of these residues, we generated alanine substitutions C43A and R89A and tested their effects on Nav1.4 channels.
Our modelling results indicated that Arg89 formed a salt bridge with Glu87 (,2.5 Å distance). Alternatively, it may interact with the carbonyl main chain of Val111. The equivalent position of R89 on the template (Arg69) also formed a salt bridge.
It has been acknowledged that low homology crystal structures conserve relevant three-dimensional information [30], which may be consistent with this case. The superposition of other Ig-like crystal structures indicated the potential to form an electrostatic bond at this position was conserved through divergent evolution. The residue at position 89 in b1 was always a basic residue (arginine or lysine, Figure 1B stick atom representation) despite very low sequence identity (,25%). In particular, two proteins showed salt bridges between arginine and aspartate or glutamate within a range of 5 Å . Finally, the Smith-Waterman 3D alignment yielded score values of 25%, 16%, and 28% homology for the three Ig-like crystal structures.
Wild type (WT) b1 accelerated the inactivation and the recovery from inactivation, and shifted the voltage dependence of activation and inactivation to more hyperpolarized potentials. Figure 2 showed representative current traces for the Nav1.4 channel stably co-expressed on frog oocytes alone and with the b1 subunit and the point mutants C43A and R89A. The effects on the current decay were disrupted by the alanine substitutions at Cys43 and Arg89 (Figure 2). b1-C43A Disrupts Modulation of the Voltage Dependence on Nav1.4 Channels Figure 3 shows the hyperpolarizing shifts that b1 induced on activation and steady-state inactivation on the Nav1.4 channel. This modulation was disrupted by the C43A mutant but not by the R89A mutant. The half voltage of activation and inactivation for C43A was significantly different from WT b1 (Table 1) constituting a loss-of-function effect on the voltage dependency.
b1-C43A and b1-R89A Delayed the Inactivation and Recovery of Nav1.4 Channels In Figure 4 A, we illustrate the recovery from inactivation. WT b1 generated a current recovery of approximately 50% after 1 ms. R89A and C43A slowed this process by 90% y 95% respectively (P,0.05). The recovered fraction of R89A was different from WT after a Dt of 10 ms (P = 0.05), however C43A was significantly different until a Dt of 100 ms (P,0.05). Data from each oocyte from each group was fitted individually and the means of the time constants were compared. Only the time constant of recovery from fast inactivation (t fast ) of C43A differed from Nav1.4+b1 (Table 2).
We analyzed the effects of the mutants on current decay and found that although both mutants decrease the b1-induced acceleration of inactivation (Figure 4 B), R89A yielded a statistically significant inactivation time constant, and C43A induced an irregular decay that did not fit an exponential equation. This behaviour was similar to the behaviour of the Nav1.4 a subunit alone (Table 2).

b1-C43A Induced a Use-dependent Decay
The presence of b1 also allows some isoforms of the sodium channel to recover from high frequency depolarisations. As shown in Figure 4 C, the current from the a subunit of Nav1.4 alone declined 30% after the second pulse at 1 Hz and exhibited a cumulative rundown up to 60% from control after repetitive peak voltage stimulations. Mutant C43A induced a 10% decline after the second pulse and had a cumulative rundown of 20% from control as opposed to R89A and WT b1. At 5 Hz, we observed the same 20% rundown for the C43A mutant and no effects for WT b1 or the R89A mutant (data not shown).

Computed Physicochemical Effects of b1 Mutants
To investigate the changes in potential energy produced by each mutant, we calculated the ESP of b1. In Figure 5A, the local increase in acidity produced by the neutralizing substitution of  Arg89 for alanine and the resultant loss of the proposed intramolecular salt bridge with Glu87 were shown. Importantly, these electrostatic changes did not affect the voltage dependence of either activation or inactivation. Our model predicted that the C43A mutant produced a local increase in polarity caused by the release of the thiol group of C21 (Figure 5 B, sphere atom representation), however the net effect of this mutant reflected a greater increase in entropy and reduction of rigidity at the N-terminal loop, which seemed to be a critical participant in the modulation of voltage dependence on channel activation and inactivation.
To further evaluate the role of this loop in the b1-induced modulation of the channel, we designed two point mutants in the first turn in the C-terminus direction. We neutralized the aromatic polar side chain of Y32 by an alanine substitution (Y32A), and tested the introduction of a large hydrophobic side chain with the mutant G33M. Each single point mutation selectively disrupted the voltage dependence of inactivation. Y32A and G33M significantly shifted the half voltages of inactivation to more depolarized potentials compared with WT (Table 1). Importantly these mutants did not modify the voltage dependence of activation, the time constants of inactivation, or the recovery from inactivation (Table 1 and 2). Our model showed that the alanine substitution of the Tyr43 aromatic ring uncovered a nearby acidic residue (Asp148) involved in a hydrogen bond network with the Tyr43 hydroxyl group. The substitution of the next residue by methionine (G33M) introduced a bulky hydrophobic side chain producing the same selective loss-of-function effect as Y32A ( Figure S2).

Discussion
The main goals of combining site-directed mutagenesis, electrophysiological assays, and modelling studies were to identify residues in the b1 subunit that play roles in VGSC current modulation and acquire hints about the location of the interface between the subunits. Voltage-gated ion channels are highly flexible complexes that undergo conformational changes to transition between resting, activated, inactivated, and closed states to respond to membrane depolarization [31]. Therefore, we considered the possibility of multiple binding modes between the a and b1 subunits.
In this study, the selective effects of the b1 point mutants suggested the existence of multiple functional binding modes that affected channel processes differently. Our interpretation is consistent with the state model for reversible protein-protein associations, in which conformational rearrangements involving side chain movements optimize binding [32].
A number of electrophysiological studies, including single channel recordings and simulations using Eyring reaction rate theory models has supported the assumption that the b1 subunit modulates VGSC currents by shifting the fraction of channels in fast and slow gating modes [6,33,34]. This implies that a subunits alone are able to inactivate and recover in a fast gating mode, but the probability is low because the process is thermodynamically unfavourable [34]. Consequently, the binding of b1 to the a   The time constants of inactivation were obtained by fitting with single exponential equations, and the time constants of recovery were obtained from fitting with double exponential equations. Each cell was fitted individually and the mean time constants from each group were compared with the Nav1.4+b1 wild type, statistically significant differences were determined with a one-way ANOVA test followed by a Holm-Sidak test for multiple comparisons (n = 4-9, *P,0.05, **P,0.01, ***P,0.001). doi:10.1371/journal.pone.0081995.t002 subunit may reduce the energy barrier and increase the probability of fast mode gating. Our experimental findings and computed predictions are consistent with these assumptions, they also suggest that the mechanisms by which b1 accelerates the kinetics of sodium channels and modulates voltage dependency are independent. Here, the electrophysiological results revealed that R89A decreased the b1-induced acceleration of inactivation and recovery from inactivation without affecting the voltage dependence or the availability after repetitive depolarisations. Conversely, Y32A and G33M, point mutants in the first loop of the b1 extracellular domain, selectively affected the voltage dependence of inactivation without affecting the kinetics. Qu et al. discovered similar site-specific effects on b1-induced modulation of sodium channel a subunits [9]. Their results with chimeric Nav1.2/ Nav1.5 channels showed that the substitution of a single extracellular loop in domain IV decreased the ability of b1 to modulate the kinetics of inactivation and recovery without major effects on the voltage dependence [9].
The extracellular domain of the VGSC b1subunit and the extracellular vestibule of the pore-forming a subunit were predicted to be highly acidic. Based on in silico studies, we hypothesized that R89 acted within its side chain to form a solvent exposed salt bridge, which could be broken as it approached the external vestibule of the a subunit. After depolarization, R89 may directly bind the a subunit. The simulated side-flipping contribution from the intracatenar to intercatenar protein-protein association supports the findings indicating that a shift in the fraction of VGSCs that undergo rapid inactivation and recovery could occur without modifying the midpoint voltage of inactivation.
In contrast, the C43A mutant produced general loss-of-function effects, possibly due to an increase in translational entropy caused by the loss of the theoretical disulphide bond with C21. This may allow the N-terminal loop in which this cysteine is located more flexibility (Figure 5 B). Based on the comparisons with crystal structures, we argue that this structural change could not catastrophically affect folding. Indeed, none of the Ig-like proteins form disulphide bonds at this site although they conserved antiparallel b-sheet motifs (Figure 1 B). However, we cannot rule out the possibility that the cytoskeletal trafficking or the membrane insertion of the mutants was diminished, therefore current amplitudes between groups cannot be compared. Other studies have obtained similar results on the VGSC b3 subunit. The C24A point mutant in b3, which corresponded to C43A in b1, shifted the voltage of activation and modified the half voltage of inactivation of the Nav1.5 channel [35]. In another series of studies, C24A prevented homophilic binding as shown by immunoprecipitation assays in HEK293 cells. The interpretation was that the disulphide bond formed by C24 does not directly participate in the binding site for dimerization but that it contributes to orienting the interacting residues [29].
The reductions in the fractions of mutant channels recovered after the first millisecond were 95% for C43A and 90% for R89A, and the changes remained significant after a Dt of 10 ms for R89A and after a Dt of 100 ms for C43A (figure 1 A). These effects seem consequential according to the literature. Importantly, the evidence for the principles ruling b1 association to the a subunit have not advanced in over 14 years [5,10]. These studies indicate that an acid triad domain (E23, D25, E27) in the extracellular domain of b1 located between C43 and Y32 plays a key role. Expression of the b1 triple mutant E23Q D25Q E27Q with the Nav 1.2 channel reduced the fraction of channels that underwent rapid inactivation by 22% compared to WT Nav1.2+b1 [10]. Point mutations and larger sequence manipulations caused loss-of-function effects with comparable magnitudes to changes elicited by C43A and R89A [10].
The biophysical response to the R89A substitution seemed relevant compared to R85C or R85H, two mutations on the b1 subunit associated with generalized epilepsy with febrile seizures plus (GEFS+) [36]. Xu et al. (2007) found that R85C but not R85H affected current inactivation of the Nav1.2 channels [36]. After b1 WT model refinement, we found that R85 had the potential to form electrostatic interactions with nearby acidic residues, similar to R89. The reported discrepancy was computationally simulated. The results indicated that substituting histidine for arginine could rearrange the non-covalent bond network and that cysteine replacement disrupts it.
In addition, our model predicts that E87 forms a salt bridge with R89 and thus may be of functional relevance. Indeed the E87Q mutation is associated with Brugada syndrome, a channelopathy identified in a Turkish family [37]. The C121W mutation was also associated with GEFS+, and biochemical studies have indicated that C121 forms an extracellular intramolecular disulphide bridge that was disrupted by this mutation [38]. Our b1 model suggested that C40 would be its most probable binding partner ( Figure 1A).
The deceleration of recovery and the increase in rundown on high frequency stimulation may be important factors that restrain hyper-excitability in pain conduction and conditions such as epilepsy [11]. Interestingly, R89A delayed the recovery from inactivation without increasing the rundown, and C43A disrupted both b1-induced effects.
Our results were consistent with the Monod-Wyman-Changeux model in that: (I) mutations can shift the spontaneous allosteric equilibrium between relaxed and tense conformational states of the complex (in this case fast versus slow gating and voltage sensitive to voltage hypersensitive states), and (II) in a system of interacting protein subunits (or protomers) within a membrane lattice, various classes of responses to specific regulatory signals can exist, from graded to all-or-none phase transitions, depending on the isomerisation and the free energy of the interaction between protomers in large and periodic protein assemblies [39].
In this study, we suggest the existence of 2 or more discrete a-b1 binding states, one that controls rapid gating and another that dictates voltage sensing modulation. Taken together, these data may enhance the understanding of pathologies associated with the interactions of a and b1 in VGSCs and enable the design of molecular compounds that selectively interfere with a particular electrophysiological process. Figure S1 Multiple sequence alignment of b1 sequences. Vertebrate b1 sequences in gray (Danio rerio, Sternopygus macrurus, Osmerus mordax, Takifugu rubripes, Homo sapiens, Macaca mulatta, Canis lupus familiaris, Mus musculus, Oryctolagus cuniculus) aligned with 3 related proteins including the rat b2 subunit, the frog b3 subunit, and the rat myelin protein P0 (lines 10-12 in brown). Highly conserved residues were highlighted in cyan. The acid triad domain (ExDxD) was implicated in the modulation of the VGSC brain isoform (first panel in red). Substitutionss of either Cys or His for Arg85 (second panel in red) were related to familial cases of GEFS+. Regions highlighted in yellow and green form b strands and a helixes, respectively. Consensus secondary structure predictions for the b sequence (performed by Jpred) were presented on line 14. (TIF) Figure S2 Effects of C43A and R89A b1 mutants on the electrostatic surface potential (ESP) of the extracellular domain. ESP represented by Coulombic colouring within a radius of 5 Å from residue 32 (in the left panel). The substitution of tyrosine (Y32A) increased the polarity by destabilizing a networked coordination of intermolecular H-bonds. This was caused by the loss of the hydroxyl group, which shared a hydrogen with Asp148. The hydrophobicity was reduced by the loss the aromatic phenyl ring. The ESP within a radius of 5 Å from position 33 (in the right panel) showed the introduction of a long hydrophobic side chain in the G33M mutant. Calculations were made under the AMBER-Gasteiger force field in Chimera 1.3.5 after energy minimization under AMMP in VEGA ZZ 2.3.2. (TIF)