Structure based analysis of KATP channel with a DEND syndrome mutation in murine skeletal muscle

Developmental delay, epilepsy, and neonatal diabetes (DEND) syndrome, the most severe end of neonatal diabetes mellitus, is caused by mutation in the ATP-sensitive potassium (KATP) channel. In addition to diabetes, DEND patients present muscle weakness as one of the symptoms, and although the muscle weakness is considered to originate in the brain, the pathological effects of mutated KATP channels in skeletal muscle remain elusive. Here, we describe the local effects of the KATP channel on muscle by expressing the mutation present in the KATP channels of the DEND syndrome in the murine skeletal muscle cell line C2C12 in combination with computer simulation. The present study revealed that the DEND mutation can lead to a hyperpolarized state of the muscle cell membrane, and molecular dynamics simulations based on a recently reported high-resolution structure provide an explanation as to why the mutation reduces ATP sensitivity and reveal the changes in the local interactions between ATP molecules and the channel.

. Membrane potential and MgATP sensitivity of Kir6.2/R50P myotubes. (a-c) Mean membrane potential (± SEM) of C2C12 myotube (a), before (white) and after (black) application of glibenclamide (b), and ACh response (c) without infection (N = 5-6), with infection of the Kir6.2/WT AAV vector (N = 5-6), the Kir6.2/R50P AAV vector (N = 5-6). Multiple groups were compared using an ordinary One-Way ANOVA test, followed by Tukey's test, and two groups were compared using Paired t test. www.nature.com/scientificreports/ results indicate that R50P mutation is responsible for the hyperpolarized shift of membrane potential observed in the present study. Also, the MgATP dose-response curves show that Kir6.2/R50P-expressing myotubes have less ATP sensitivity compared to Kir6.2/WT-expressing or non-transfected myotubes (Fig. 1d). Taken together, the R50P mutation induces a reduction of the membrane potential within the physiological range of intracellular MgATP concentrations, and this in turn may result in a reduction of the acetylcholine (ACh) stimuli. Therefore, we subsequently examined whether depolarization induced by 10 mM acetylcholine (ACh) differ in Kir6.2/R50P-expressing myotubes (Fig. 1c). We concluded that there were no significant differences when examining on the cell line basis in response to ACh (− 12.2 ± 2.1 mV for non-transfected cells; − 12.2 ± 3.8 mV for Kir6.2/WT-expressing cells; − 10.9 ± 4.8 mV for Kir6.2/R50P-expressing cells). In addition, the results of the Ca 2+ release from the intracellular Ca 2+ store in the endoplasmic reticulum in Kir6.2/R50P-expressing myotubes showed that there is not a significant difference between Kir6.2/WT-and Kir6.2/R50P-expressing myotubes in response to the same concentration of ACh (Fig. 2). As downstream factor, the qRT-PCR of Baf60c, Deptor and AKT2 revealed that no significant difference was found in these factors in our cell line system (Fig. S1).

Structural analyses and molecular dynamics simulations of R50P polymorphism in Kir6.2. The
R50P-mutated K ATP channel exhibits gain-of-function properties with reduced ATP sensitivity 17,18 . Structural analyses from recently resolved structures 37,40 (PDB IDs: 6BAA and 6JB1) show that the Arg50 residue is located at a solvent exposed loop in the intracellular region of the Kir6.2 subunit, and that this residue facilitates wrapping around the nearby Kir6.2 subunit (Fig. 3a). Due to the spatial proximity, the side chain of Arg50 establishes electrostatic interactions with the phosphate groups in ATP, facilitating ATP sensing in the cytosol. To further characterize the ATP binding site and the interactions between ATP and the channel residues, molecular dynamics simulations were performed on Kir6.2/WT and Kir6.2/R50P inserted in a model membrane of POPC (1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine). A simple measure of the stability of the structure during the dynamics is given by the drift from its initial conformation provided by the root-mean-square deviation (RMSD) between the channel structure at a given time and the initial structure. Analysis of the alpha carbon atom RMSD for the entire tetramer, and the M1 (residues 67-97) and M2 (residues 142-172) helices versus time revealed a similar pattern in each of the simulations (Fig. S5). The initial structural drift reflects the fact that the restraints imposed during the equilibration period have been switched off and that the protein is free to relax within the lipid bilayer and the solvent. The overall RMSD values are the same for all the simulations suggesting overall conformational stability of the protein at these timescales.
While Martin et al. reported that Arg50 interacts with only the gamma phosphate group (PDB ID: 6BAA; 3.6 Å-resolution) 37 , Chen et al.reported that the residue extensively interacts with the beta and gamma phosphate (e) Mean firing frequency (± SEM) before and after 10 mM ACh stimulation at 5 min (N = 6-12). Multiple groups were compared by an ordinary one-way ANOVA test, followed by Tukey's test (**p < 0.01; ****p < 0.001). The figures were analyzed and created with GraphPad Prism software (GraphPad Prism 8.2.1; https:// www. graph pad. com/ scien tific-softw are/ prism/).  40 (Fig. 3a). In line with Chen's reports, the hydrogen atoms of the sidechain of Arg50 interact with the beta and gamma phosphate groups of ATP in our simulations of the wild type ( Fig. 3d,f,g, Fig. S6). In addition, interactions with the alpha phosphate group of ATP are also recorded (Fig. 3b,c,e). Replacement of Arg50 with a Pro residue results in the loss of the extensive electrostatic interactions established between Arg50 and the ATP molecules. Furthermore, the relatively rigid five-membered ring of the Pro replacement drastically limits the range of the values of the torsion angle φ (rotation about N-C α bond of the peptide backbone) to approximately − 60°, leading to loss of flexibility.  www.nature.com/scientificreports/ During the simulations, interactions between ATP and nearby residues were recorded as a function of time. It was found that residues Lys185, Lys39 and Arg54 as well as Arg50 in the wild-type and Gln52 in the mutant are in close contact with the ATP molecules at the beginning of the simulations. Some of these interactions are maintained during the entire simulation time and in some other instances, there is an exchanged between the pair of atoms of a given residue that interact with the ATP molecule ( Fig. 3b-d). While four ATP molecules remain bound to the wild-type Kir6.2 through their robust links with residues Lys185 and Arg50, and to some extend Arg54 and Lys39 in the absence of SUR1 subunits, the majority of the ATP-protein interactions diminishes as the simulation progresses in the case of the R50P mutant simulations (Table 1). Overall, the association between ATP and the protein is weaker in the R50P mutant compared to the wild type as the interactions described between the oxygen atoms from the phosphate groups of ATP and the hydrogen atoms of the guanidinium side chain of arginine are impaired (Fig. 4, Fig. S6). During the simulations, the side chain of Lys39 interacts with the phosphate groups of ATP, facilitating the association of ATP with the protein although the final orientation of ATP is far from its original state. Notably, residues Arg50 and Arg54 of the same monomer are in direct contact with a SUR subunit and likewise, residues Lys39 and Lys185 each in a nearby monomer also interact directly with the same SUR subunit. Only the interactions between Lys185 and the ATP molecules remain in place in both the wild type and the mutant in the four monomers of the protein during the whole time in all the replicas Table 1. Summary of simulations performed and the interactions established between ATP molecules and the Kir6.2 channel residues after 500 ns of simulation time in the wild type and mutant systems considered. Arg50 is used for ATP binding in all four monomers in wild type.

System Replica Number Protein-ATP interactions (presence in # monomers)
Wild type Figure 4. Comparison between before and after molecular dynamics simulation in R50P mutant. Before simulation, four ATP molecules were accommodated in the cavity (a), whereas after 500 ns, each adenosine moiety in the ATP is exposed to the solvent (b). The ATP molecule is weakly bound to the channel with three residues Lys39, Arg54, Lys185. The simulation was calculated by gromacs v2020. www.nature.com/scientificreports/ computed, although there are signs of weakening as the time progresses. The simulation results agree with the observation that the mutation K185Q is also involved in the DEND syndrome 44 . Although the simulations were carried out in the absence of the SUR subunits, altering these SUR-Kir interactions could potentially affect the function of the complex.

Effect of myogenesis and glucose uptake in Kir6.2/R50P-expressing myotubes.
To examine the physiological effects of the R50P mutation of the K ATP channel on skeletal muscle, myogenesis and glucose uptake measurements were performed. As shown in Fig. 5, the myogenic index values between non-transfected Kir6.2/WT-and Kir6.2/R50P-expressing murine skeletal muscles were 0.168 ± 0.029, 0.164 ± 0.045, 0.163 ± 0.029, respectively, which suggests that K ATP channels are not involved in the myogenic index. At the same time, myotube width analyses show that both Kir6.2 (WT)-expressing and Kir6.2 (R50P)-expressing myotubes exhibit almost similar width distribution with a maximum peak at 27-30 μm. Taken together, these results indicate that K ATP channels are not involved in myogenesis on a cell line basis. Subsequently, glucose uptake was evaluated with confocal laser-scanning microscopy. The mean average values of the green fluorescence were recorded for the glucose uptake in Kir6.2/WT-and Kir6.2/R50P-expressing murine skeletal myotubes. Figure 6 shows that, in the presence of insulin, the mean average value of the green fluorescence of the Kir6.2/R50P-expressing myotubes is almost similar to that of the Kir6.2/WT-expressing myotubes, indicating that K ATP channels are not involved in glucose uptake on a cell line basis.

Discussion
A previous study demonstrated that muscle dysfunction was of neuronal origin by comparing muscle-and nerve-specific expression of the Kir6.2/V59M mutation in transgenic mice (named m-V59M and n-V59M mice, respectively), the most commonly known mutation inducing intermediate DEND syndrome 19 . Interestingly, there was no significant difference in muscle membrane potential between wild-type and m-V59M mice or wild-type and n-V59M mice, indicating that the gain-of-function mutation in the K ATP channel does not directly affect muscle function itself. However, the present study shows that the R50P mutation, DEND syndrome inducing mutation with more severe symptoms than the V59M, induces a hyperpolarized membrane potential. Although the muscle weakness is mainly originated in the brain, our present data indicate that severely reduced ATP sensitivity in the muscle K ATP channels may take part to develop muscle weakness in DEND syndrome patients. From the recent K ATP channel structures solved at atomic resolution, it is inferred that Arg50 can play pivotal roles in the recognition of nucleotides by establishing electrostatic interactions with ATP, and it is indeed part of the ATP binding pocket in Kir6.2 40 . It is well established that this mutation decreases ATP sensitivity without   17 . Our present simulation results further confirmed that R50P mutation indeed contributes to reduce the ATP sensitivity. It is of interest because the simulation enabled us to see the serial interaction between ATP and its binding residues while released from the binding pocket (such as the K185Q mutation described above), which facilitates the understanding of molecular bases of ATP binding mechanism and development of neonatal diabetes/DEND syndrome. In contrast, Val59, which has been reported to affect ATP sensitivity with increased open probability (denoted "gating mutation") 19 is located at a distal place, > 18 Å away from the ATP location (PDB ID: 6JB1) and not contributing directly to ATP binding. Consequently, the effect on the activation mechanism of each of these two mutations seems to be rather different from an atomistic perspective. Inside-out patch clamp experiments demonstrated that the C-terminal tail comprised by 42 a.a. residues in SUR2A possess an inhibitory effect on NBD2-mediated ADP-induced channel activation 45 . In line with these observations, higher concentrations of ADP were required to activate K ATP channels that contained more SUR2A than SUR1 or SUR2B regulatory domains 46 . As a consequence, the K ATP channel with the SUR2A subunit tended to be found in a closed state. Therefore, it is possible that a "binding mutation", such as R50P, has higher impact on ATP sensitivity in SUR2A-containing K ATP channels compared to "gating mutations" such as V59M, which showed no difference in membrane potential in skeletal muscle. As suggested by the recently reported K ATP channel atomic structure [35][36][37][38][39][40] , SUR1 adopts an inward-facing conformation which is an inactive form, with a folded lasso motif that facilitates ATP binding in Kir6.2. In contrast, the lasso motif transforms to an extended form when the complex adopts an outward-facing structure. Considering these observations, given that SUR2A adopts an inward-facing conformation more likely than SUR1, the R50P mutation may disrupt SUR2A regulation more drastically than the SUR1 regulation, which is consistent with the previous observations 46 . Consequently, the R50P "binding mutation" in Kir6.2 may have a higher impact on the membrane potential than the "gating mutation" in skeletal muscles, or even in cardiomyocytes [47][48][49][50] .
In this study, we used C2C12 cells with transfected WT and R50P-mutated Kir6.2. Since Kir6.2 and SUR are both required to be expressed in cell surface, transfected Kir6.2 channel surface expression should be under the regulation of endogenous SUR. Because except for the rare case, most of the reported neonatal diabetes/DEND syndrome patients are heterozygous with mutated channel subunit 51 , we consider that our system is reflecting the heterozygous state of the patients without excessive expression of K ATP channels in cell surface.
In previous study, ATP sensitivity experiments show that in pancreatic Kir6.2/R50P-SUR1 heteromeric K ATP channel, approximately 33% channel current is unblockable 17 . However, high dose of MgATP was able to fully www.nature.com/scientificreports/ block K ATP channel current in Kir6.2/R50P transfected C2C12 cells (Fig. 1d). Therefore, although we need to keep in mind that the current we have recorded in this study should reflect the mixture of Kir6.2/SUR1 and Kir6.2/ SUR2A (Fig. S1), it can be considered that R50P mutation has less impact on skeletal muscle compared to that of pancreas and brain, and that the R50P mutated neuronal K ATP channel may mainly contribute to muscle weakness. Nevertheless, this study shows that R50P mutated sarcoK ATP channel does contribute to hyperpolarization unlike V59M mutated sarcoK ATP channel. Furthermore, this study revealed that the application of glibenclamide shows no effects on non-transfected and Kir6.2 (WT)-expressing myotubes, whereas the application partly restores the hyperporalized state in Kir6.2 (R50P)-expressing myotubes. These results indicate that in skeletal muscle, wild-type K ATP channels are closed in physiological state with very small contribution on membrane potential and also indicating the possible effects of glibenclamide on DEND syndrome patients. Also, we need to keep in mind that glibenclamide is reported to act on activate exchange protein directly activated by cAMP 2 (Epac2A) or carnitine palmitoyltransferase 1 (CPT-1) 52-55 . Epac2A is not expressed in skeletal muscle but CPT-1 is known to be expressed in skeletal muscle and functions as the key regulatory enzyme of mitochondrial longchain fatty acid oxidation. Therefore, further study is required to clarify the targets of glibenclamide other than K ATP channels and its function in skeletal muscle. It is possible that the hyperpolarization caused by the R50P mutation may affect downstream targets of muscle cells, and Meng et al. has shown that Baf60c-Deptor-AKT2 axis is affected by the activities of K ATP channel in skeletal muscle 56 . Therefore, we have investigated whether opening of K ATP channel induced by R50P mutation has any effect on the downstream signaling of the channel. The Baf60c-Deptor-AKT2 axis was found not to be affected by the R50P mutation (Fig. S1). Furthermore, the intracellular Ca 2+ levels and glucose uptake were also not affected by this mutation (Figs. 2, 6). This was of surprise since previous studies using Kir6.2 or SUR2 KO mice showed changes in intracellular Ca 2+ or glucose uptake 57,58 .
As for the limitation, data from the transgenic mice with Kir6.2 (R50P) or any other neonatal diabetes inducing mutation are not available in this study. However, the present study provided insight on the links between muscle weakness and sarcoK ATP channel. Further study using a transgenic model of the muscle-specific Kir6.2(R50P)-expressing mouse is required to clarify the mechanism of muscle weakness in the DEND syndrome patients.
Structural analyses and molecular dynamics simulations. The K ATP channel structures were selected from the Protein Data Bank (https:// www. rcsb. org), PDB IDs: 6C3O; 6C3P; 6BAA; 5WUA; 5TWV; 5YKE; 5YKF; 5YKG; 5YWA; 5YWB; 5YWC; 5YW8; 5YW9; 6JB1, all of which were obtained by cryogenic electron microscopy [35][36][37][38][39][40] . Structures with resolution higher than 4 Å were selected: 6C3O (3.9 Å), 6BAA (3.6 Å), 6JB1 (3.3 Å). The Arg50 residue in Kir6.2 is only assigned in PDBs 6BAA (3.6 Å) and 6JB1 (3.3 Å), and thus, these two structures were selected for structural analyses using PyMol (http:// pymol. sourc eforge. net/). To begin www.nature.com/scientificreports/ to address the structural mechanisms by which the mutations considered in this study alter the action of ATP on the Kir6.2 + -channel, molecular dynamics simulations were carried out using the wild type channel and the mutant R50P. The crystal structure of Kir6.2 closed at its inner gate was retrieved from the Protein Data Bank (PDB ID 6JB1) 40 and resolved residues 26-114 were used for subsequent modelling. This particular crystal structure of Kir6.2 was chosen for this study over other available structures of the channel because it has a good atomic resolution, the two residues focused of the study are resolved, and it was obtained in the presence of the SUR subunits. R50P mutant channels were generated using the Mutator plugin of VMD. N-and C-termini were acetylated and methylated, respectively. The VMD solvate plugin was used to solvate internal cavities of the protein. The structures were aligned perpendicular to the bilayer and inserted into a membrane, a neutral membrane containing 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC), with x and y dimensions of 96 Å. The VMD solvate plugin 60 was then used to create a rectangular water box around the membrane-protein system. The overlapping water and lipid molecules around the ion channel structure were removed with the cutoff distance (1.2 Å). Potassium and chloride ions were added using Autoionize Plugin of VMD to neutralise the systems and obtain a concentration of 150 mM. The final system size was approximately 183,000 atoms. MD simulations were performed with the gromacs v2020.1 software (https:// doi. org/ 10. 5281/ zenodo. 35625 12). Four replicates of 500 ns were run for the wild type and R50P mutant amounting to 4 µs in total. CHARMM36 61 parameters were used for the protein, ions, ATP and lipids, and the TIP3P model was used for water. The particle mesh Ewald method was used for the treatment of periodic electrostatic interactions 62 , with an upper threshold of 1 Å for grid spacing. Electrostatic and van der Waals forces were calculated every time step. A cut-off distance of 12 10 Å was used for van der Waals forces. A switching distance of 10 Å was chosen to smoothly truncate the non-bonded interactions. Only atoms in a Verlet pair list with a cut-off distance of 13.5 Å reassigned every 20 steps were considered. The SETTLE LINCS algorithm 63 was used to constrain all bonds involving hydrogen atoms, to allow the use of a 2 fs time step throughout the simulation. MD simulations were performed in the NPT ensemble. The Parrinello-Rahman barostat 64 was employed to control the pressure with a 50 fs period, and the desired value of 1 atmosphere. The system was coupled to a Nose-Hoover thermostat 65,66 to sustain a temperature of 310.15 K throughout.
The same equilibration protocol was applied to all the systems. The systems were subjected to 10,000 steps of minimization, with harmonic constraints on protein backbone atoms and protein sidechains (20 kcal/mol/ Å 2 ), the ATP (force constant 40 kcal/mol/Å 2 ), and lipid headgroups (10 kcal/mol/Å 2 ). Harmonic restraints were gradually reduced to a force constant of 2 kcal/mol/Å 2 and were removed in consecutive steps from the lipid headgroups, protein sidechains and protein backbone over the course of a 10 ns trajectory. Production runs followed and analysis of the simulations was performed with python v3.7.3 scripts using the MDTraj v1.9.4 python package.
Electrophysiology. Electrophysiological analyses were performed according to the previous study 67 . Following incubation for 48 h after medium change to differentiation medium, 10 μL of AAV vectors carrying Kir6.2/WT or Kir6.2/R50P were added to the cultured dishes. Subsequently, after infection for 72 h, the medium was replaced with fresh medium, and incubated for an additional 5-6 days. All electrophysiological measurements were performed at room temperature (22-25 °C) using an EPC-800 patch clamp amplifier (HEKA Electronics, Lambrecht/Pfalz, Germany) and pCLAMP 10 software (Molecular Devices, CA, USA). The membrane potentials were recorded in the current clamp mode in standard whole-cell technique 68  where [MgATP] is the MgATP concentration, IC 50 is the [MgATP] at which inhibition is half-maximal and h is the slope factor (Hill coefficient). Data were analyzed using Origin Graphing and Data Analysis Software (Version 9, OriginLab Corporation, MA; https:// www. origi nlab. com).

Criteria for [Ca 2+
] i responses. [Ca 2+ ] i was measured as previously reported 69 . Mouse skeletal myotubes were incubated with 2 μmol/L Fura-2/AM (Dojin Chemical, Kumamoto, Japan) for 30 min at room temperature, mounted in a chamber and superfused with KRB at 1 mL/min at 37 °C. Fluorescence images due to excitation at 340 and 380 nm were captured and the ratio (F340/F380) images were produced by an Argus-50 system (Hamamatsu Photonics, Hamamatsu, Japan www.nature.com/scientificreports/ medium was replaced with fresh medium, and incubated for an additional 24 h. After removing the incubation media, the C2C12 cells were washed in cold PBS fixed in 100% methanol. The DAPI-containing media were mounted on a coating dish for confocal laser-scanning microscope (Fluoview FV10i; Olympus, Tokyo, Japan) and used to allow clear visualization of nuclei and myotube structures for quantitative measurements; the dishes were then covered with a cover glass. The myogenic index was used as a morphological parameter of muscle differentiation 70 . The viral infected myoblasts and myotobes were also confirmed by RFP fluorescence automatically cleaved from Kir6.2/WT or Kir6.2/R50P inside the cytosol. The number of nuclei in each myotube containing ≥ 3 nuclei and the total number of nuclei were counted in 30 randomly selected fields per well. The myogenic index (in %) was then calculated as: ([number of nuclei in myotubes]/[total number of nuclei in 30 randomly selected fields]) using ImageJ software version 1.51 (National Institutes of Health, Bethesda, MD, USA; https:// imagej. nih. gov/ ij/ downl oad. html). Myotube widths were measured by using a confocal laser-scanning microscope or the ImageJ software. In brief, cells were evaluated in 30 randomly selected fields per well. The width of each myotube containing ≥ 3 nuclei was measured on the cell, and the average width per myotube was calculated. A total of 60-90 myotubes was evaluated for each treatment group.
qRT-PCR and immunohistochemical analyses. Polymerase chain reaction amplification was performed by using a StepOne Real-time PCR system (Life Technologies, CA, USA) under the following conditions: one cycle of 95 °C for 30 s; 45 cycles of 95 °C for 5 s and 56 °C for 10 s and 72 °C for 15 s. The amplified fragments were applied to 2% (w/v) agarose gel, and visualized by staining with SYBR Safe DNA staining reagent (Life Technologies, CA, USA). In immunostaining experiments, the glass bottom dishes, which were cultured transfected and non-transfected myotubes, were washed in phosphate buffered saline (PBS) and incubated in blocking solution containing 2% bovine serum albumin (BSA) and 5% normal goat serum for 1 h. Then, the dishes were incubated with rabbit polyclonal Kir6.2 antibody (1:1000, ab81121, abcam, Cambridge, UK) for 48 h at 4 °C. Dishes were washed in PBS and incubated with Alexa 488-labelled anti-rabbit IgG (Invitrogen, CA, USA) for 30 min and mounted with a DAPI-containing mount medium and covered with a cover glass.

Statistical analysis.
As previously reported 71 , the data are presented as the means ± standard error of the mean (SEM) unless otherwise indicated. One-way analysis of variance (ANOVA) was used, and the treatment groups were compared with Tukey's post hoc test for honest significant difference. All statistical analyses were performed using GraphPad Prism 8.2.1 (GraphPad Software, Inc., San Diego, CA, USA; https:// www. graph pad. com/ scien tific-softw are/ prism/) for MacOS, version 10.14.6. A p value of < 0.05 was considered statistically significant. www.nature.com/scientificreports/