Deregulated ion channels contribute to RHOBTB2-associated developmental and epileptic encephalopathy

Abstract

While de novo missense variants in the BTB domains of atypical RhoGTPase RHOBTB2 cause a severe developmental and epileptic encephalopathy, de novo missense variants in the GTPase domain or bi-allelic truncating variants are associated with more variable neurodevelopmental and seizure phenotypes. Apart from the observation of RHOBTB2 abundance resulting from BTB-domain variants and increased seizure susceptibility in Drosophila overexpressing RhoBTB, our knowledge on RHOBTB2-related pathomechanisms is limited. We now found enrichment for ion channels among the differentially expressed genes from RNA-Seq on fly heads overexpressing RhoBTB. Subsequent genetic interaction experiments confirmed a functional link between RhoBTB and paralytic, the orthologue of human sodium channels, including epilepsy associated SCN1A, in vivo. We then performed patch-clamp recordings on mature neurons differentiated from human induced pluripotent stem cells with either homozygous frameshifts or patient-specific heterozygous missense variants in the GTPase or the BTB domains. This revealed significantly altered neuronal activity and excitability resulting from BTB domain variants but not from GTPase domain variants or upon complete loss of RHOBTB2. Our study indicates a role of deregulated ion channels in the pathogenesis of RHOBTB2-related developmental and epileptic encephalopathy and points to specific pathomechanisms underlying the observed genotype–phenotype correlations regarding variant zygosity, location and nature.

developmental and epileptic encephalopathy, iPSC derived neurons, RHOBTB2, Drosophila, genetic interaction

Introduction

Recently, we and others identified variants in RHOBTB2 (RHO-related BTB domain-containing protein 2) as causative for variable neurodevelopmental phenotypes. RHOBTB2 belongs to an atypical subfamily of the evolutionary conserved Rho family of small GTPases [1–3]. All three members (RHOBTB1-3) share a distinct architecture with a RHO GTPase domain, a proline rich region, two BTB domains and a C-terminal region [1, 2]. In contrast to typical Rho GTPases, there is no specific guanine nucleotide exchange or activating factor known to interact with and activate RHOBTB2. The BTB domains are involved in intra- and intermolecular protein–protein interaction [1, 4].

RHOBTB2 interacts with the ubiquitin ligase scaffold protein Cullin3, allowing the assembly into an ubiquitin ligase complex and initiating auto-ubiquitination of RHOBTB2 and ubiquitination of substrates recruited by RHOBTB2 to the complex and directing them to the proteasome for degradation [1, 5].

In accordance with its high expression in brain [2, 6], RHOBTB2 has been implicated in various neurodevelopmental and epilepsy phenotypes during the last years. De novo, heterozygous missense variants in the BTB domain region were identified to cause a developmental and epileptic encephalopathy (DEE64, MIM#618004) with consistent presentation of early-onset seizures, severe to profound intellectual disability, microcephaly and movement disorders [7–17]. De novo, heterozygous missense variants in the GTPase domain and outside the domains are associated with more variable neurodevelopmental and neurological phenotypes with or without seizures [12, 16, 18, 19]. Furthermore, bi-allelic truncating variants were reported in several families with variable intellectual disability and commonly seizures [19, 20].

Missense variants in the BTB domain region resulted in increased abundance of RHOBTB2 in vitro, probably due to impaired proteasomal degradation [7]. In accordance, flies overexpressing the Drosophila ortholog RhoBTB pan-neuronally, presented with seizure susceptibility [7]. As missense variants outside the BTB domain region did not result in such protein abundance, and as for truncating variants a loss-of-function effect is assumed, these observations indicated a genotype–phenotype correlation based on location, zygosity and nature of variants and potentially different functional consequences [19].

However, so far, little is known about the role of RHOBTB2 and its particular variant classes in neuronal function and dysfunction and about the pathomechanism of RHOBTB2-related seizures.

We now found enrichment of genes encoding ion channels among the differentially expressed genes in transcriptome analysis from flies overexpressing RhoBTB. By genetic interaction experiments in Drosophila we confirmed a functional link between RhoBTB and paralytic (sodium channel) as well as possibly ir76a (ionotropic glutamate receptor). In a human model, we performed patch clamp recordings and observed altered electrophysiological behavior of hiPSC- (human induced pluripotent stem cells) derived neurons carrying heterozygous missense variants in the BTB domain region of RHOBTB2, but not in neurons carrying missense variants in the GTPase domain or bi-allelic truncating variants. Our study therefore implicates deregulated ion channels specifically in the pathomechanism of RHOBTB2-related DEE due to heterozygous missense variants in the BTB domain region and suggests different pathomechanisms for other RHOBTB2-related neurological and neurodevelopmental phenotypes associated with missense variants in the GTPase domain or bi-allelic truncating variants.

Results

Differentially expressed genes are enriched for ion channels in flies overexpressing RhoBTB

As increased RHOBTB2 levels were observed in association with variants in the BTB domain region, and as accordingly flies overexpressing RhoBTB recapitulate the seizure phenotype [7], we utilized Drosophila as a model to further study RHOBTB2-related pathomechanisms. As substrates of ubiquitin ligases have often been found to be involved in transcriptional regulation [21, 22], alterations of ubiquitin ligase RhoBTB may influence transcriptional profiles as well [23]. To test this hypothesis, we performed RNA sequencing on heads of flies with pan-neuronal overexpression of RhoBTB (n = 6, 3 wildtype (WT) and 3 OE samples).

Gene expression in RNA from heads of control and RhoBTB overexpressing flies could be clearly separated along the first principal component in the principal component analysis (PCA), albeit with some variation among the RhoBTB samples (Supplementary Fig. S1A). Nevertheless, we found a broad deregulation of gene expression in RhoBTB overexpressing flies. Of 1082 differentially expressed genes, 331 were downregulated, and 751 were upregulated with an adjusted p-value of 0.05 (Fig. 1A, Supplementary Fig. S1B). For identification of commonly deregulated pathways, we performed gene enrichment analysis for down—and upregulated genes using the PANTHER knowledge base [24, 25]. All significantly enriched molecular functions among the downregulated genes are linked to neurotransmitter receptor activity, ion channel activity or transmembrane transporter activity. Biological processes enriched are linked to transmembrane transport, signaling in general and signal transmission more specifically. Correspondingly, cell junctions in general and synapse specifically were the only enriched cellular components (Fig. 1B, Supplementary Fig. S1C and D, and Supplementary Table S7 and Supplementary Tables S1 and S2). In contrast, no clear pattern could be defined for ontology terms enriched among upregulated genes (Supplementary Fig. S1E and F). As individuals with variants in the BTB domain region of RHOBTB2 present with a severe DEE phenotype, and as DEEs have frequently been linked to defects in various ion channels [26–29], we focused on deregulated genes encoding for ion channels. Of the 331 downregulated genes, 20 were involved in ion channel function, which represents a significant enrichment compared to ion channels among all expressed genes (OE = 2.7, P < 5.9 × 10⁻⁵). To better link differential gene expression in Drosophila heads to human disorders, we first identified human orthologs for all genes. Almost all downregulated ion channel genes had at least one ortholog in humans (19/20 genes). When looking at human orthologs that have been linked to a neurodevelopmental disorder (NDD) according to the SysNDD database [30], we found this to be the case for more than 40% of the downregulated ion channel encoding genes (Supplementary Table S1). From these we selected three genes that represented different ion channel classes, were already associated with NDDs and seizures [31–36] and for which transgenic fly lines were available (paralytic, slowpoke, ir76a).

Figure 1

Functional link between RhoBTB and ion channels in Drosophila melanogaster. (A) Heatmap displays the 1082 differentially expressed genes (padj < 0.05) between fly head samples overexpressing RhoBTB pan-neuronally (RhoBTB_OE, n = 3) and matching, wildtype control flies (WT, n = 3). (B) Enriched gene ontology (GO) terms among downregulated genes were grouped based on function and show enrichment of molecular functions (MF), biological processes (BP) and cellular components (CC) all relating to signaling and transmembrane transport of various molecules at cell junctions. Detailed results on individual enriched GO terms can be found in Supplementary Fig. S1. (C) Summary of genetic interaction screen between RhoBTB and three selected ion channel genes. Effect of additional manipulation of different genes on the RhoBTB associated phenotype in three assays is categorized as antagonistic, synergistic or additive. (X/Y) number of fly lines with phenotype or alteration of the RhoBTB associated phenotype out of the number of tested combinations, L: Not assessable due to lethality. For slowpoke, only an overexpression line of the human orthologue KCNMA1 was available and used in the screen. (D) Simultaneous manipulation of paralytic resulted in an amelioration of the RhoBTB overexpression associated bang sensitivity phenotype quantified by the number of flies with spasms after being vortexed. Asterisks indicate statistical significance (^*P < 0.05, ^**P < 0.01, ^***P < 0.001, Wilcoxon test with Bonferroni correction) for the different timepoints. Significance is only shown for comparison between RhoBTB_OE and the pairwise manipulations. (E) Simultaneous manipulation of paralytic resulted in an amelioration of the RhoBTB overexpression associated wing morphology phenotype, quantified by presence or absence of cross veins, exemplary pictures of the wings are shown. Asterisks indicate statistical significance (^*P < 0.05, ^**P < 0.01, ^***P < 0.001, chi-squared test with Bonferroni correction). Significance is shown for comparisons between RhoBTB_OE and the pairwise manipulations as well as wildtype and effects of single interactor gene manipulation.

Genetic interaction between RhoBTB and ion channel genes in Drosophila

To further investigate a functional link between RhoBTB and three selected ion channels (paralytic, slowpoke, ir76a), we performed genetic interaction experiments in Drosophila.

We first generated lines combining the transgenic RhoBTB element with either a pan-neuronal (elav-Gal4) or a wing (ms1096-Gal4) driver, thus stably overexpressing RhoBTB in either the nervous system (qPCR validation in Supplementary Fig. S2) or the wing. We confirmed the previously reported [7] bang sensitivity and negative geotaxis phenotypes (Supplementary Fig. S3A and B) and newly observed a wing phenotype with a significant number of one or both short cross veins lacking in male flies (Supplementary Fig. S3C and D). While wing phenotypes are not modelling a neurological phenotype, they are a suitable, quantifiable parameter to investigate genetic interactions or molecular commonalities between genes [30, 37]. We confirmed expression of the ion channel genes in wing tissue by reverse transcriptase (RT) PCR (Supplementary Fig. S2C). Upon crossing the stable overexpression line with a line containing a UAS-element from scrambled c. elegans cDNA, we did not observe a significant ameliorating effect on the phenotypes in any of the assays (Supplementary Fig. S3). This excludes dilution effects by distribution of GAL4 between two UAS elements.

Ubiquitous or—if lethal—pan-neuronal knockdown or overexpression of the three selected ion channel genes by one to four different transgenic RNAi- or overexpression lines was confirmed by quantitative RT-PCR (Supplementary Fig. S2). RNAi or overexpression lines not inducing knockdown or overexpression or lines with early lethality were excluded from further experiments.

Knockdown or overexpression of either of the three selected genes resulted in a variable phenotypic spectrum ranging from no effect to neuronal or wing phenotypes or lethality (Supplementary Table S3). By crossing any of these lines with the stable RhoBTB overexpression line, we could obtain between three to four pairwise combinations (Supplementary Table S3) that were assessed for phenotype alterations. The detailed results of this screen are listed in Supplementary Table S3 and are summarized below.

For the pairwise combination of RhoBTB and paralytic four different combinations (two RNAi and two overexpression lines for paralytic) were available. Not considering combinations with phenotypes that could result from an additive effect (Supplementary Table S3), we observed eight situations (out of 12) in which the simultaneous dosage modification of paralytic resulted in an altered bang or wing phenotype compared to overexpression of RhoBTB alone (six antagonistic and two synergistic) (Supplementary Table S3, Fig. 1C–E). While simultaneous knockdown of paralytic resulted in an ameliorated phenotype and overexpression in a more severe RhoBTB wing phenotype, in the case of pan-neuronal modification, both simultaneous overexpression or knockdown of paralytic resulted in an ameliorated RhoBTB bang phenotype. These results do not allow to establish an exact regulatory mechanism, however, it is possible that the dosage alteration of paralytic itself might be more relevant than the direction of this alteration. This would be compatible with the observation that both loss-of-function and gain-of-function variants in SCN1A or other sodium channels result in seizure phenotypes [29, 31, 34]. In summary, our findings indicate a high evidence for genetic interaction between RhoBTB and paralytic.

Combined dosage manipulation of RhoBTB and ir76a was assessed in three combinations (two RNAi and one overexpression line for ir76a) in three assays, each. Out of the nine tested situations, phenotypic modification of the RhoBTB overexpression associated phenotype was observed in five, four times antagonistic and once synergistic (Fig. 1C, Supplementary Table S3, Supplementary Fig. S4). These findings indicate genetic interaction between RhoBTB and ir76a with moderate evidence as alterations could only be observed for one combination in two assays.

Combined dosage manipulation of RhoBTB and slowpoke was assessed in three combinations (two RNAi lines targeting slowpoke and one hKCNMA1 overexpression line) in three assays, each. Eight combinations were viable and could be tested. Phenotypic modification of the RhoBTB overexpression phenotype beyond an additive effect was observed in two situations. Simultaneous knockdown of slowpoke or overexpression of hKCNMA1 resulted in a milder geotaxis phenotype compared to overexpression of RhoBTB alone (Fig. 1C, Supplementary Table S3, Supplementary Fig. S4). For this pairwise combination, we obtained low evidence for genetic interaction.

Our genetic interaction screen therefore indicates a functional link between RhoBTB and the sodium channel paralytic and possibly the glutamate receptor encoding ir76a in vivo.

Establishing human neuronal models for RHOBTB2 related disorders

To follow up our observations in a human neuronal model, we induced various representative RHOBTB2 aberrations (Fig. 2A) in a standard hiPSC line from a healthy donor and differentiated them to neuronal progenitor cells (NPCs) and neurons (Fig. 2B). As overexpression of RhoBTB in flies was rather used as a proxy due to limitations in introducing specific variants, manipulating hiPSC lines by CRISPR/Cas9 allows a better reflection of the patient situation. We introduced three different missense variants from the GTPase domain and three different missense variants from the BTB domain region (Fig. 2A) in a heterozygous state to model the situation of de novo missense variants associated with either a severe DEE (BTB domain region) or a more variable NDD with or without epilepsy (GTPase domain) [19]. Based on the clear genotype–phenotype correlation regarding variant location and clinical manifestation, missense variants in the same domain region were grouped and analyzed together. Additionally, we induced homozygous knockout to reflect the situation in individuals with bi-allelic truncating variants reported previously as causing a variable neurodevelopmental and seizure phenotype [19] (Fig. 2A).

Figure 2

HiPSC-derived neuronal model for RHOBTB2-related disorders. (A) Schematic overview of RHOBTB2 (GenBank: NM_001160036.2) with domains and the selected missense variants located either in the GTPase domain or BTB domain region above. Domains were identified and re-colored based on SMART prediction [38]. Variants resulting in RHOBTB2 knockout are indicated below the scheme. (B) Schematic overview on generation and differentiation of hiPSC-derived neurons with either complete loss of RHOBTB2 or harboring specific heterozygous missense variants. Schema was generated with BioRender.com (C) exemplary immunofluorescence stainings of hiPSC with NANOG (green) and OCT4 (red). Images were taken with the AxioImager Z2 with the Apotome 3 with a 20× objective. Scale bar 100 μm. (D) Exemplary immunofluorescence stainings of NPC with PAX6 (green) and NES (red). Images were taken with the AxioImager Z2 with the Apotome 3 with a 40× objective. Scale bar 20 μm. (E) All NPCs showed a similar PAX6 positivity rate, calculated as the quotient of PAX6 positive cells/DAPI positive cells. (F) Exemplary immunofluorescence stainings of neurons after 5 weeks of differentiation with MAP2 (red) and TUBB3 (green). Images were taken with the AxioImager Z2 with the Apotome 3 with a 40× objective. Scale bar 20 μm. Images of all used cell lines can be found in Supplementary Fig. S7.

Presence of variants was confirmed by Sanger sequencing (Supplementary Fig. S5), and integrity of hiPSC lines was confirmed by staining with NANOG and OCT4 (Fig. 2C, Supplementary Figs S6, S7, and S8). In case of knockout, reduction of RHOBTB2 was confirmed by quantitative RT-PCR (Supplementary Fig. S6B) but could not be analyzed on protein level due to lack of an antibody recognizing endogenous RHOBTB2. After differentiation to NPCs, variants were confirmed by sequencing again (data not shown), and predicted off-target variants were excluded (Supplementary Table S4). Differentiation to NPCs was confirmed with specific markers NES and PAX6 (Fig. 2D, Supplementary Figs S6, S7, and S8), and all lines showed a similar PAX6 positivity rate of around 75% (Fig. 2E, Supplementary Figs S6, S7, and S8). Further differentiation to neurons at week five was confirmed by staining with specific markers MAP2 and TUBB3 (Fig. 2F, Supplementary Figs S6, S7, and S8) Of note, we could not observe any differences in expression of different ion channel genes whose orthologues have been found to be deregulated in the fly transcriptomes, in ten weeks old hiPSC-derived neurons (Supplementary Fig. S9). However, it has to be considered that fly transcriptomes were performed upon RhoBTB overexpression, while expression analysis in hiPSC-derived neurons was performed upon introducing specific missense variants.

Different variant classes in RHOBTB2 result in a different electrophysiological phenotype

To test if the reported variants result in altered electrophysiological properties, whole-cell patch clamp recordings were performed on 10 weeks old, hiPSC-derived neurons (Fig. 3A). We tested three lines, each, with heterozygous missense variants either located in the GTPase domain (Arg116Cys, Arg154Gln, Arg183Met) or in the BTB domain region (Arg483His, Arg507Cys, Arg511Gln) plus four different wildtype clones from the same isogenic control and two homozygous knockout lines plus three different wildtype clones from the same isogenic control.

Figure 3

Electrophysiological characteristics of mature hiPSC-derived neurons. (A) Exemplary hiPSC-derived neuron after 10 weeks of differentiation and its action potential (AP) responses to increasing current injections. Image was taken on a Nikon Ts2-FL microscope with an Axiocam 305 color camera (Zeiss). Scale bar 20 μm. (B) Neurons exhibited either non-repetitive (1 AP, small peaks) or repetitive (some AP, multiple AP) firing. Example traces shown for neurons responding to 35 pA current injections. All tested neurons (WT or harboring either domain specific missense variants or frameshifting variants) presented with repetitive/mature firing patterns, no differences between the groups were observed. (C) Neurons with variants in the BTB domain region show increased AP firing frequency compared to the wildtype (WT). Asterisks indicate statistical significance (^*P < 0.05, ^**P < 0.01, ^***P < 0.001, Wilcoxon test with Bonferroni correction). Traces of individual cell lines can be found in Supplementary Fig. S11. (D) Example responses to 25 pA or 100 pA current injections in the WT, GTP-variants and BTB-variant groups. (E–I) Characterization of AP and membrane electrophysiological properties obtained from a single AP (E–H) or multiple AP (I). While neurons with variants in the GTPase domain do not show any differences compared to the wildtype, neurons with variants in the BTB domain region present with decreased depolarization speed (G) and increased AP half width (F). Asterisks indicate statistical significance (^*P < 0.05, ^**P < 0.01, ^***P < 0.001, Wilcoxon test with Bonferroni correction).

Neurons from all lines fired with regular and repetitive action potentials (AP) upon current injections with comparable rheobase between wildtype and genetically altered cells. In all tested lines, more than 75% of cells possessed a mature firing pattern without significant differences in the distribution of firing pattern categories (Fig. 3B and Supplementary Figs S10A and S11A), therefore indicating similar states of maturation.

When assessing the firing frequency upon increasing current injections, we observed no significantly altered pattern for the homozygous knockout cells compared to the wildtype (Supplementary Fig. S10). Current clamp recordings from cells with heterozygous missense variants in the GTPase domain resulted in a similar firing frequency as the wildtype, however, with consistently but not significantly increased firing responses at higher current injections (Fig. 3C and D). In contrast, missense variants in the BTB domain region resulted in a different pattern with increased firing frequency at all current injections (Fig. 3C and D) as well as increased AP half width, decreased depolarization speed (Fig. 3E–I). This indicates that neurons with missense variants in the BTB domain region react to current injections with enhanced excitability and increased capacity to fire action potentials. The altered half width and depolarization speed suggested alterations in voltage-gated sodium or potassium channel composition.

Overall, neurons with variants in the BTB domain region display an altered electrophysiological behavior compared to wildtype neurons, while both knockout neurons as well as neurons with missense variants in the GTPase domain show no significantly altered response.

Discussion

Only recently, we and others identified different variants in ubiquitin ligase complex subunit RHOBTB2 to cause variable epilepsy and other neurological and neurodevelopmental phenotypes [7–20, 39]. While de novo missense variants in the BTB domain region are associated with a distinct, severe DEE, de novo missense variants in the GTPase domain or bi-allelic truncating variants are associated with milder and more variable neurodevelopmental and seizure phenotypes [19]. This adds RHOBTB2 to the growing list of UPS-related genes implicated in DEEs and NDDs. However, little is known about the underlying pathomechanisms, neither for RHOBTB2-related phenotypes in general, nor for the observed genotype–phenotype correlations. We now found several lines of evidence that deregulation of ion channels might play a role in the pathogenesis of seizures in RHOBTB2-related DEE due to de novo missense variants clustering in the BTB domain region, while at least in vitro such a role is not obvious for heterozygous missense variants in the GTPase domain or complete loss.

First, transcriptome analysis in flies overexpressing RhoBTB, which models the RHOBTB2 abundance associated with BTB domain variants, revealed enrichment of ion channel genes among the differentially expressed genes. One of the downregulated genes, paralytic, encodes the only voltage gated sodium channel α-subunit in Drosophila. A gain-of-function mutation in paralytic represents one of the most acknowledged fly models for seizure susceptibility [40, 41]. Interestingly, paralytic is the fly ortholog of several epilepsy related genes in humans, including SCN1A, SCN2A, SCN3A and SCN8A. Variants in SCN1A are a relatively frequent cause of various epilepsy-related and other neurodevelopmental and neurological disorders (including DEE6B, non-Dravet: MIM#619317, Dravet syndrome: MIM#607208, Febrile seizures, familial, 3A and Generalized epilepsy with febrile seizures plus, type 2: MIM#60443, Migraine, familial hemiplegic, 3: MIM#609634, early-onset DEE associated with gain-of-function variants [42]). SCN2A is associated with a developmental and epileptic encephalopathy (DEE11 (MIM#613721)) as well as episodic ataxia, type 9 (MIM#618924) and benign familial infantile seizures type 3 (MIM#607745) with a variable phenotypic spectrum [43]. Variants in SCN3A are associated with a developmental and epileptic encephalopathy (DEE62, MIM#617938) and with familial focal epilepsy with variable foci 4 (MIM#617935) and other neurodevelopmental disorders typically including epilepsy [44]. Variants in SCN8A are associated with different neurodevelopmental disorders [45, 46] (cognitive impairment with or without cerebellar ataxia: MIM#614306) and seizure related disorders (DEE13: MIM#614558, Seizures, benign familial infantile, 5: MIM#617080).

For all four of these sodium channels, both pathogenic variants with gain-of-function effects as well as variants resulting in loss of function have been reported, resulting in various neurodevelopmental and/or neurological phenotypes, but often with seizures as an overlapping aspect [29, 31, 34]. This indicates that in most cases more than a single pathomechanism can contribute to altered neuronal excitability and uncontrolled discharge presenting as seizures.

Second, we validated a functional link between RhoBTB and paralytic in vivo by performing genetic interaction experiments in Drosophila. While knockdown or overexpression of paralytic alone had no effect on bang sensitivity, pairwise combination with RhoBTB overexpression resulted in modulation of the RhoBTB-overexpression associated phenotypes. This suggests that the family of voltage gated sodium channels might play a role in the pathomechanisms of RHOBTB2-related DEE. As also other classes of ion channels were among the deregulated genes in fly transcriptomes they might have an additional contributory role. Similar genetic interaction experiments with ionotropic-receptor ir76a also indicated a potential functional link with RhoBTB. Ir76a belongs to a specific subgroup of ionotropic glutamate receptors (iGluRs) which play an important role in interneuron communication [47, 48]. In humans, many genes from this family of (ionotropic) glutamate receptors, including the subfamilies GRIA, GRID, GRIK and GRIN have been associated with various neurological and neurodevelopmental disorders, including spinocerebellar ataxia, intellectual disability and epilepsy [33, 35, 36, 49–52]. Genetic interaction was less evident for potassium channel encoding slowpoke. Variants in the human ortholog KCNMA1 have been associated with different forms of epilepsy and seizures (EIG16: MIM#618596), developmental delay (Cerebellar atrophy, developmental delay, and seizures: MIM#617643) and neurologic dysfunction, dyskinesia and polymalformation (Liang-Wang syndrome: MIM#618729, Paroxysmal nonkinesigenic dyskinesia, 3, with or without generalized epilepsy: MIM#609446).

Third, we could confirm the role of ion channels in a neuronal model of RHOBTB2-related DEE. HiPSC-derived neurons carrying heterozygous missense variants in the BTB domain region presented with higher excitability. Whether the observed effect is based on a reduced number of ion channels (as suggested by fly transcriptomes) or an altered channel activity or opening cannot be concluded from our analysis.

Our electrophysiological observations confirm the specific effect of BTB domain variants, previously reflected in a more severe and more specific clinical phenotype in affected individuals and different in vivo molecular behavior [7, 19], Distinct, yet unidentified pathomechanisms might underlie the milder and more variable neurodevelopmental and seizure phenotypes caused by de novo missense variants in the GTPase domain or by bi-allelic truncating variants since neurons carrying these variant types did not show significantly electrophysiological alterations.

Such genotype–phenotype correlations and underlying different pathomechanisms challenge not only diagnostics, variant interpretation and prognosis, but also the understanding of specific disease mechanisms that might require particular therapeutic approaches. Variant specific therapies e.g. based on loss- or gain-of-function effects have been applied for several epilepsy genes already [53–55].

While we provide cumulative evidence that increased abundance of RhoBTB/RHOBTB2 either upon overexpression or in association with particular missense variants results in an ion channel deficiency or dysfunction, the exact nature of this link remains elusive. Therefore, we propose two possible scenarios: Firstly, RHOBTB2 might contribute to the transcriptional regulation of ion channel genes, as suggested by the transcriptome data from flies. As RHOBTB2 is not known to act as a transcriptional regulator itself due to lack of any obvious DNA-binding motif [56, 57], an indirect mechanism has to be postulated, possibly via deregulation of one or more of its substrates which in turn might be direct regulators of ion channel expression. In accordance, several proteins with BTB domains have been shown to act as transcription regulators by recruiting co-repressor complexes [58]. Related to its proposed role as a tumor suppressor, altered expression of RHOBTB2 has also been associated with (indirect) downregulation of genes involved in membrane trafficking and actin cytoskeleton organization [23]. However, as altered expression of ion channel encoding genes was only observed upon overexpression of RhoBTB in flies, but not in hiPSC-derived neurons carrying specific missense variants, incongruent mechanisms between organisms and/or dosage and variant situation have to be considered. As discussed previously [7, 19], missense variants in RHOBTB2 seem to result in rather specific effects without being clearly classifiable as gain of function, loss of function or neomorphic.

Alternatively, there may be a more direct effect of the ubiquitin proteasome system on the amount of ion channels, e.g. by ubiquitination and/or deubiquitination [59, 60]. For example, different voltage-gated ion channels have been shown to be a direct substrate of the E3 ubiquitin ligase Nedd4-2 in excitable cells [61, 62]. So far, no direct interaction between RHOBTB2 and ion channels has been documented (http://www.thebiogrid.org), and Musashi-2 is the only identified direct substrate of RHOBTB2 to date [63].

To conclude, our findings shed further light into RHOBTB2-related disorders by supporting a role of deregulated ion channels in the pathomechanism of RHOBTB2-related DEE. Furthermore, they confirm distinct mechanisms underlying the observed genotype–phenotype correlations regarding zygosity, location and nature of RHOBTB2 variants.

Materials and methods

RNA isolation from flies or cells

For RNA extraction from flies, whole larvae, adult flies, wings or heads were collected and frozen at −80°C. Fly tissues were shredded (QIAshredder columns, Qiagen), and RNA was extracted according to the RNeasy Lipid Tissue Mini Kit (Qiagen) with substitution of QIAzol by TRIzol and DNAse digestion to remove DNA contaminations. RNA was isolated from cultured (human) neurons after ten weeks of differentiation. Cells were harvested in RLT buffer, and RNA was then isolated with the RNeasy Plus Mini Kit (Qiagen) with an additional DNAse digest step to remove potential DNA remains.

RNA was reversely transcribed into cDNA using the SuperScript II (ThermoFisher) with random hexamer primers.

Transcriptome sequencing and analysis

RNA quality was checked using the Bioanalyzer 2100 system (Total RNA Analysis ng sensitivity (Eukaryote) Kit, Agilent) and confirmed to be of high quality and not degraded. Library preparation for RNA-sequencing was done using the TruSeq Stranded mRNA kit (Illumina), and sequencing was performed on a HiSeq2500 platform (Illumina) with single-end sequencing (101 bp). Conversion of reads to.fastq format and demultiplexing was performed with bcl2fastq v2.17. For adapter trimming and quality filtering cutadapt v1.18 [64] was used, and quality control was done with fastqc v0.11.8. Mapping of processed reads to the Drosophila reference genome BDGP6.11 with the Ensembl Annotation v98 was done using STAR aligner v.2.6.1c [65]. Read counting of mapped reads (number of reads per gene = sum of all exons) was performed using samtools v1.8 and subreadv1.6.1 [66]. The following analyses were performed using R v3.6.1. Differential gene expression analysis was performed using DESeq2 v1.24.0 [67] with the apeglm package [68] for noise removal. For further data analysis, only annotated protein-coding genes were considered.

Gene ontology enrichment analysis was performed in Panther [24, 25] with the slim GO annotation and with a Fisher’s exact test and FDR (false discovery rate) correction for multiple testing. As a background, a list of all protein-coding genes for which expression was detected in RNA-Seq was used. For comparison with human genes, Drosophila genes were mapped to its human orthologs using DIOPT ortholog finder retaining only one hit (best score) per gene [69, 70]. For annotation of human disease genes, the SysNDD database (https://sysndd.dbmr.unibe.ch/) was used [30].

Quantitative RT-PCR for expression analysis

Quantitative RT-PCR was performed to confirm knockdown or overexpression of genes of interest in the flies as well as RHOBTB2 expression in hiPSCs and SCN1A, SCN2A, KCNMA1, GRID2 and MAP2 expression in human neurons. Expression analysis was performed on a QuantStudio 3 or a QUantStudio 12 K Flex (Life Technologies) with exon spanning primers for the different genes (Supplementary Table S5) as well as the corresponding endogenous controls (dmTubulin (Drosophila), hB2m (hiPSC) or hACTB, hGAPDH, hRPL13a, hC1orf43 (neurons)) and the PowerUP SYBR Green Master Mix (Thermo Scientific). To compare the expression levels, RQ values were calculated based on the ∆∆Ct method and quadruplicates or triplicates.

Drosophila lines and conditions

All flies were kept on standard food containing sugar, corn flour and yeast at room temperature. RNAi-mediated knockdown or overexpression was achieved using the Gal4-UAS system [71], crosses were performed at 28°C and confirmed by quantitative RT-PCR (method see above, Supplementary Fig. S2). If lethality was observed at this temperature, crosses were performed at room temperature [72].

Tissue specific driver lines (elav (pan-neuronal), actin (ubiquitous), ms1096 (wing)), RNAi lines (for paralytic, ir76a, slowpoke) and transgenic lines (for ir76a) as well as the matching controls were obtained from the Bloomington Stock Center, the Vienna Drosophila Resource center (VDRC) or have been created previously (UAS-RhoBTB_1, UAS-RhoBTB_2) [7]. For slowpoke, only an overexpression line of the human orthologe KCNMA1 was available from the Bloomington Stock Center. Transgenic paralytic lines for overexpression were a gift from Beverly Pigott, University of California at San Francisco, and the line containing an UAS element with c.elegans cDNA and the double balancer line Kr/CyO;D/Tm6CSbTb were a gift from Annette Schenck, Radboud University, Nijmegen; the elav-Gal4/CyO;MKRS/Tm6BTb was a gift from Beat Suter, University Bern.

Three double balancer lines (Kr/CyO;D/Tm6CSbTb, elav-Gal4/CyO;MKRS/Tm6BTb and FM6(B)/Y;Sco/CyO) were used to generate double transgenic fly lines with either stable pan-neuronal (elav-Gal4/CyO;UAS-RhoBTB_2/(Tm6BTb) or wing specific (ms1096_Gal4/(FM6B);UAS-RhoBTB_1/(CyO)) overexpression of RhoBTB. Additionally, a fly line with stable wing specific overexpression of an UAS-element with scrambled cDNA from c.elegans was created to serve as a control for a potential distribution effect of the Gal4 element between two UAS-elements. These lines were then crossed with fly lines inducing either knockdown or overexpression of the potential interactor genes. All lines except paralytic_KD 1 and 2 (RT) were crossed at 28°C.

A list of all used fly lines and their genotype is provided in Supplementary Table S6. In the text and figures, RNAi-lines are referred to as gene_KD, and transgenic overexpression lines as gene_OE.

Geotaxis and bang sensitivity assay

Climbing and bang sensitivity assays upon pan-neuronal manipulation with the Gal4-UAS system were performed as described previously [7]. For the bang sensitivity assay, the flies were vortexed for 10 s and filmed for 2 min. The fraction of flies that were paralyzed or having spasms were determined every 10 s for 2 min. A minimum of 3 × 10 flies (average 20 × 10) were tested per condition. In case of male lethality only female flies were analyzed, and the control was adapted accordingly. Analysis of both assays was performed blindfolded.

Assessment of wing morphology

To analyze wing morphology, total numbers of male and female flies, respectively, were counted, and morphological aberrations such as shape, size, curling, bristles as well as the presence of additional wing veins and the number of intact cross wing veins (a-cv and p-cv) were assessed. For each fly, the two wings were analyzed independently. Pictures were taken on a Zeiss Discovery V8 stereo microscope. P-values were calculated using the chi-squared test followed by Bonferroni correction.

Genetic interaction screens

For the genetic interaction studies, RhoBTB and a potential interactor were modified simultaneously using the Gal4/UAS system. Flies with RhoBTB overexpression and either RNAi induced knockdown or overexpression of the potential interactor were analyzed. If possible, multiple RNAi or overexpression lines were used per interactor, resulting in different numbers of combinations (Supplementary Table S3). Phenotypic modifications/alterations were assessed in the geotaxis and bang sensitivity assay and regarding wing morphology. Phenotypic modifications were categorized as antagonistic (phenotype milder), additive (combination of the single phenotypes) or synergistic (phenotype more severe than just addition of the single phenotypes). As controls, flies with RhoBTB overexpression alone as well as with overexpression of scrambled cDNA were used.

Experiments were performed in two rounds, where first wildtype flies and flies with modification of the interactor alone were compared, and then flies with RhoBTB overexpression alone and together with manipulation of the interactor. P-values were calculated with the Wilcoxon-Mann–Whitney test, followed by Bonferroni correction for the two rounds separately. For the representing figures, all four conditions were displayed in one figure.

We categorized the potential interaction into different groups regarding the evidence for a functional link between the tested gene and RhoBTB: High evidence, when alteration of the RhoBTB associated phenotype was observed in at least two out of three assays in two combinations each. Moderate evidence, when alteration of the RhoBTB associated phenotype was observed in two out of three assays in one combination. Low evidence, when alteration of the RhoBTB associated phenotype was observed in one assay for one combination.

Cell lines and maintenance

All cells were cultured at 37°C and 5% CO₂ under sterile conditions and passed according to their cell type. Human induced pluripotent stem cells (hiPSCs) were obtained from Thermo Scientific (A-18944) and are derived from cord blood derived CD34+ lymphocytes of a healthy individual and cultured in mTeSR Plus (STEMCELL Technologies) with 1% penicillin/streptomycin (ThermoFisher). To passage them, cells were treated with Gentle Cell Dissociation Reagent (STEMCELL Technologies) for 5 min, resuspended in medium supplemented with 1:1000 Rock Inhibitor (RI, Y-27632, Lucerna-Chem) and seeded on plates coated with Geltrex (ThermoFisher).

Genome editing with CRISPR/Cas9

To generate RHOBTB2 knockout cells or to introduce point variants with the CRISPR/Cas9 system, different guide RNAs targeting the specific region were designed using CRISPOR (http://crispor.tefor.net/) and cloned into the pX330-GFP plasmid. The efficiency of the guides was confirmed with the T7 assay in HEK293 cells. Single stranded donor templates were designed to introduce the specific point mutations as well as mutations to prevent recutting of the position (https://eu.idtdna.com/pages). A full list of the different guide RNAs and templates is provided in Supplementary Table S4.

For the genome editing, two wells of a 6 well plate of confluent hiPSCs were dissociated with Accutase (Millipore) for 10 min, washed in DMEM/F12 and pooled. Cells were nucleofected with 2.5 μg plasmid DNA (pX330 containing the different sgRNA, Supplementary Table S4) and in case for the missense variants additionally the matching template (IDT) (Supplementary Table S4) and alt-R HDR Enhancer V2 (IDT) using the Nucleofector 2b kit (program B16) and device (Lonza) and transferred to one well of a 12 well plate containing media and RI. After 72 h, cells were dissociated with Accutase, washed and resuspended in mTeSR Plus supplemented with 1:10 CloneR (STEMCELL Technologies). Single GFP positive cells were sorted into 96 well plates, half media changes were performed until cells reached 70% confluency in the well. Cells were dissociated with Gentle Cell Dissociation Reagent and transferred to two 48 well plates for expansion and for genotyping. When confluent, cells for expansion were dissociated with Gentle Cell Dissociation Reagent and frozen in 100 μl Bambanker freezing medium (LubioScience).

Validation of hiPSC lines

For genotyping, DNA was extracted by dissociating cells with Gentle Cell Dissociation Reagent and resuspending in 50 mM NaOH. After incubation at 90°C for 20 min, 1/10 1 M Tris–HCl was added and 2 μl were used for genotyping PCR (primers: Supplementary Table S5). Here, the edited region was amplified using specific primers and analyzed using Sanger sequencing. Additionally, the clones were screened for potential exonic off target effects if the guide RNA sequence was predicted to bind at additional regions with three mismatches. All used primers are listed in Supplementary Table S7.

Differentiation of hiPSC into NPCs and neurons

Differentiation of hIPSC to NPCs was performed using a modified version of the STEMdiff SMADi neural induction kit (STEMCELL Technologies).

To initiate differentiation of hiPSCs into NPCs, the medium was replaced with DMEM/F12 supplemented with 1× N2 (ThermoFisher), 1 μM Dorsomorphin (Tocris), 2 μM A-83-01 (Tocris) and 1% penicillin/streptomycin. After 24 h, cells were dissociated with Accutase, resuspended in STEMdiff Neural Induction medium (STEMCELL Technologies) supplemented with RI and seeded into one AggreWell800 (STEMCELL Technologies) (3.0 × 10⁶ cells). After five days, the neuronal aggregates were harvested and plated on poly-ornithine (Sigma-Aldrich) and laminin (Sigma-Aldrich) (POL) coated plates and cultured for seven days. After that, the neural progenitor cells were isolated using the STEMdiff Rosette selection reagent (STEMCELL Technologies) and plated onto POL coated plates and cultured in STEMdiff Neural Induction medium until they were confluent. Afterwards, medium was changed to DMEM/F12 supplemented with 0.5× N2 (ThermoFisher), 0.5× B27 (ThermoFisher), 20 ng/μl FGF2 (Peprotech, ThermoFisher) and 1% penicillin/streptomycin. Before passing or freezing, cells were treated with Accutase (Millipore) for 10 min, resuspended in media and seeded to POL coated plates.

To further differentiate NPCs into neurons, the medium was changed to STEMdiff Forebrain Neuron Differentiation medium (STEMCELL Technologies). After one week, the cells were dissociated with Accutase, washed and resuspended in STEMdiff Forebrain Neuron Maturation Medium and counted. For immunofluorescence staining or patch clamp recordings, cells were then seeded onto POL coated coverslips at a density of 4.5 × 10⁵ cells per well in a 24 well plate. Neurons were cultivated for at least 10 weeks with medium changes every 3 days.

Immunofluorescence

For immunofluorescence, cells were grown on coated coverslips (NPC or neurons: POL, hiPSC: Geltrex) and fixated with 4% paraformaldehyde in PBS. For the staining, primary antibodies against NANOG, PAX6, NES, MAP2, TUBB3 and the corresponding secondary antibodies as well as DAPI (1:50000, Serva) were used. All antibodies used are listed in Supplementary Table S8. Images were taken with a Zeiss Axio Imager Z2 microscope (Zeiss) with a 20×, 40× or 63× lens and Apotome 3 and the Zen software v3.4 and analyzed using ImageJ (v1.53) [73] and CellProfiler (v4.2.5) [74]. The ratio of PAX6 positive cells was quantified from at least five different pictures taken from random positions of NPCs using CellProfiler.

Electrophysiological recordings

For electrophysiological recordings, neurons were grown and differentiated on cover slips for 10 to 11 weeks. Cover slips were transferred into the recording chamber at room temperature and superfused with artificial cerebrospinal fluid (ACSF) containing 125 mM NaCl, 2.5 mM KCl, 25 mM NaHCO₃, 1.25 mM NaH₂PO₄, 1 mM MgCl₂, 2 mM CaCl₂ and 25 mM glucose, saturated with 95% O₂/5% CO₂ at a flow rate of 0.35 ml/min. Heat-pulled borosilicate glass pipettes with a tip resistance of 5–8 MΩ were used for whole cell patch clamp recordings. All recordings were performed with an intracellular solution containing 130 mM K-gluconate, 5 mM KCl, 10 mM Na-phosphocreatine, 4 mM Mg-ATP, 0.3 mM Na-GTP, and 10 mM HEPES. Signals were recorded with a patch-clamp amplifier (BCV-700A, Dagan), filtered at 5 kHz and digitized at 10 kHz (ITC-16, Instrutech).

In current clamp experiments, current was injected to hold neurons at −60 mV (±5 mV), and the input–output relationship between a series of 500 ms-long current steps (−25 pA—100 pA) and the action potential (AP) frequency response was evaluated as a measure of intrinsic neuronal excitability.

Only cells with repetitive firing were considered mature and used for the following analysis. The minimum current required to elicit a single action potential was extrapolated from the input–output curve and expressed as the AP rheobase. AP shape, depolarization speed and repolarization speed were analysed from individual APs elicited by a single 5 ms, 600 pA current step.

In a subset of experiments, the contribution of voltage gated sodium currents was analyzed by bath-applying the Na + channel blocker tetrodotoxin (TTX, 500 nM) for 10 min prior to measurements.

All data acquisition and analyses were performed using custom-made procedures in Igor Pro (WaveMetrics) and MatLab (MathWorks) software.

Statistics

A Wilcoxon signed rank test was performed for the climbing and bang sensitivity assay in Drosophila as well as most electrophysiological recordings except or distribution of different respond types. For the wing morphology phenotype as well as the distribution of different respond types in the electrophysiological recordings, the qui-squared t-test was used. In all cases, Bonferroni correction for multiple testing was applied.

Acknowledgements

RNA-Seq was performed at the NGS core unit of the Medical Faculty of the FAU Erlangen-Nürnberg. We thank staff of the IMSV Bern statistical counseling. We thank Beverly Pigott, Annette Schenk and Beat Suter for providing transgenic fly lines. Fly stocks obtained from the Bloomington Drosophila Stock Center (NIH P40OD018537) and the Vienna Drosophila Research Center (VDRC, www.vdrc.at) were used in this study.

Conflict of interest statement: The authors report no conflict of interest.

Funding

C.Z. has been supported by the German Research Foundation (DFG, ZW184/6-1, 270949263/GRK2162), the Novartis Foundation for medical-biological research (#22C194) and the IZKF Erlangen (E31).

Data availability

All data that support the findings of this study and are not included in manuscript or the supplement are available from the corresponding author upon request. All materials generated in this study (cell lines, vectors, fly lines) are available from the corresponding author upon request with a completed material transfer agreement.

The datasets generated during and/or analysed during the current study are available from the corresponding author on reasonable request.