Computational Design of High-Affinity Blockers for Sodium Channel NaV1.2 from μ-Conotoxin KIIIA

The voltage-gated sodium channel subtype 1.2 (NaV1.2) is instrumental in the initiation of action potentials in the nervous system, making it a natural drug target for neurological diseases. Therefore, there is much pharmacological interest in finding blockers of NaV1.2 and improving their affinity and selectivity properties. An extensive family of peptide toxins from cone snails (conotoxins) block NaV channels, thus they provide natural templates for the design of drugs targeting NaV channels. Unfortunately, progress was hampered due to the absence of any NaV structures. The recent determination of cryo-EM structures for NaV channels has finally broken this impasse. Here, we use the NaV1.2 structure in complex with μ-conotoxin KIIIA (KIIIA) in computational studies with the aim of improving KIIIA’s affinity and blocking capacity for NaV1.2. Only three KIIIA amino acid residues are available for mutation (S5, S6, and S13). After performing molecular modeling and simulations on NaV1.2–KIIIA complex, we have identified the S5R, S6D, and S13K mutations as the most promising for additional contacts. We estimate these contacts to boost the affinity of KIIIA for NaV1.2 from nanomole to picomole domain. Moreover, the KIIIA[S5R, S6D, S13K] analogue makes contacts with all four channel domains, thus enabling the complete blocking of the channel (KIIIA partially blocks as it has contacts with three domains). The proposed KIIIA analogue, once confirmed experimentally, may lead to novel anti-epileptic drugs.


Introduction
The voltage-gated sodium (Na V ) channels are responsible for the initiation and propagation of action potentials, and thus play a key role in membrane excitability and electrical signaling [1][2][3]. The nine subtypes of Na V channel α-subunits (named Na V 1.1-1.9) are differentially encoded by SCNxA (x = 1-5 corresponding to Na V 1.1-Na V 1.5, and x = 8-11 corresponding to Na V 1.6-Na V 1.9, respectively) and expressed throughout tissues [4,5]. Na V 1.2 was first purified from the rat brain, hence it is also called brain type-II Na V channel. It functions mainly in the central nervous system [6,7]. Genetic mutations in Na V 1.2 have been linked to seizures, infantile spasm, and autism spectrum disorders. Na V 1.2, therefore, is studied as a potential target for the development of therapeutics for epilepsy [8].
A class of conotoxins from cone snails (called µ-conotoxins) have been found to target the extracellular vestibule of the Na V channel pore [9] and thus offer promising scaffolds for the rational design of novel peptide-based therapeutic drugs [10,11]. µ-conotoxin KIIIA (KIIIA) contains only 16 amino acids, but retains the typical µ-conotoxins disulfide bond pattern [12], and has shown potential inhibition of TTX-sensitive Na V channels [9,13,14]. KIIIA binds to different Na V channel subtypes with variable degrees of affinity and block, e.g., with 5 nM affinity for rat Na V 1.2, 37 nM for rat Na V 1.4, and 97 nM for human Na V 1.7 [9,14,15]. Experimental structure-activity analysis indicates that the KIIIA residues K7, R10, H12, and R14 are important for binding to various Na V channel subtypes [9,14]. The relative contribution of KIIIA residues in binding to Na V channels differ between the Na V channel subtypes. For example, the neutral replacement of R14 by A14 in KIIIA group of R10 is coordinated by D1426. In addition, W8 engages in an interaction with Y362 on the extracellular loop of domain I, and the guanidinium group of R14 is H-bonded with the carbonyl group of L920 [17]. The lack of any interactions with any domain IV residues provides an explanation for the partial block of Na V 1.2 by KIIIA.
Once a toxin is selected as a potential drug lead for a target ion channel, it is important to improve its affinity and selectivity for the target to reduce the dosage and avoid side effects. For toxin peptides, this can be achieved by designing analogues of the toxin through mutations of selected residues, which are predicted to make further contacts with the channel residues. In the case of KIIIA, the residues S5, S6, and S13 are not involved in any interactions with the Na V 1.2 residues, hence they provide potential mutation sites to improve its affinity for Na V 1.2.
In this paper, we use computational methods to analyze the molecular mechanism of KIIIA binding to Na V 1.2. Starting with the cryo-EM structure of Na V 1.2-KIIIA, we first performed molecular dynamics (MD) simulations to check the robustness of the observed interaction network. Most of the key contacts were found to be preserved during the MD simulations. Then we used the VMD plugin Mutator [18] to generate structures of the KIIIA analogues with the mutations S5R, S6D, and S13K individually and in various combinations. We used the MD simulations to interpret the binding modes of these KIIIA analogues in complex with Na V 1.2. Compared to the wildtype peptide, we have identified several unique contacts between the KIIIA analogues and Na V 1.2, which indicate that the mutations could considerably enhance their affinity for Na V 1.2. Moreover, after the S5R mutation, R5 makes contact with a domain IV residue, resulting in a complete block of the pore. These results suggest that the proposed KIIIA analogues could provide valuable leads in the development of therapeutics for epilepsy.

Critical Interaction Networks Identified from the Na V 1.2-KIIIA Complex
We start with a comparison of the Na V 1.2-KIIIA binding mode that was obtained from the MD simulations with that of the cryo-EM structure ( Table 1). The average distances between the nearest charged N and O atoms that were computed from the MD simulations were compared to the corresponding experimental values.
Most of the N-O distances that were obtained from MD were in good agreement with the experimental values that are in Table 1, especially when the pair of atoms are at contact distances of~3 Å when the interaction is strongest. These include the key pairwise interactions, E945-K7, D1426-R10, D949-H12, and L920-R14. A notable difference occurs for the interaction of R14 with E919 and E1444. In the cryo-EM structure, the R14 side chain is in between the side chains of E919 and E1444, and shows a slight preference for the former, whereas in the MD results, the R14 side chain has a distinct preference for the E1444 side chain. Most of the N-O distances that were obtained from MD were in good agreement with the experimental values that are in Table 1, especially when the pair of atoms are at contact distances of~3 Å when the interaction is strongest. These include the key pairwise interactions, E945-K7, D1426-R10, D949-H12, and L920-R14. A notable difference occurs for the interaction of R14 with E919 and E1444. In the cryo-EM structure, the R14 side chain is in between the side chains of E919 and E1444, and shows a slight preference for the former, whereas in the MD results, the R14 side chain has a distinct preference for the E1444 side chain.
Mar. Drugs 2022, 20, 154 4 of 11 R14 is located within the C-terminal tail region, with the side chain extending between the extracellular P1 helix of domain II and the extracellular P2 helix of domain III. The guanidinium group of R14 is wedged between the side chains of E919 and E1444, and it also makes an H-bond with the carbonyl group of L920 ( Figure 1). An inspection of the time series of the R14 interactions with L920, E919, and E1444 (Supplementary Materials Figure S1) shows that initially R14 interacts with E919 (which is the bias of the cryo-EM structure), but after a short while it switches to E1444 and forms a stable contact with it. The R14A mutation reduces the KIIIA affinity by two orders of magnitude, which can be explained by the stronger coupling of R14 to the channel residues that were observed in the MD simulations. A detailed view of R14 interacting with the channel residues L920 and E1444.

MD Simulations of the KIIIA[S5R] Analogue
R5 is positioned near the N-terminal, and its positively charged side chain is in between the channel residues D1692 and V1689 in domain IV ( Figure 2A). The guanidinium group of R5 forms a salt bridge with D1692 and makes an H-bond with the main chain carbonyl of V1689. A time series of the pair distances show that R5 makes strong and stable contacts with D1692 and V1689 ( Figure 2B), which is expected to make a significant contribution to the binding free energy of KIIIA[S5R]. In addition, the S5R mutation introduces a positively charged side chain near the K7

MD Simulations of the KIIIA[S5R] Analogue
R5 is positioned near the N-terminal, and its positively charged side chain is in between the channel residues D1692 and V1689 in domain IV ( Figure 2A). The guanidinium group of R5 forms a salt bridge with D1692 and makes an H-bond with the main chain carbonyl of V1689. A time series of the pair distances show that R5 makes strong and stable contacts with D1692 and V1689 ( Figure 2B), which is expected to make a significant contribution to the binding free energy of KIIIA[S5R].
In addition, the S5R mutation introduces a positively charged side chain near the K7 and R10 residues, which perturbs their positions in the Na V 1.2-KIIIA structure, while other KIIIA residues (W8, D11, H12, and R14) have maintained their contacts after the mutation. The K7 side chain interacts with E945 in the EEDD ring in Na V 1.2-KIIIA, but after the S5R mutation, it slightly shifts downwards and starts interacting with E942 in the DEKA ring as well, thus gaining extra binding affinity ( Figure 3).
Similarly, the R10 side chain, which interacts with D1426 in the EEDD ring in Na V 1.2-KIIIA, shifts towards E942 in the DEKA ring ( Figure 4). The R10 side chain is seen to make a more favorable contact with E942 in the Na V 1.2-KIIIA[S5R] complex compared to its contact with D1426 in the Na V 1.2-KIIIA complex. Thus, we expect the S5R mutation to indirectly enhance the contribution of the R10 interaction to the affinity.
R5 is positioned near the N-terminal, and its positively charged side chain is in between the channel residues D1692 and V1689 in domain IV (Figure 2A). The guanidinium group of R5 forms a salt bridge with D1692 and makes an H-bond with the main chain carbonyl of V1689. A time series of the pair distances show that R5 makes strong and stable contacts with D1692 and V1689 ( Figure 2B), which is expected to make a significant contribution to the binding free energy of KIIIA[S5R]. In addition, the S5R mutation introduces a positively charged side chain near the K7 and R10 residues, which perturbs their positions in the NaV1.2-KIIIA structure, while other KIIIA residues (W8, D11, H12, and R14) have maintained their contacts after the mutation. The K7 side chain interacts with E945 in the EEDD ring in NaV1.2-KIIIA, but after the S5R mutation, it slightly shifts downwards and starts interacting with E942 in the DEKA ring as well, thus gaining extra binding affinity ( Figure 3). Similarly, the R10 side chain, which interacts with D1426 in the EEDD ring in NaV1.2-KIIIA, shifts towards E942 in the DEKA ring ( Figure 4). The R10 side chain is seen to make a more favorable contact with E942 in the NaV1.2-KIIIA[S5R] complex compared to its contact with D1426 in the NaV1.2-KIIIA complex. Thus, we expect the S5R mutation to indirectly enhance the contribution of the R10 interaction to the affinity.   Similarly, the R10 side chain, which interacts with D1426 in the EEDD ring in NaV1. KIIIA, shifts towards E942 in the DEKA ring ( Figure 4). The R10 side chain is seen to ma a more favorable contact with E942 in the NaV1.2-KIIIA[S5R] complex compared to contact with D1426 in the NaV1.2-KIIIA complex. Thus, we expect the S5R mutation indirectly enhance the contribution of the R10 interaction to the affinity.

MD Simulations of the KIIIA[S5R, S13K] Analogue
We next include the S13K mutation in the KIIIA[S5R] analogue. The side chain of K is fully extended and forms ionic bonds with the side chains of D917 and E919 ( Figu  5A). An inspection of the contact distances in the NaV1.2-KIIIA[S5R, S13K] complex show that K13 makes ionic bonds with D917 and E919 ( Figure 5B). Other pairwise interactio have been maintained as in the NaV1.2-KIIIA[S5R] complex, so there is no penalty to th

MD Simulations of the KIIIA[S5R, S13K] Analogue
We next include the S13K mutation in the KIIIA[S5R] analogue. The side chain of K13 is fully extended and forms ionic bonds with the side chains of D917 and E919 ( Figure 5A). An inspection of the contact distances in the Na V 1.2-KIIIA[S5R, S13K] complex shows that K13 makes ionic bonds with D917 and E919 ( Figure 5B). Other pairwise interactions have been maintained as in the Na V 1.2-KIIIA[S5R] complex, so there is no penalty to the binding affinity due to breaking of the existing contacts. Thus, we expect the S13K mutation to increase the affinity of the KIIIA analogue to Na V 1.2 further.

MD Simulations of the KIIIA[S5R, S6D, S13K] Analogue
The S6D mutation creates a polar interaction between D6 and N333, which further reinforces the toxin binding in domain I ( Figure 6A). The time series of the N-O distances for the D6-N333 interaction indicate a slightly longer contact distance and larger fluctuations compared to those of the R5 and K13 interactions ( Figure 6B). Thus, we don't expect the S6D mutation to contribute to the KIIIA affinity as much as the S5R and S13K mutations. However, for drugs targeting the central nervous system, the peptide absorption is affected by the blood-brain barrier. The S6D mutation introduces one negative charge, which partially neutralizes the surplus positive charges in the KIIIA analogue. After the three mutations, the toxin analogue gains a net charge of +e, and the total charge of KIIIA[S5R, S6D, S13K] becomes +5e. The binding modes of the NaV1.2-KIIIA[S5R, S6D, S13K] and NaV1.2-KIIIA complexes are compared in Table 2. It is seen that the pairwise interactions that exist in NaV1.2-KIIIA have mostly been maintained after the three mutations which add three more contacts. Therefore, the KIIIA[S5R, S6D, S13K] analogue is expected to yield a substantially higher binding affinity to NaV1.2.

MD Simulations of the KIIIA[S5R, S6D, S13K] Analogue
The S6D mutation creates a polar interaction between D6 and N333, which further reinforces the toxin binding in domain I ( Figure 6A). The time series of the N-O distances for the D6-N333 interaction indicate a slightly longer contact distance and larger fluctuations compared to those of the R5 and K13 interactions ( Figure 6B). Thus, we don't expect the S6D mutation to contribute to the KIIIA affinity as much as the S5R and S13K mutations. However, for drugs targeting the central nervous system, the peptide absorption is affected by the blood-brain barrier. The S6D mutation introduces one negative charge, which partially neutralizes the surplus positive charges in the KIIIA analogue. After the three mutations, the toxin analogue gains a net charge of +e, and the total charge of KIIIA[S5R, S6D, S13K] becomes +5e. The binding modes of the Na V 1.2-KIIIA[S5R, S6D, S13K] and Na V 1.2-KIIIA complexes are compared in Table 2. It is seen that the pairwise interactions that exist in Na V 1.2-KIIIA have mostly been maintained after the three mutations which add three more contacts. Therefore, the KIIIA[S5R, S6D, S13K] analogue is expected to yield a substantially higher binding affinity to Na V 1.2.

Blocking of Sodium Channel Na V 1.2
Wildtype KIIIA residues make contact with the extracellular segments in domain I to III of Na V 1.2, but they make no contacts with domain IV, which leaves a narrow gap on the side of domain IV for entrance of a sodium ion into the selectivity filter. The observation of entry of an Na + ion into the selectivity filter ( Figure 7A) is consistent with the incomplete block of Na V 1.2 with~95% blockage. After the S5R mutation, the exposed ion permeation pathway is blocked completely by the new contacts between R5 and domain IV residues ( Figure 7B), which prevents the permeation of Na + across the toxin analogue. three mutations, the toxin analogue gains a net charge of +e, and the total charge of KIIIA[S5R, S6D, S13K] becomes +5e. The binding modes of the NaV1.2-KIIIA[S5R, S6D, S13K] and NaV1.2-KIIIA complexes are compared in Table 2. It is seen that the pairwise interactions that exist in NaV1.2-KIIIA have mostly been maintained after the three mutations which add three more contacts. Therefore, the KIIIA[S5R, S6D, S13K] analogue is expected to yield a substantially higher binding affinity to NaV1.2.

Modelling of NaV1.2 Channel
We use the cryo-EM structure of human NaV1.2-KIIIA (PDB ID: 6J8E) [17] in the MD simulations. As KIIIA binds at the pore domain of NaV1.2, we included only the S5-S6 helices and the connecting P-loops and helices in the computations, which correspond to the residues 250-435 in domain I, 884-983 in domain II, 1341-1479 in domain III, and 1655-1780 in domain IV. In the cryo-EM structure of NaV1.2, there is a missing piece in the extracellular loop of domain I (residues 285-313) that needs to be included in the structure (Figure 8). The homology model of the NaV1.2 pore domain was created with SWISS-MODEL [19]. The sequence and coordinates of the NaV1.2 pore domain were submitted as the initial parameters. The target-template alignments were obtained with built-in BLAST [20] and HHblits [21] (Supplementary Material, Table S1). Based on the alignment, the model was generated using ProMod3 [22], following an energy minimization process ( Figure 8).

Modelling of Na V 1.2 Channel
We use the cryo-EM structure of human Na V 1.2-KIIIA (PDB ID: 6J8E) [17] in the MD simulations. As KIIIA binds at the pore domain of Na V 1.2, we included only the S5-S6 helices and the connecting P-loops and helices in the computations, which correspond to the residues 250-435 in domain I, 884-983 in domain II, 1341-1479 in domain III, and 1655-1780 in domain IV. In the cryo-EM structure of Na V 1.2, there is a missing piece in the extracellular loop of domain I (residues 285-313) that needs to be included in the structure (Figure 8). The homology model of the Na V 1.2 pore domain was created with SWISS-MODEL [19]. The sequence and coordinates of the Na V 1.2 pore domain were submitted as the initial parameters. The target-template alignments were obtained with built-in BLAST [20] and HHblits [21] (Supplementary Materials, Table S1). Based on the alignment, the model was generated using ProMod3 [22], following an energy minimization process (Figure 8).

Modelling of NaV1.2 Channel
We use the cryo-EM structure of human NaV1.2-KIIIA (PDB ID: 6J8E) [17] in the MD simulations. As KIIIA binds at the pore domain of NaV1.2, we included only the S5-S6 helices and the connecting P-loops and helices in the computations, which correspond to the residues 250-435 in domain I, 884-983 in domain II, 1341-1479 in domain III, and 1655-1780 in domain IV. In the cryo-EM structure of NaV1.2, there is a missing piece in the extracellular loop of domain I (residues 285-313) that needs to be included in the structure (Figure 8). The homology model of the NaV1.2 pore domain was created with SWISS-MODEL [19]. The sequence and coordinates of the NaV1.2 pore domain were submitted as the initial parameters. The target-template alignments were obtained with built-in BLAST [20] and HHblits [21] (Supplementary Material, Table S1). Based on the alignment, the model was generated using ProMod3 [22], following an energy minimization process ( Figure 8).  In order to refine the Na V 1.2 channel model and check its stability, we performed the MD simulations using the NAMD code [23], employing the protocols that were established in previous MD simulations of Na V channels [24]. Briefly, the channel-toxin complex was embedded in a lipid bilayer and hydrated with a 0.15 M NaCl solution. The simulation system was then relaxed in several steps, gradually relaxing all the restraints on the atoms of the complex system. All the MD simulations were run until equilibration of the system was ensured (typically several hundred nanoseconds).

3.2.
Search for KIIIA Mutations to Improve Its Affinity for Na V 1.2 Excluding the cysteines that are involved in disulfide bonding and those KIIIA residues that already make contact with Na V 1.2 leaves only S5, S6, and S13 as potential mutation sites to improve KIIIA's affinity for Na V 1.2. We inspected the channel residues in the vicinity of S5, S6, and S13 in the Na V 1.2-KIIIA model, looking for potential charge interactions. We also would like the analogue to achieve a complete block of the channel. S5 is located within the N-terminal region, facing the extracellular loop in domain IV, and pointing towards the residues E1688 and D1692 ( Figure 9A). We considered a mutation of S5 to K5 and R5 but ultimately selected R5, which provides more stable contacts with domain IV. For S6, N333 is the only channel residue nearby ( Figure 9B), and its mutation to D6 is promising to provide a contact interaction between D6-N333. For the E6 mutation, the side chain is too long to make a good fit with N333. S13 is in the C-terminal region of KIIIA, and its side chain is located between D917 and E919 in domain II ( Figure 9C). A mutation of S13 to K13 could provide a perfect fit for the interaction of K13 with both D917 and E919. Therefore, we first considered a mutation of S5R, and then include the S13K and S6D mutations. The Na V 1.2-KIIIA analogue complexes were constructed from that of Na V 1.2-KIIIA using the mutator plugin in VMD and refined with the MD simulations.
Mar. Drugs 2022, 20, x FOR PEER REVIEW 9 of 11 In order to refine the NaV1.2 channel model and check its stability, we performed the MD simulations using the NAMD code [23], employing the protocols that were established in previous MD simulations of NaV channels [24]. Briefly, the channel-toxin complex was embedded in a lipid bilayer and hydrated with a 0.15 M NaCl solution. The simulation system was then relaxed in several steps, gradually relaxing all the restraints on the atoms of the complex system. All the MD simulations were run until equilibration of the system was ensured (typically several hundred nanoseconds).

Search for KIIIA Mutations to Improve Its Affinity for NaV1.2
Excluding the cysteines that are involved in disulfide bonding and those KIIIA residues that already make contact with NaV1.2 leaves only S5, S6, and S13 as potential mutation sites to improve KIIIA's affinity for NaV1.2. We inspected the channel residues in the vicinity of S5, S6, and S13 in the NaV1.2-KIIIA model, looking for potential charge interactions. We also would like the analogue to achieve a complete block of the channel. S5 is located within the N-terminal region, facing the extracellular loop in domain IV, and pointing towards the residues E1688 and D1692 ( Figure 9A). We considered a mutation of S5 to K5 and R5 but ultimately selected R5, which provides more stable contacts with domain IV. For S6, N333 is the only channel residue nearby ( Figure 9B), and its mutation to D6 is promising to provide a contact interaction between D6-N333. For the E6 mutation, the side chain is too long to make a good fit with N333. S13 is in the C-terminal region of KIIIA, and its side chain is located between D917 and E919 in domain II ( Figure 9C). A mutation of S13 to K13 could provide a perfect fit for the interaction of K13 with both D917 and E919. Therefore, we first considered a mutation of S5R, and then include the S13K and S6D mutations. The NaV1.2-KIIIA analogue complexes were constructed from that of NaV1.2-KIIIA using the mutator plugin in VMD and refined with the MD simulations.

Conclusions
In conclusion, we used MD simulations to refine the Na V 1.2-KIIIA complex and compared the pairwise interactions with that in the cryo-EM structure. The simulations provide a comprehensive understanding of the structural dynamics of KIIIA binding to Na V 1.2, which are overall in agreement with the reported cryo-EM structure. One feature of the computational modeling are the differences in the R14 binding site. The cryo-EM structure shows that the R14 side chain stably interacts with E919 on the extracellular P1 helix in domain II [17], while our MD simulations indicate that R14 forms a stronger ionic bond with E1444 in domain III. As R14 is identified as a critical residue in KIIIA binding [15], its stronger binding that was observed in our simulations is more consistent with the R14A mutation data.
The Na V 1.2-KIIIA complex that was constructed here has been used in designing KIIIA analogues with improved binding affinity and complete channel blockage. The S5R mutation forms stable interactions with D1692 in domain IV, which block the ion permeation pathway. The S13K side chain extends towards the extracellular loops in domain II and makes strong ionic bonds with D917 and E919. Both interactions are expected to enhance the binding affinity considerably. However, two positive charges are introduced by the [S5R, S13K] double mutation which may affect the drug absorption. The S6D mutation partly neutralizes the extra positive charges and also its side chain interacts with the adjacent N333 in domain I. The final KIIIA[S5R, S6D, S13K] analogue binds to all four domains of the Na V 1.2, which results in a complete block of the pore and an increased binding affinity. Our simulation results suggest new potent KIIIA analogues for the blocking of Na V 1.2. Once tested and confirmed experimentally, they could provide promising leads for anti-epileptic drugs.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/md20020154/s1, Table S1: Sequence alignment used in homology modeling, Figure S1: Time series of the R14 contacts with the channel residues.