Na_v1.4-CHO cells and Na_v1.7-CHO cells were maintained in Iscove’s modified Dulbecco’s medium (Gibco) supplemented with 10% (v/v) fetal bovine serum (FBS) and 200 μg/mL G418 (Sigma). All cells were cultured at 37 °C and 5% CO₂ atmosphere. Mouse myoblasts C2C12, obtained from the Stem Cell Bank of The Chinese Academy of Sciences, were grown in Dulbecco’s modified Eagle‘s medium (DMEM) enriched with 10% heat-inactivated FBS, 2 mM L glutamine, 100 U/mL of penicillin, and 100 μg/mL of streptomycin. The C2C12 cells typically reached confluence in 3–5 days in 75-T flasks. For passaging, monolayers were treated with 0.25% trypsin for 3–5 min, rinsed twice with DMEM, and sub-cultured at a ratio of 1:4. Similarly, the Chinese hamster ovary (CHO) cells were cultured in IMDM medium supplemented with 10% fetal bovine serum, 1% penicillin/streptomycin, 100 μM hypoxanthine, and 16 μM thymidine. When 90% confluence was achieved, cells were separated by 0.25% trypsin EDTA, and patch clamp electrophysiology experiments were carried out.

The ProTx2 gene and its mutants (101 bp) were cloned into the pET-32a vector, incorporating TrxA, 6xHis, and enterokinase cleavage sites. Agarose gel electrophoresis confirmed the presence of the expected 100-bp bands, indicating successful plasmid construction (Supplementary Figure S1). Sequencing confirmed 100% amino acid fidelity with the original sequence, affirming the successful creation of the active peptides (Supplementary Figure S2). IPTG (isopropyl β-D-1-thiogalactopyranoside) induction yielded prominent expression bands for both ProTx2 and its mutants, aligning with the expected molecular weight of approximately 21 kDa (Figure 2). Subsequent analysis revealed over 90% expression of the target protein in the supernatant post-cell lysis, confirming effective soluble expression in the pET-32a system. Purification via nickel chelate affinity chromatography (buffer A: Tris-HCl, pH 8.0; buffer B: Tris-HCl with 50 mM imidazole, pH 8.0; buffer C: Tris-HCl with 150 mM imidazole, pH 8.0; and buffer D: 0.05M EDTA) and subsequent SDS-PAGE analysis identified a clear band corresponding to the ProTx2 fusion protein at 21 kDa, indicative of successful elution (Figures 2C,D). Further purification was achieved by Q Sepharose Fast Flow ion-exchange chromatography, employing a gradient elution strategy (gradient elution, the high-concentration eluent was Tris-HCl with 150 mM imidazole and 1M NaCl, pH 8.0; the low-concentration eluent was Tris-HCl with 150 mM imidazole, no salt) to isolate the target protein effectively. Enterokinase was finally utilized to cleave the fusion of His tags. Protein purity was assessed using Tricine–SDS-PAGE and quantified by HPLC analysis with area normalization (Supplementary Figures S3A, S3B).

Male KM mice (18–22 g) were obtained from Liaoning Changsheng Biotechnology Co., Ltd. The mice were acclimated for 1 week at 18°C–22°C and 50%–70% humidity, following approval by the Animal Ethics Committee of Shenyang Pharmaceutical University, China. The pain response was assessed using a formalin-induced paw-licking test. The saline solution, 5 mg/kg morphine, commonly used analgesic medications in clinical practice, or 0.75 mg/kg peptides, ProTx II, and its mutants, were injected via the tail vein to relief the pain. After half an hour, a 20-μL 5% (v/v) formalin solution was administered as an activator to induce pain in mice (n = 10) (Miranda et al., 2007). The cumulative pain behavior that occurred 0–5 min (phase I, an immediate nociceptive phase) and 20–60 min (phase II, a later inflammatory phase) after injection was recorded for 60 min at 5 min intervals in individual cages to measure the total licking time in phases I and II, respectively.

Additionally, a forced swimming test measured endurance following 7 days of treatment with ProTx2 or its mutants (Porsolt et al., 1977; Galeotti et al., 2002; Xu et al., 2017). Fifty male KM mice were randomly divided into five groups with 10 mice in each group. The groups were classified as normal saline control, ProTx2, R13A, R13D, and K27Y.

C2C12 cell viability in response to ProTx2 or it mutants was evaluated using the CCK8 assay (Pan et al., 2019; Ge et al., 2020). A set of concentration gradients ranged from 0, 20, 40, 60, 80, 100, 120, 140, and 160 μg/mL, using a serum-free medium as a diluent. Cells were seeded in 96-well plates with 5*10³/well in five duplicates. After 18 h of treatment, cell viability was assessed by measuring the absorbance at 450 nm, following the addition of the CCK-8 reagent.

The simulation utilized human-derived PDB structures 6AGF (hNa_v1.4) and 6N4I (hNa_v1.7), with lipid bilayers modeled using POPC in Charmm-GUI (https://charmm-gui.org/) (Chen and Chung, 2012). Simulations were performed using Amber 18, employing the Amber-ff14SB force field (Maier et al., 2015), with the particle mesh Ewald (PME) algorithm for electrostatics (Pederson and McDaniel, 2022) and the SHAKE algorithm (Kräutler et al., 2001) for hydrogen movement. Langevin software maintained the temperature, with a 2-fs time step over a 100-ns simulation. Energy decomposition was conducted with MM/GBSA to analyze binding energies and key residue contributions (Izaguirre et al., 2001; Walton and Vanvliet, 2006).

After 100-ns MD, the root mean square deviation (RMSD) and cluster analysis were performed to determine the quality of simulation by the CPPTRAJ (Roe and Cheatham III, 2013). RMSD curves showed the volatility of atoms in the athletic systems whose stationary band represented a stable simulation trajectory. Cluster analysis was used to identify the similarity of poses in the simulation from the stable trajectory. hieragglo epsilon was set as 1.5 Å, and the backbone atoms were calculated. For the analysis, the structures would be steady if the largest cluster occupied a proportion greater than 50% in all clusters.

Energy decomposition analysis (EDA) was carried out to analyze the binding energy between ProTx2 and Na_v1.7/Na_v1.4 in order to detect the ability of the combination of different ligands and receptors and show the contribution of every residue in ligands. The molecular mechanics generalized Born surface area (MM/GBSA) method (Hou et al., 2011) was performed on the basis of the snapshots extracted from the stable period of the trajectories. The total binding energy would be calculated, which was the overall evaluation for a combination including the sum of average molecular mechanical gas phase energies (EMM) and salvation-free energies (Gsolv). In more detail, EMM was divided into the summation of van der Waals (Evdw) and electrostatics (Eele) interaction energies, and Gsolv was described as a polar and non-polar solvation interaction due to an implicit solvent model and being solvent-accessible. Meanwhile, the decomposition of binding free energy of each ligand–residue pair emerged from which the key residues with high energy contribution (lower than −1 kcal/mol) were selected, and their dominant interactions were detected.

Dynamic cross-correlation analysis (Hünenberger et al., 1995) (DCCM) was performed to the last stable trajectories of MD to speculate the relative region in the dynamic combination between ligands and receptors using CPPTRAJ of Amber20. The coordinate matrix was calculated based on the coordinates of Cα on the complexes which were set as the input data. The relevancy in the dynamic simulation will be analyzed between the residues from the construction of the former matrix results.

In addition, the salt-bridge, one kind of a strong electrostatic force, was calculated using VMD 1.9.2 (Humphrey et al., 1996), which was regarded as existent when the distance was less than 4 Å with a proportion higher than 50%. The hydrogen bond was analyzed using CPPTRAJ in Amber20, whose criterion was also 50% in all snapshots. Other parameters in the methods were kept as default.

Automatic electrophysiological platforms (NPC-1 and Port-a-Patch systems, Nanning Technologies GmbH, Munich, Germany) and EPC10 amplifiers (HEKA Elektronik, Lambrecht, Germany) assessed the effects of recombinant ProTx2 and its mutants on CHO-Na_v1.4 and CHO-Na_v1.7. PatchControl HT (Nanion Technologies GmbH, Germany) and PatchMaster software applications (HEKA Elektronik, Lambrecht, Germany) were used to execute the experimental process.

Cells underwent separation using 0.25% trypsin, followed by centrifugation at 800 rpm for 5 min, and were subsequently resuspended in extracellular fluid. CHO cells displaying characteristic morphology and vitality under microscopic examination were selected, and 10 cells were grouped per set. Microelectrodes were fashioned from blank glass tubes through a two-step pulling process, resulting in tip diameters of 0.8–1 μm and resistance levels between 1 and 3 MΩ. The composition of the extracellular solution was as follows (in mol/L): 140 NaCl, 4 KCl, 5 D-glucose monohydrate, 1 MgCl₂, 2 CaCl₂, and 10 HEPES (N-2-hydroxyethylpiperazine-N-2-ethane sulfonic acid)/NaOH, adjusted to a pH of 7.4 and maintained at 4°C. The intracellular electrode solution contained the following (in mol/L): 50 CsCl, 10 NaCl, 60 CsF, 20 EGTA (Xu et al., 2021), and 10 HEPES/CsOH, with a pH of 7.2. Both solutions were filtered through a 0.22-μm filter and stored at −20°C (Xiangxue et al., 2015; Xu et al., 2020).

In our whole-cell recordings, potassium currents were inhibited using Cs and a relatively high Mg²⁺ concentration (2 mmol/L). Cells exceeding 35 μm in diameter were exclusively examined, with a holding potential set at −80 mV. Sodium currents were elicited with a test pulse from −80 to +80 mV, in increments of 10 mV over 50-ms intervals (Xu et al., 2017). The transient sodium current amplitudes in CHO cells both before and subsequent to the addition of ProTx2 and R13D were documented. All experimental procedures were conducted at an ambient temperature of 22°C.

Statistical analyses were conducted using SPSS 22.0, employing one-way ANOVA and Tukey’s post hoc tests for group comparisons (Shenoy et al., 2021). Results were considered significant at p < 0.05. GraphPad Prism 5, Pymol (Schiffrin et al., 2020) and Gunplot (Latendresse et al., 2022) generated kinetic plots from molecular dynamic simulation data.

ProTx2 gene and its mutants (101 bp) were successfully cloned into the pET-32a vector, incorporating TrxA, 6xHis, and enterokinase cleavage sites (Figure 2).

Agarose gel electrophoresis confirmed the presence of the expected 100-bp bands, indicating successful plasmid construction (Figure 2) (Supplementary Figure S2).

IPTG induction yielded distinct expression bands of approximately 21 kDa, matching the anticipated sizes of the fusion proteins (Figure 2).

Subsequent purification via nickel chelate affinity chromatography and Q Sepharose Fast Flow ion-exchange chromatography resulted in a prominent band corresponding to the ProTx2 fusion protein at 21 kDa. Enterokinase cleavage and ultrafiltration provided a purified ProTx2 sample, with a final purity of approximately 96%, as determined by HPLC (Figures 2C,D) (Supplementary Figure S3B).

Three potential hNa_v1.7 subtype-specific ProTx2 mutants (R13D, R13A, and K27Y) were identified through a virtual screen. The formalin-induced paw-licking assay in mice was established, and the cumulative pain behavior was recorded for 60 min at 5-min intervals (Figure 3A). The total licking time in phase I (0–5 min) and phase II (20–60 min) was measured, respectively, (Figure 3B). In phase I, the durations of total licking time of the control, morphine, wild ProTx2, R13D, R13A, and K27Y groups were 112.3 ± 3.9 s, 7.5 ± 1.6 s, 53.3 ± 1.6 s, 33.5 ± 1.6 s, 45.5 ± 1.1 s, and 48.7 ± 2.4 s, respectively, while in phase II, they were 502.6 ± 9.3 s, 42.1 ± 6.1 s, 203.5 ± 5.3 s, 128.0 ± 5.2 s, 173.5 ± 2.8 s, and 199.1 ± 2.9 s, respectively. These results revealed that the R13D mutant significantly alleviated pain in both the acute and inflammatory phases, outperforming the wild-type ProTx2. Specifically, R13D exhibited a marked reduction in licking time during both pain phases, suggesting an improved analgesic profile.

Following a 7-day treatment regimen, mice received injections of ProTx2 or its mutants and underwent forced swimming tests 20 min after the final treatment. Notably, R13D-treated mice exhibited significantly prolonged endurance compared to the wild-type ProTx2 group, suggesting reduced muscular toxicity (p < 0.05). In contrast, the performance of mice treated with R13A and K27Y did not significantly differ from those treated with wild-type ProTx2 (Figure 4).

Additionally, to elucidate the interaction between the peptide’s toxicity and Nav1.4, the viability of C2C12 myoblasts, which express Nav1.4 abundantly, after 18-h exposure to varying concentrations of ProTx2 and its mutants were assessed. Remarkably, R13D demonstrated significantly improved cellular viability in a concentration-dependent manner, particularly between 100 and 160 μg/mL, underscoring its reduced toxicity relative to the wild-type peptide. No significant differences in cellular viability were observed with R13A and K27Y, aligning with the forced swimming test outcomes.

Molecular dynamic (MD) simulations over 100 ns provided detailed insights into the interactions between ProTx2 and its mutants (R13D, R13A, and K27Y) with hNav1.7 and hNav1.4 (Figure 5). Root mean square deviation (RMSD) analyses indicated that with the exception of K27Y, both wild-type ProTx2 and the other mutants bound similarly to hNav1.7, suggesting minimal alterations in their binding modes. However, notable conformational changes were observed when these peptides interacted with hNav1.4, particularly for R13D.

Dynamic cross-correlation matrix (DCCM) analyses of stable simulation trajectories (90–100 ns) further highlighted the specific interactions driving these effects (Piao et al., 2019). Enhanced interactions between the peptides’ C-terminal (residues 160–170) and hNav1.7’s extracellular loop (DII S3–S4, residues 99–108) were evident post-mutation, as marked by intensified correlation signals (red box in Supplementary Figure S4). Similarly, changes were noted at the peptides’ N-terminal and the extracellular loop connecting S1–S2 at DII (residues 39–48) on hNav1.7 (brown box in Supplementary Figure S4).

In contrast, the interaction dynamics between R13D and hNav1.4’s extracellular loops showed significantly strengthened negative correlations, suggesting a distinct modulation of binding post-mutation, which is not observed with R13A or K27Y. These findings align with the structural data from a published chimera model (Xu et al., 2019), confirming that specific residue interactions, rather than broad conformational shifts, underpin the altered specificity of ProTx2 mutants toward hNav1.7 over hNav1.4. Crucially, R13D’s enhanced analgesic efficacy, paired with reduced muscle toxicity, underscores the potential of targeting specific residues to improve therapeutic profiles. This nuanced understanding of peptide–channel interactions points a new direction for refining the specificity and safety of Nav1.7-targeted analgesics.

Free energy decomposition analysis (EDA) revealed that the significant decrease in total binding free energy (ΔEtotal) for R13D when bound to hNav1.7 highlighted a notable improvement in affinity primarily due to enhanced van der Waals and electrostatic interactions. Intriguingly, this enhancement stemmed from a cascade of effects (Table 1) initiated by the single-point mutation at R13D, which did not directly interact with hNav1.7 but significantly amplified the peptide’s affinity.

The effect of the single-point mutation was first propagated through the intramolecular hydrogen bond network as D13-D10-S11-K14-Q3-C16 on mutant R13D (Supplementary Figure S6), which was called point-to-surface. This network formed an extremely stable intramolecular 3D architecture, resulting in the formation of a strong steric hindrance between C16 and K27. Predictably, K27 was the point “with wide-ranging implications.”

Further analysis pinpointed K27’s pivotal role in this affinity enhancement. Adjustments in K27’s spatial positioning dramatically increased K26’s proximity to ASP96 on hNav1.7, forming a new, potent salt bridge absent in the wild-type ProTx2. This strategic mutation at site 13 thus effectively concentrated the binding energy on K26, manifesting as a surface-to-point amplification in affinity. To decipher the role of K27 in the combination of R13D and hNa_v1.7, ΔEtotal was decomposed into each residue on the receptor and ligand in order to identify the individuals with the greatest contributions in the combination. The similar fluctuation regions (Figure 6) indicated that the binding modes of ProTx2 and its mutants to hNa_v1.7 were almost identical, which was consistent with the analysis of DCCM and 3D models. However, it was evident that the contribution of M6 in the N-terminal and R22/W24/K26 in the C-terminal on the peptide was significantly distinct when they were located on different ligands (Table 2). Among them, the affinity increase in K26 in the mutant R13D was the most obvious as five times in contrast with wild-type ProTx2 versus two times with other residues. With the further EDA of the interaction type of the critical residues on the ligand in each system (Figures 7A, 8A), it was found that the major reason was that the spatial position adjustment of K27 affected K26 to be closer to ASP96 on hNa_v1.7 and formed a new powerful salt bridge between them, which did not exist in the wild-type ProTx2 system. So far, the energy caused by the mutation at site 13 was concentrated on K26, which was called surface-to-point. Another notable difference was the newly developed interactions between W24 and two highly conserved negatively charged residues ARG104/ARG107, which were identified as the direct factors related to the activation of Na_vs (Hampl et al., 2016). In addition, the presence of the salt bridge formed between R22 and GLU98 brought out the improvement in the contribution of R22.

While R13D emerged as a formidable variant in enhancing hNav1.7 affinity, R13A also displayed a notable increase in energy contribution from R22, leading to a stronger hNav1.7 binding compared to wild-type ProTx2 (Table 2). This was attributed to a newly formed salt bridge between R22 and GLU98 (Figures 7B, 8B). The collective data underscore R13D’s standout role in augmenting ProTx2’s affinity for hNav1.7, attributed to a synergistic interplay among key residues rather than mere conformational changes, setting a new paradigm in the specificity enhancement of peptide-based therapeutics.

Obviously, for the binding affinity between the peptide and hNa_v1.7, the allied heroes in the mutant R13D were more intimidating than the lone hero in mutants R13A and K27Y.

Contrary to its interaction with hNav1.7, the total binding free energy (ΔEtotal) of ProTx2 mutants remained similar to the wild-type when bound to hNav1.4 (Table 1). However, our analysis revealed significant shifts in their binding dynamics (Figure 9). The N-terminal region, initially the primary binding site in wild-type ProTx2 (Figures 10A, 11A), underwent notable changes in its interaction pattern with hNav1.4 post-mutation, especially for R13D, where a pronounced shift toward C-terminal interactions was observed (Table 2; Figure 11C).

The alteration in R13D’s electrical properties at residue 13 eliminated its electrostatic interaction with TYR42 on hNav1.4, transforming the N-terminal’s attractive forces into repulsive forces (Figure 11C). This shift led to an increased affinity of the peptide’s C-terminal for the receptor, a pattern also seen in R13A, albeit with the N-terminal still engaging due to the neutrality of glycine at position 13 (Figures 10B, 11B). Remarkably, the mutation at residue 27 significantly influenced the binding dynamics. For instance, substituting lysine with tyrosine at this position increased the hydrophobicity of the peptide’s C-terminal, drawing K26 closer to ASP49 on hNav1.4 and forming a critical salt bridge. This interaction underscores the direct impact of mutations at site 27 on the peptide’s binding affinity (Figure 11D).

Overall, our findings highlight the N-terminal as a pivotal domain in the ProTx2’s interaction with hNav1.4, with mutations at site 13 inducing indirect effects that modulate the peptide receptor affinity. The alterations in the peptide’s binding modes, particularly the enhanced role of the C-terminal in the mutant R13D, mark a significant departure from the wild-type interactions, emphasizing the nuanced influence of specific residues on the binding efficacy.

Energy decomposition analysis (EDA) (Figures 6A, 9A) highlighted the critical role of the extracellular loops connecting S1–S2 and S3–S4 on the voltage-sensor domain 2 (VSD2) in Nav channels, identified as Loop1–2 and Loop3–4, in mediating ProTx2’s interaction with hNav1.4 and hNav1.7. Despite an 85.71% sequence similarity in the VSD2 region of these subtypes, Loop3–4’s composition diverges significantly (blue box in Supplementary Figure S7) with hNav1.7 featuring acidic residues that form strong salt bridges with R22 and K26 of ProTx2, an interaction absent in the neutral Loop3–4 of hNav1.4. This distinction enhances the binding affinity to hNav1.7 while constraining the ligand’s interaction modes, particularly weakening its engagement with Loop1–2 in hNav1.7 compared to hNav1.4. Notably, the R13D mutant exploits these structural disparities, disrupting its bond with Loop1–2 in hNav1.4 while bolstering connections with both loops in hNav1.7, as deduced from DCCM analysis (Supplementary Figure S5C; Supplementary Figure S6C).

Electrophysiological assays quantified the effects of ProTx2 and its R13D mutant on sodium currents within CHO cells expressing hNav1.4 and hNav1.7 (Figure 12). Although both ProTx2 and R13D inhibited currents across these subtypes, R13D exhibited a marginally diminished effect on hNav1.4, yet pronouncedly intensified inhibition of hNav1.7, aligning with predictions from molecular dynamics simulations. Specifically, at 1 μM concentration, R13D notably decreased the peak sodium current amplitude in hNav1.7 from 82% to 68% (p < 0.01) and modestly increased it in hNav1.4 from 50% to 59%. This shift underscores R13D’s enhanced analgesic potential and reduced muscle toxicity, corroborating the insights gained from computational and experimental analyses.