Biallelic missense variants in ZBTB11 can cause intellectual disability in humans

Abstract

Exploring genes and pathways underlying intellectual disability (ID) provides insight into brain development and function, clarifying the complex puzzle of how cognition develops. As part of ongoing systematic studies to identify candidate ID genes, linkage analysis and next-generation sequencing revealed Zinc Finger and BTB Domain Containing 11 (ZBTB11) as a novel candidate ID gene. ZBTB11 encodes a little-studied transcription regulator, and the two identified missense variants in this study are predicted to disrupt canonical Zn²⁺-binding residues of its C2H2 zinc finger domain, leading to possible altered DNA binding. Using HEK293T cells transfected with wild-type and mutant GFP-ZBTB11 constructs, we found the ZBTB11 mutants being excluded from the nucleolus, where the wild-type recombinant protein is predominantly localized. Pathway analysis applied to ChIP-seq data deposited in the ENCODE database supports the localization of ZBTB11 in nucleoli, highlighting associated pathways such as ribosomal RNA synthesis, ribosomal assembly, RNA modification and stress sensing, and provides a direct link between subcellular ZBTB11 location and its function. Furthermore, given the report of prominent brain and spinal cord degeneration in a zebrafish Zbtb11 mutant, we investigated ZBTB11-ortholog knockdown in Drosophila melanogaster brain by targeting RNAi using the UAS/Gal4 system. The observed approximate reduction to a third of the mushroom body size—possibly through neuronal reduction or degeneration—may affect neuronal circuits in the brain that are required for adaptive behavior, specifying the role of this gene in the nervous system. In conclusion, we report two ID families segregating ZBTB11 biallelic mutations disrupting Zn²⁺-binding motifs and provide functional evidence linking ZBTB11 dysfunction to this phenotype.

Introduction

The high levels of cognition and intelligence that have evolved in humans have provided us with extraordinary abilities to adapt, survive and thrive. It has been proposed that genes predispose 50–70% of cognition variation among people. Although high heritability is reported for intelligence, gene–environment correlations cannot be overlooked, and intelligence is considered as a complex polygenic human trait (1,2). Understanding the biological and molecular pathways underpinning the processes of learning and memory can solve the complex puzzle of how cognition develops. One of the common approaches to study the molecular basis of human cognition is to study genetic mutations in cognitively impaired patients (3).

Intellectual disability (ID), defined as significant deficit in intellectual functioning and adaptive behavior, is a common condition that affects 1–3% of people worldwide, and genetic factors are considered as one of the main contributors and to moderate and severe forms in particular (IQ <50) (4). In recent years, due to the advent of new technologies such as next-generation sequencing (NGS), large-scale systematic studies are being performed to identify genes involved in the pathology of ID, and over 700 genes have been identified for different types of this heterogeneous disorder. Although we may be approaching the identification of most or all X-linked ID genes, there are still many autosomal genes yet to be discovered (4). Candidate gene identification is undoubtedly valuable in diagnosis and appropriate genetic counseling, but it should be noted that an etiopathological connection has yet to be established for many of the candidate ID genes identified through whole exome/genome sequencing approaches. In addition, it has been shown that genetic variation underlying polygenic traits such as intelligence are enriched in genes responsible for the related mendelian disorders. So, the identification of new genes for ID introduces candidate genes that may also harbor genetic variation impacting intelligence within the non-ID population (1).

In the post-NGS era, in vivo/in vitro functional systems for validating candidate genes are rapidly being developed. These new approaches are aimed at clarifying the specific molecular effects of a candidate mutation, along with determining the biological pathways and networks in which the candidate gene is involved (4).

These functional approaches will be beneficial at the clinical level, as they attempt to authenticate the candidate genes being used in prenatal diagnosis. Moreover, they can determine major pathways and networks involved in brain development and cognitive functioning. The latter may lead to the identification of potential therapeutic targets or pinpoint where gene therapies may be applied. One of the main pathways for ID genes that has emerged over the recent years is regulation of gene transcription(4). Transcription, as a critical process for memory consolidation, is required for long-term synaptic plasticity and memory. It acts by regulating the expression of specific genes and synthesis of proteins required in activated neurons. So far, many transcription-related genes are being identified as candidate ID genes presenting new challenges to associate these candidates into their possible neuronal targets (5). Mutations in such genes can explain complex phenotypes with a broad range of symptoms in patients, as transcriptional regulators may themselves influence many different and diverse pathways (3). As part of two separate systematic studies to identify candidate ID genes, here we report likely pathogenic variants of Zinc Finger and BTB Domain Containing 11 (ZBTB11)—a little-studied transcriptional regulator—in two consanguineous ID families and report on the molecular and cellular consequences of the two missense mutations.

Results

Missense variants in highly conserved residues of ZBTB11 C2H2 zinc finger domains co-segregate in two consanguineous families with ID

Two consanguineous families featuring autosomal recessive ID were selected in this study based on the overlapping homozygosity by descent (HBD) regions and similar gene defects latterly. Family A (Iranian) was investigated in a study comprising 404 ID families with more than two affected individuals. The three affected children of this family were from a third-cousin marriage, presenting with moderate ID (IQ: 40–45), microcephaly (-2.8 S.D.), slightly delayed motor milestones, slight facial hypotonia and maxillary hyperplasia, ataxic gait, spasticity, cerebellar atrophy and mild ventriculomegaly in cranial magnetic resonance imaging (MRI; Fig. 1; Supplementary Material, Table S1)(6). Family B (Pakistani) was identified from a study of 192 multiplex ID families and included affected individuals presenting moderate ID and variable clinical features within three branches in which the parents were first cousins (Fig. 1; Supplementary Material, Table S1) (7). Whole-genome single-nucleotide polymorphism (SNP) genotyping for six family members of Family A were used for homozygosity mapping and multipoint parametric linkage analysis revealed a more prominent homozygous interval (the largest) located between heterozygous SNP markers, rs9867935 and rs3849511. This interval is located on chr3: 86538834–102651023 encompassing 122 genes with calculated parametric Logarithm of odds (LOD) score of 2.529 and nonparametric Z-mean of 6.41. Genotype analysis for 14 individuals in Family B showed a single large homozygous region on chr3: 71687846–109144543 flanked by rs12492625 and rs797078 encompassing 247 genes (homozgyosity (HM) score: 1250) (7). Subsequent sequencing of probands’ DNA revealed two homozygous missense variants within the HBD intervals in ZBTB11 gene as follows: Chr3: 101374954G>A (NM_014415.3: c.2185C>T, p.H729Y) in Family A and Chr3: 101371344A>C (NM_014415.3: c.2640 T > G, p.H880Q) in Family B.

Figure 1.

(A) Pedigree of Family A (The filled symbols indicate affected individuals). Co-segregation of c.2185C>T is shown under pedigree and the homozygous intervals obtained by linkage analysis are shown in the right. (B) The MRI images and Facial appearance of patients from Family A and B. a. T1 and T2-weighted in cross sectional (axial) view showing ventriculomegaly in individual III-2; b. T1 and T2-weighted in cross sectional (axial) view showing mild ventriculomegaly in individual III-4; c. T1-weighted in sagittal view showing cerebellar hypoplasia in individual III-2, the white matter is normal; d. Facial appearance of Family A patients from left to right; individual III-2, individual III-3, individual III-4; all showing hypotonia and maxillary hyperplasia e. T1 and T2-weighted in cross sectional (axial) view showing no prominent defect in individual IV-7 from Family B. This patient is younger from other two affected individuals who showed defect in their MRI images. Also, the clinical presentations in this individual is significantly mild compared to the others which may explain the normal MRI results. e. Facial appearance of Family B patients from left to right; individual IV-14 (showing hypotonia and maxillary hyperplasia), individual IV-7(normal) (C) Pedigrees of Family B (Pakistani). Co-segregation of c.2640 T > G is shown under pedigree. The HBD region on chromosome 3 is shown in the right.

Moreover, to evaluate whether additional pathogenic variants might play a role in complex and variable clinical features of affected individuals in both families, whole-exome sequencing (WES) data were also analyzed beyond the homozygous intervals. This led to the identification of seven and two other homozygous rare variants in Family A and B, respectively (Supplementary Material, Table S2). Among these, all were de-prioritized based on nonstrongly-related gene functions to the clinical features presented in patients and only an LOF variant (NP_001161693.1: p.Arg1779Ter) in ALS2CR11 was prioritized in Family A, as it was located in an interval previously mapped for an autosomal recessive form of juvenile amyotrophic lateral sclerosis. However, this stop-gain variant is only present in two transcripts of ALS2CR11and is located in the last exon downstream of the final exon–junction complex, which reduces the possibility of nonsense-mediated decay. Furthermore, co-segregation studies revealed that Arg1779Ter genotype is not similar in the three affected children of the family (the two other affected individuals were WT for this variant) and one of the normal siblings showing none of the clinical features observed in individual III-2 is also homozygous for the this variant. So, it seems that the complex and variable clinical features observed in patients of both families can be considered as a whole entity caused by defects in the ZBTB11 gene.

The two variants are located within HBD intervals, predicted to be damaging by bioinformatics predictions and are absent from controls in public databases as well as 800 normal Iranian individuals (www.iranome.com) (Table 1). Positive z = 1.91 of missense variants in this gene suggests that ZBTB11 is intolerant to missense variants. Also, ZBTB11 is not considered as extremely intolerant to LOF variants (haploinsufficient) but the value (pLI = 0.27) supports the recessive (tolerated heterozygous LOFs) mode (inferred from ExAC Browser). In addition, no homozygous similar missense variants at Zn2+-binding residues, can be observed in gnomAD and ExAC databases for this gene. Co-segregation studies confirmed that both variants are inherited along with the observed phenotype (Fig. 1A and C) and introduced ZBTB11 as a novel candidate ID gene. The two variants are located on conserved C2H2 residues of this classical C2H2-type zinc finger (Znf) protein. Replacement of the conserved histidine by tyrosine and glutamine respectively may lead to failure to bind the Zn²⁺ ion due to changes in amino acid spatial orientation that were predicted by in silico 3D modeling (Fig. 2). This, in turn, is likely to alter the overall conformation of the respective Znf s, which may then impact the DNA-binding specificity of the protein.

Table 1.

Description of identified variants in ZBTB11 gene

Family number	Nationality	No. of affected individuals	Transcript	Variant	Logit score	SIFT	PolyPhen2	Mutation Taster	CADD Phred Score	ExAC database	gnomAD database	Iranome database
Family A	Iranian	3	NM_014415	Chr3: 101374954 G > A (c.2185C>T; p.H729Y)	6.199	Deleterious	Probably Damaging	Disease Causing	27.5	NA	NA	NA
Family B	Pakistani	6	NM_014415	Chr3: 101371344 A > C (c.2640T>G; p.H880Q)	NA	Deleterious	Probably Damaging	Disease Causing	24.7	NA	NA	NA

Figure 2.

In-silico three-dimensional (3D) models of the two missense variants are shown. (A) ZBTB11 WT C2H2-type zinc finger domain sequences are shown for the 5th and 10th zinc fingers (from residue 707 to 730 and from residue 858 to 881, respectively). The H729 and H880 residues are shown by arrow. WT (upper) and variant (lower) 3D in-silico models are shown for (B) H729Y and (C) H880Q, indicating the C2H2-type zinc fingers structure. WT residues H729 and H880 are shown versus the variant residues Y729 and Q880. The in-silico modeling shows their different 3D structure and spatial orientation, which may lead to dysfunctional zinc atom stabilization.

ZBTB11 is predominantly localized in nucleolus, and the two missense variants affect proper protein localization

We established an in vitro localization study to evaluate the pathogenicity of the observed missense variants. To investigate whether p.H729Y and p.H880Q variants disturb protein localization, HEK293T cells were transfected with WT and mutant Green Fluorescent Protein (GFP)-ZBTB11 constructs. Visualization of recombinant WT ZBTB11 protein showed that the intense GFP signals were accumulated in regions with low 4′,6-diamidino-2-phenylindole (DAPI) fluorescent signals predicted as nucleoli regions in the nucleus. The subsequent immunofluorescent staining (IF) of the transfected cells with a nucleolar marker protein, nucleolar phosphoprotein B23 (or nucleophosmin NPM), confirmed the nucleolar localization of the recombinant WT ZBTB11 protein, while concurrent immunofluorescent staining using anti-ZBTB11 antibody showed undetectable levels of endogenous ZBTB11 protein in nontransfected HEK293 cells. We also observed that the mutant fusion proteins of both missense variants failed to localize within the nucleolus, showing a scattered protein distribution in the nucleus rather than a concentrated clustering in nucleoli (Fig. 3A). Pearson correlation coefficients (PCCs) were generated for ZBTB11 and NPM colocalization and indicated significant alteration in the localization pattern of the mutant fusion ZBTB11 proteins in comparison to recombinant WT (one-way analysis of variance; P = 0.011) (Fig. 3B). It should be mentioned that transfected HEK293 cells using GFP-only constructs as a control in our experiments supports that mislocalization is due to the defect in ZBTB11 protein and not because of the GFP itself (data not shown). In addition, in a small number of cells transfected with H880Q mutant protein, localization in both nucleus and cytoplasm was observed at the same time (data not shown). So, we conclude that the two missense variants result in cellular mislocalization of the ZBTB11 protein, which is consistent with ZBTB11 protein dysfunction or loss of function.

Figure 3.

Cellular localization of WT recombinant ZBTB11 and the mutants (H729Y and H880Q) in HEK293T cells. (A) Confocal image stacks of immunofluorescent (IF) staining with nucleolar marker; nucleolar phosphoprotein B23 or nucleophosmin (NPM) is represented in this image. From left to right: DAPI (blue) in Column 1, GFP/FITC (green) in Column 2, Nucleolar phosphoprotein B23 (red) in Column 3 and the merged images are shown in Columns 4 and 5 for each of the WT and mutants. (B) Pearson correlation coefficients (PCCs) were generated for ZBTB11 and NPM and indicate significant alteration for the variant ZBTB11 in comparison to WT (one-way analysis of variance; P = 0.011). Bars represent mean ± standard error.

The enriched pathways and molecular/biological functions of ZBTB11 are associated with its localization in nucleolus, potentially linking this transcription factor to other known ID genes

As ZBTB11 is thought to be a transcription regulator, we applied Genomic Regions Enrichment of Annotations Tool (GREAT) (8) to investigate the enriched molecular/biological processes of available ZBTB11 ChIP-seq data deposited in the ENCODE database (9); ENCSR706BJO, ENCSR882ZTS, ENCSR155VDK and ENCSR331GDC. These data were obtained from three different cell lines, K562, HEK293 and MCF7, including ChIPed DNA derived from anti-ZBTB11 antibody and also from cells stably expressing eGFP-ZBTB11 fusion protein. The results showed that enriched pathways and molecular/biological functions of ZBTB11 targets can be categorized into main nucleoli-related functions such as ribosomal RNA (rRNA) synthesis, ribosomal assembly, RNA modification, protein translation and stress sensing, linking ZBTB11 function to its subcellular location identified in this study (Table 2).

Table 2.

Common enriched GO Cellular Component, molecular functions, biological process and MSigDB Pathways ascertained from pathway analysis of ENCODE ZBTB11-ChIP-seq data; eGFP-HEK293 (ENCSR882ZTS), eGFP-K562 (ENCSR706BJO), ZBTB11-K562 (ENCSR331GDC), ZBTB11-MCF-7 (ENCSR155VDK) applying Genomic Regions Enrichment of Annotations Tool (GREAT)

Nucleolus-related function	Ontology	Term	eGFP-HEK293, eGFP-K562, ZBTB11-K562, ZBTB11-MCF7
			Binom FDR Q-value	Binomial fold
RNA synthesis, Ribosomal assembly, RNA modification, mRNA biogenesis, RNA metabolism	GO Cellular Component	Small nuclear ribonucleoprotein complex	2.78e-10, NA, 9.00e-07, 1.10e-06	2.0765, NA, 3.0064, 5.1141
		Spliceosomal complex	NA, NA, 3.23e-16, 3.29e-10	NA, NA, 2.1841, 2.7557
		Mitochondrial ribosome	NA, NA, 4.74e-10, 4.37e-06	NA, NA, 3.1634, 4.213
	GO molecular Functions	Ribonucleoprotein complex binding and structural constituent of ribosome	5.54e-12, NA, 9.43e-12, 4.52e-06	2.0079, NA, 3.1805, 2.718
	GO Biological Process	rRNA processing	NA, NA, 2.28e-11, 1.06e-07	NA, NA, 2.3802, 3.1603
		rRNA metabolic process	NA, NA, 4.63e-11, 5.05e-07	2.3115, 2.9741
		ribosome biogenesis	NA, NA, 4.16e-17, 3.65e-10	NA, NA, 2.3625, 2.9671
		Nuclear-transcribed mRNA catabolic process, nonsense-mediated decay	2.56e-33, 2.18e-35, 1.01eE-13, 1.62e-06	2.2079, 2.1732, 2.2099, 2.529
	MSigDB Pathway	Genes involved in mRNA Splicing	1.89e-46, 3.56e-51, 3.91e-15, 3.03e-08	2.7969, 2.7713, 3.0311, 3.8061
		Genes involved in metabolism of mRNA	9.62e-45, 1.91e-59, 8.28e-16, 1.32e-12	2.0105, 2.1147, 2.1971, 3.0562
		Genes involved in processing of capped intron-containing pre-mRNA	2.62e-39, 1.41e-54, 5.48e-17, 4.36e-09	2.2917, 2.4646, 2.7374, 3.3486
		Genes involved in nonsense mediated decay enhanced by the exon junction complex	1.99e-30, 2.00Ee-37, NA, NA	2.29, 2.3604, NA, NA
		Spliceosome	9.09e-30, 2.07e-37, 1.08e-21, 1.00e-10	2.0913, 2.1647, 2.998, 3.6044
Protein translation	GO cellular component	Proteasome accessory complex	NA, NA, 7.10e-08, 1.84e-08	NA, NA, 4.0561, 7.7741
		Cytosolic large ribosomal subunit	6.58e-14, 5.16e-13, NA, NA	2.1632, 2.0372, NA, NA
	GO molecular Functions	Unfolded protein binding, translation initiation factor activity	6.00e-20, NA, 3.18e-05, NA	2.0753, NA, 2.2486, NA
	GO biological process	Translation (initiation, elongation and termination)	1.07e-26, 1.72e-30, 7.85e-15, 4.43e-09	2.3264, 2.3374, 2.4197, 2.158
		Protein targeting and localization to ER	NA, 1.07e-36, 2.11e-10, NA	NA, 2.1751, 2.3296, NA
	MSigDB pathway	Genes involved in translation	2.04e-22, 1.04e-49, 2.48e-14, NA	2.0368, 2.2345, 2.3862, NA
		Genes involved in peptide chain elongation	9.73e-26, 4.30e-32, NA, NA	2.3209, 2.3998, NA, NA
		Genes involved in unfolded protein response	NA, 2.16e-31, 4.08e-12, NA	NA, 2.4171, 2.8932, NA
		Genes involved in 3' -UTR-mediated translational regulation	2.04e-22, 7.34e-28, 1.96e-10, NA	2.0368, 2.0965, 2.4211, NA
Stress sensing	GO cellular component	Cytoplasmic stress granule	9.42e-17, 5.63e-16, NA, NA	2.4134, 2.2766, NA, NA
	GO molecular functions	NA	NA	NA
	GO biological process	Regulation of DNA-dependent transcription in response to stress, response to endoplasmic reticulum stress	NA, 1.20e-19, 2.19e-09, NA	NA, 2.0656, 2.0018, NA
	MSigDB pathway	NA	NA	NA

In order to associate ZBTB11 to its possible neuronal targets and also linking this gene to the network of ID-related genes, the final ChIP IDR peaks were confined based on their location on ID-related genes as described in methods. This strategy limited the possible ZBTB11 targets to a list of 34 ID genes (Supplementary Material, Table S3). Next, the presence of cardinal clinical features observed in Family A in all 34 ID gene phenotypes highlighted AP4M1, AAAS and SCO1 target genes among the others (Table 3). Interestingly, we found AP4M1, sharing all the 6 cardinal features in accordance with clinical manifestations of Family A. This highlights AP4M1 among the other ZBTB11 targets and brings up a hypothesis that the observed phenotypes in Family A patients might be somehow related to impaired regulation of this gene at some level which remains to be clarified experimentally.

Table 3.

Clinical comparison that highlighted AP4M1, AAAS and SCO1 target ID genes among the possible ZBTB11 targets obtained from ChIP-seq data

Gene names	Syndrome	Cognitive impairment	Microcephaly	Spasticity	Ataxia	Drooling	MRI
ZBTB11	Our patient (Family A), AR	+	+	+	+	+	Cerebellar hypoplasia, mild ventriculomegaly
AP4M1	Spastic paraplegia 50, AR	+	+	+	+	+	Cerebellar atrophy, thin corpus callosum, ventriculomegaly, white matter lesions
AAAS	Achalasia-addisonianism-alacrimia syndrome, AR	+	+	+ (in one patient)	+	–	NA
SCO1	Mitochondrial complex IV deficiency, AR	+	+	+	+	–	Symmetric lesions in the basal ganglia

Prominent reduction in mushroom body size of CkIIα-i1 (ZBTB11 orthologue) knockdown flies

To gain some insights into ZBTB11 function in the brain, we utilized Drosophila melanogaster as an animal model. The role of Drosophila BTB-ZF proteins in regulation of neuronal development is known, and several ZBTB proteins have been identified to be involved in neurological functions, such as fruitless, abrupt, lola and chinmo (10,11). CKII-α subunit interactor-1 (CkIIα-i1) is the known ZBTBT11 ortholog in Drosophila melanogaster which has homology to both human and zebrafish zbtb11, with 38% overall similarity (25% identity) to the human ZBTB11 gene. This gene has higher expression in the larval central nervous system and adult Drosophila brain compared to the other genes with good matches to ZBTB11. Also, it lacks the BTB domain and could be considered as a good candidate to evaluate the possible effects of disruption in the Znf domain of the protein, where the missense variants identified in this study are located. We employed RNAi-mediated knockdown of CkIIα-i1 by applying the UAS-Gal4 system. Relative quantification of CkIIα-i1 (CG6215) knockdown using the pan-neuronal driver elav-Gal4 (Elav^c155-Gal4/+; UAS-RNAi^CkIIα-i1/+) and analysis of RNA isolated from brains of F₁ progeny showed a reduction in CkIIα-i1 expression of ∼70% among three biological replicates (P value: 0.0036; Fig. 4B).

Figure 4.

(A) Confocal image stacks of whole brain mount in Control flies: UAS-RNAi^CkIIα-i1 (first row) and OK107-Gal4/+; UAS-RNAi^CkIIα-i1/+ offsprings (second row) are shown. From left to right: anti-elav (blue) in Column 1, anti-Fasciclin II (green) in Column 2 and the merged image is shown in Column 3. (A’) The zoomed images of MBs are shown (B) Relative quantification of RNAi-mediated CkIIα-i1 (CG6215) knockdown in F1 offspring (Elav^c155-Gal4/+; UAS-RNAi^CkIIα-i1/+). The black bar shows the CkIIα-i1 expression in UAS-RNAi^CkIIα-i1 (Control) and the gray bar shows CkIIα-i1expression in Elav^c155-Gal4/+; UAS-RNAi^CkIIα-i1/+ offsprings (F1). The ∼70% CkIIα-i1 knockdown is obtained from three biological replicates (P-value: 0.0036). (C) The approximate reduction to a third of the corrected total fluorescence intensity (CTCF) of Mushroom Bodies in the F1 generation compared to Control flies (total n = 13 tissue sections) is shown in this chart (obtained from two biological replicates). The black bar shows the UAS-RNAi^CkIIα-i1 (control) mean CTCF; the light gray bar shows the mean CTCF of OK107-Gal4/+; UAS-RNAi^CkIIα-i1/+ offsprings eclosed from crosses of UAS-RNAi^CkIIα-i1 virgin females and OK107-Gal4 single male (P-value: 0.0040); the darker gray bar shows the mean CTCF of OK107-Gal4/+; UAS-RNAi^CkIIα-i1/+ offsprings eclosed from crosses of UAS-RNAi^CkIIα-i1 single male and OK107-Gal4 virgin female (P-value: 0.0004).

RNAi-mediated knockdown in the mushroom body (MB) followed by behavioral assays or studies on brain organization and architecture is one of the common approaches to study the effect of candidate ID genes at whole organism level (4,12). MB is a noticeable neuropil structure in the Drosophila brain and is mainly involved in different forms of learning, short- and long-term memory (12–14). Additionally, the structure and anatomy of the MB is very well defined (15). In this study, the MB morphology of OK107-Gal4/+; UAS-RNAi^CkIIα-i1/+ adult fly brains were examined with α-FasII antibody to visualize α, β and γ neurons of the MB and α-elav to visualize all neuronal nuclei under a confocal microscope. UAS-RNAi^CkIIα-i1 flies were used as controls. A recognizable shrinking of the MB was observed in transgenic flies compared to controls (Fig. 4A and A’), and the approximate reduction to a third of the corrected total fluorescence intensity of the MB staining supports this observation (Fig. 4C). These results were achieved and confirmed through two biological replicates (total n = 13 tissue sections). This observation might be explained by the reduction or degeneration in the number of neurons due to the CkIIα-i1 knockdown in MBs, which may lead to cognitive–behavioral perturbation.

To investigate whether the observed reduction is due to a reduction in the number of neurons, we performed further analysis using anti-dachshund antibody as a marker for Kenyon cell nuclei. Nuclei of Kenyon cells of CKIIa-i1 knockdown lines (n = 4) and controls (n = 5) were quantified using FIJI software, and a significant reduction in the number of neurons has been observed in case of the knockdown (data not shown). While these preliminary observations are exciting, further experiments including the generation of mutant flies that carry the identified mutations and transgenic lines that overexpress the mutant human Zbtb11 constructs are needed to verify these initial data.

Discussion

ZBTB11 gene encodes a classical zinc finger protein with a consensus sequence of Cys2His2 and is predicted to be involved in transcriptional regulation (11). ZBTB proteins are thought to have a critical role in development based on their diverse functions in cell such as transcriptional regulation, chromatin remodeling, cytoskeleton dynamics, protein degradation and several developmental processes such as lymphocyte development, axon guidance and morphogenesis of the gonads. Additionally, these proteins are enriched across vertebrates, suggesting that they might have evolved relatively recently and therefore be involved in key specific features of vertebrate development such as the immune or nervous system (10,11). In human, ID-related diseases are reported for ZBTB16 (Skeletal defects, genital hypoplasia and mental retardation; MIM# 612447), ZBTB18 (Mental retardation, autosomal dominant 22; MIM# 612337), ZBTB20 (Primrose syndrome; MIM#259050), ZBTB24 (immunodeficiency-centromeric instability-facial anomalies syndrome-2; MIM# 614069), ZBTB7A (16), HIC1 (17) and ZBTB40 (18).

ZBTB11 is a little-studied member of the ZBTB protein family, formerly identified as regulator of metallothionein 2A (19) and recently has been reported to be critical in basal neutrophil development and emergency granulopoiesis within PU.1-ZBTB11-TP53 pathway. In a recent study, a Zbtb11 zebrafish mutant named marsanne (mne) with a Cys116Ser change in an N-terminal HHCC domain, located ahead of the BTB and Znf domains, was investigated. Interestingly, this mutant showed central nervous system (CNS) degeneration, impaired craniofacial development and hydrocephalus in addition to a defect in myeloid development. Also, they indicated that Zbtb11 is widely expressed early in development, and later in Zebrafish life its expression is retained only in the nervous system (20). ZBTB11 is predicted to be a Znf transcriptional regulator, and transcription regulation is one of the main cellular processes affected in ID along with neurogenesis, neuronal migration and synaptic plasticity (4).

In this study, we report likely pathogenic variants of the ZBTB11 gene in two consanguineous families and provide evidence to validate this gene as the cause of autosomal recessive ID in these families. The two missense variants (p.H729Y and p.H880Q) change the second conserved histidine residues of the classical Cys2His2 motif: X2-Cys-X2, 4-Cys-X12-His-X3, 4, 5-His. The conserved Cys2His2 residues are responsible for coordination of the zinc ion within the zinc finger. It is known that this coordination is essential for the correct folding of Znf proteins (20), and therefore these missense variants are likely to disrupt the stabilization and folding of the Znf in ZBTB11, leading to altered DNA-binding specificity of the protein.

In light of supporting evidences for considering ZBTB11 as an ID gene, we established an in vitro localization study to evaluate the possible disrupting effects of the two missense variants experimentally. Although in-silico analysis predicts that the missense mutations would affect the Znf structure and thus DNA binding, our results showed a localization abnormality for the H729Y and H880Q mutant protein, with the mutants being excluded from the nucleolus, where the WT recombinant protein was predominantly localized, and supporting the pathogenicity of the identified missense variants. Nucleolar localization of the WT protein was noteworthy itself, as nucleoli are involved in RNA and protein synthesis which is considered as one of the main biological processes required for long-term memory formation (21). In addition, the nucleolus is involved in neuronal growth and injury that works through the “nucleolar stress” signaling pathway and also nucleolar dysfunction followed by rRNA dysregulation may be a mechanism for some neurological disorders including some forms of ataxia (also one of the clinical features observed in some of our patients) (22).

Although ZBTB11 accumulation in the nucleolus might be within a specific phase of the cell cycle, the result of our pathway analysis ascertained from the ENCODE ChIP-Seq data confirmed that the ZBTB11-enriched pathways and biological functions are associated with its localization within the nucleolus (Table 2). This result is in concordance with a study by Liu et al. (23) who showed that six ribosomal protein genes (RPL8, RPL14, RPL41, RPS5, RPS17 and RPS24) were upregulated, whereas ZBTB11 was underexpressed in hepatocellular carcinoma (HCC), marking inverse correlation of ZBTB11 and nucleolus-related genes. In addition, Keightley et al. (20) showed the involvement of Zbtb11 in RNA processing, cell death and survival. Furthermore, they introduced Tp53 as one of the Zbtb11 downstream targets and concluded that the Zbtb11–Tp53-dependent pathway is essential for regulating hemopoietic development and neuronal cell death in response to cell stress such as nucleolar stress, which is consistent with ZBTB11-enriched biological processes in response to stress recognized through ChIP-seq analysis in this study (Table 2).

In vivo functional approaches using animal models can provide strong supportive evidence to validate a candidate gene. Drosophila melanogaster is a good model to study candidate ID genes through different approaches of investigating neuronal morphology or behavioral assays to evaluate learning and memory. The CkIIα-i1 gene similarly encodes an understudied Znf C2H2-type transcription factor and is known as interactor of casein kinase IIα (CkIIα)—a catalytic subunit of Drosophila melanogaster protein kinase CK2 (DmCK2) that is a Ser/Thr protein kinase (24). DmCK2 is known to be involved in neuronal development and synapse maintenance, and the activity of its alpha subunit has been identified in presynaptic motoneurons (25). Likewise, hZBTB11 is identified as a brain-specific protein phosphatase 1(PP1) interacting protein that regulates PP1 targeting and enzymatic specificity via interacting with its catalytic subunit, PPP1CA (26). PP1 is a Ser/Thr protein phosphatase and a well-recognized suppressor of learning and memory. Interestingly, PP1 is present in nucleoli in its enzymatically active form (27), where we localized ZBTB11 protein. The involvement of these two interacting partners (PPP1CA in human and CkIIα in Drosophila) in protein phosphorylation—a key biochemical process underlying synaptic plasticity—is remarkable, since synaptic plasticity is considered to be an essential component of learning and memory (28).

In this study, the observed approximate reduction to a third of the MB size and number of Kenyon cells of OK107-Gal4/+; UAS-RNAi^CkIIα-i1/+ offspring as a possible outcome of neuronal reduction during development or degeneration provides further in vivo evidence for ZBTB11 involvement in our patients’ ID phenotype. Different MB segments composed of specific neurons are shown to be required for different aspects in learning and memory formation, and any alteration in MB structure may lead to disturbance of behavioral functions. It is proposed that the reduction in number of neurons may affect functional neuronal circuits in the brain, which are required for adaptive behavior. There are some gene reports that correlate reduced MB size in mutant flies to defective olfactory learning. One of these Drosophila genes is miniature (mbm) which is interestingly a C2H2 Znf protein localized in the nucleolus and is involved in ribosomal biogenesis. It has been shown that Mbm can be phosphorylated by DmCK2 and is involved in sufficient protein synthesis in neuroblasts controlling proper cell growth. The reduced MB size in Mbm mutants is thus a direct consequence of loss of Mbm function. (14,29–31). Interestingly the zbtb11 orthologe CkIIα-i1 is an interactor of CKII-α (catalytic subunit of DmCK2) and could be associated to the abovementioned network explaining the reduced MB size observed in our transgenic flies. Alternatively, it may act independently via a similar pathway due to the comparable nucleolar functions and localizations identified for hZBTB11 in this study. Therefore, we propose that the reduction in MB size as a result of CkIIα-i1 knockdown may disturb the cognitive and adaptive behavior in our Drosophila model which is comparable to observation of CNS degeneration in the zebrafish mutant.

In conclusion, we report on the first ZBTB11 mutant patients in human (6,7) presenting autosomal recessive ID. Although ZBTB11 has been implicated in myeloid development and malignancy, here we provided some functional evidences of its role in the nervous system. Similarly, there are examples of other ID genes that were well studied for their involvement in hematopoiesis and malignancy along with ID such as BCL11A (C2H2 Znf transcription factor) causing Dias-Logan syndrome (MIM# 617101) or PLZF (Zinc Finger and BTB Domain Containing 16) causing skeletal defects, genital hypoplasia and mental retardation (MIM# 612447). In addition, by investigating the target ZBTB11 genes, having a deeper look at ID-related genes, we have shown that at least 34 ID genes might be considered as possible targets of this transcription factor, highlighting AP4M1 with the most clinical similarity to Family A patients. So, we report the first ZBTB11 mutant patients presenting ID and provide functional evidence linking ZBTB11 dysfunction to this phenotype.

Materials and Methods

SNP genotyping, homozygosity mapping and linkage analysis

Genomic DNA was extracted according to the standard salting out method, and whole-genome SNP genotyping was performed using Genome-Wide Human SNP array 6.0 and Human Mapping 250K NspI array (Affymetrix, Santa Clara, CA, USA) for Families A and B, respectively. For Family A, multipoint linkage analysis with the help of Alohomora software (32) and the Merlin program (33) was performed in order to calculate LOD scores considering autosomal recessive mode of inheritance, complete penetrance and disease allele frequency of 10⁻³. For Family B, HBD regions were determined using Homozygosity mapper (34) and FSuite (35).

Targeted NGS, whole-exome sequencing

Proband’s DNA of Family A (III-2) underwent target enrichment using a custom-made Agilent SureSelect DNA capture array (Agilent Technologies Inc, Santa Clara, CA, USA) targeting HBD regions and was then sequenced on Illumina Genome analyzer II (Illumina, San Diego, CA, USA) with 96.9% coverage at 20× and average depth of 73×. Data analysis was performed by Medical Resequencing Analysis Pipeline (MERAP) (6,36). In addition, to evaluate additional pathogenic variants beyond the homozygous intervals, proband’s DNA of Family A (III-2) was later subjected for WES, underwent enrichment by Agilent SureSelectXT Human All Exon V6 (Agilent Technologies Inc, Santa Clara, CA, USA) and then sequenced on NovaSeq 6000 (Illumina, San Diego, CA, USA) with 94.4% coverage at 20× and average depth of 104×. Data analysis was performed using GATK pipeline. For Family B, WES was performed on one affected member but with analysis of just the single (37Mb) HBD locus identified (7). The potential causal variants were validated by Sanger sequencing, and co-segregation studies were performed for all the available and informative members of each family.

Cloning, mutagenesis and cell culture

The complete human ZBTB11 open reading frame was subcloned from pCS2mt-GFP into expression vector pcDNA^TMDEST53 (N-terminal GFP) using Gateway Technology system (Life Technologies, Carlsbad, CA). Site-directed mutagenesis was accomplished using the Quick-Change Lightning Site-Directed Mutagenesis Kit (Agilent Technologies, Santa Clara, CA), to generate constructs for the two mutations, p.H729Y and p.H880Q. The Human embryonic kidney cells 293T (HEK293T) (ATCC, Manassas, VA) were transiently transfected using PolyFect Transfection Reagent (Qiagen, Hilden, Germany) with WT and mutant GFP-ZBTB11 constructs according to the manufacturer’s protocol. To investigate hZBTB11 subcellular localization, the transfected cells were fixed with 4% paraformaldehyde followed by DAPI (4', 6-diamidino-2-phenylindole) staining to mark nuclear DNA and then visualized by confocal microscopy. Images were acquired using Olympus FluoView FV1200 confocal microscope (z-stacks with a frame size of 1024 × 1024 pixels).

Immunofluorescent staining

IF was accomplished using Monoclonal Anti-B23 (Nucleophosmin/NPM) (#B0556; Sigma-Aldrich, St. Louis, Missouri) with Goat Anti-Mouse IgG Antibody, Cy3 conjugate (#AP124C; EMD Millipore, Billerica, Massachusetts) secondary antibody. Transfected HEK293T cells were fixed with 4% paraformaldehyde, permeabilized with phosphate-buffered saline (PBS), Tween-20, 0.05% and blocked in 10% goat serum for 2 hours. Cells were incubated in primary (1: 500) and secondary antibodies (1: 1000) dilutions and then counterstained with DAPI and visualized using confocal microscopy as described above. In addition, anti-ZBTB11 Antibody (#HPA015328; Atlas Antibodies, Stockholm, Sweden) with Goat Anti-Rabbit IgG H&L (FITC) (#ab6717; Abcam, Cambridge, USA), secondary antibodies were used to assess expression of endogenous ZBTB11 in HEK293 cell lines.

ENCODE ChIP-seq data analysis

We applied GREAT (8) to investigate the enriched pathways, molecular functions and biological processes of available ZBTB11 ChIP-seq data deposited in the ENCODE database (9). In addition, the final ChIP IDR peaks for each of the above cell lines were annotated using the Peak Annotation and Visualization (PAVIS) web-based tool (37), and Galaxy web-based platform was applied for the comparisons (38).

In order to associate ZBTB11 into its possible neuronal targets and also linking this gene to the network of ID-related genes, the final ChIP IDR peaks were compared to the lists of 1095 ID-related genes and 2284 candidate genes introduced through Deciphering Developmental Disorders (DDD) Study (Supplementary Material, Table S3) (39). Then, we confined the analysis to the common target peaks identified in Upstream, 5’-UTR and Exon/CD regions and identified 75 and 46 common target genes among DDD and ID gene lists, respectively (Supplementary Material, Table S3). Additionally, the 46 target ID gene list was refined to a list of 34 ID genes by choosing the target peaks located within a known CpG islands.

Drosophila melanogaster stocks and maintenance

The Drosophila stocks were maintained on standard media under 12h:12h light–dark cycle at 25°C. The following stocks were purchased from Bloomington Stock Center: #60102, UAS-RNAi stock expressing dsRNA for RNAi of CkIIα-i1 (CG6215) under UAS control; #854, mushroom body driver (OK107-Gal4) and #458, pan-neural driver (elav^c155-Gal4).

Evaluation of RNAi-mediated knockdown by real-time quantitative polymerase chain reaction

To evaluate RNAi-mediated CkIIα-i1 knockdown in the F₁ progeny (Elav^c155-Gal4/+; UAS-RNAi^CkIIα-i1/+) eclosed from crosses between UAS-RNAi single males and elav^c155-Gal4 virgin females, total RNA was extracted from ∼50 fly heads using TRIzol-Chloroform protocol followed by complementary DNA synthesis from 1 µg of RNA with QuantiNova™ Reverse Transcription Kit (Qiagen, Hilden, Germany). The SYBR Green-based quantitative polymerase chain reaction (qPCR) was performed using an ABI7500 real-time PCR instrument (Applied Biosystems, Foster City, CA, USA). CkIIα-i1 knockdown level was assessed via the comparative ddCt method analyzed by ABI7500 Software v2.0.6. The expression level was normalized to the housekeeping gene rp49 (CG7939).

Immunohistochemistry and confocal microscopy of Drosophila melanogaster brain

Adult Drosophila melanogaster brains were dissected, fixed in 4% paraformaldehyde, washed in PBS, 0.3% Triton X-100 and then blocked in PAXD (PBS containing 5%BSA, 0.3% TritonX-100 and 0.3% Sodium Deoxycholate) for 2 hours. Brains were incubated in 1D4 anti-Fasciclin II (DSHB, University of Iowa, #AB_528235) (1: 200 dilution) and Rat-Elav-7E8A10 anti-elav (DSHB, University of Iowa, #AB_528218) (1: 20 Dilution) or anti-dachshund (DSHB, University of Iowa, #AB_579733) (1: 100 Dilution) (overnight at 4°C, washed with PBS, 0.3% Triton X-10 and then incubated in goat anti-Mouse IgG (H + L), Alexa Fluor^® 647 conjugate (Invitrogen, #A-21236) (1: 800 dilution) and goat anti-Rat IgG (H + L), Alexa Fluor^® 555 conjugate (Invitrogen, #A-21434) (1: 800 dilution) at room temperature. These brains were mounted using Vectashield Mounting Medium (Vector Laboratories, Burlingame, CA, USA), and images were taken using a Leica TCS SP5 confocal microscope (Leica Microsystems GmbH, Wetzlar, Germany) with the same settings for all samples: z-stack images of 1.5 µm and 1024 × 1024-pixel resolution. The captured images were analyzed by Fiji distribution of Image J software 1.51 (40), and fluorescent signals of MBs were quantified according to the following formula: CTCF (Corrected Total Cell Fluorescence) = Integrated Density − (Area of selected cell × Mean fluorescence of background readings) (41). Also, confocal stacks through Kenyon cell layers of the MB calyx were processed in FIJI 1.52b (40) by correcting for background (rolling disk radius: 8 px), filtering for noise and smoothing before segmentation of stacks by the 3D watershed plugin of the 3D Suite (42) using empirically determined seed and image thresholds and a seed size of 12 px. Label-positive cells were quantified using the 3D Object Counter v2.0 (43), with a minimum size filter of 100.

Acknowledgements

This work could not be accomplished without the participation of our patients and their families. We are grateful for their patience and collaboration during the project. We acknowledge Prof. Graham Lieschke (Lieschke Lab, Australian Regenerative Medicine Institute) for providing the complete human ZBTB11 open reading frame cloned to pCS2mt-GFP vector as a gift. We declare no potential conflict of interest in relation to this work described here.

Conflict of Interest statement. None declared.

Funding

This research project was supported by Iran National Science Foundation (INSF); grant No.920384578, the Canadian Institutes for Health Research MOP-102758, EU-FP7 project GENCODYS (grant no. Health-241995) and Higher Education Commission (HEC) Pakistan Grant (NRPU-1118).