GLI1 inactivation is associated with developmental phenotypes overlapping with Ellis–van Creveld syndrome

Abstract

GLI1, GLI2 and GLI3 form a family of transcription factors which regulate development by mediating the action of Hedgehog (Hh) morphogens. Accordingly, inactivating variants in GLI2 and GLI3 are found in several developmental disorders. In contrast, loss-of-function mutations in GLI1 have remained elusive, maintaining enigmatic the role of this gene in the human embryo. We describe eight patients from three independent families having biallelic truncating variants in GLI1 and developmental defects overlapping with Ellis–van Creveld syndrome (EvC), a disease caused by diminished Hh signaling. Two families had mutations in the last exon of the gene and a third family was identified with an N-terminal stop gain variant predicted to be degraded by the NMD-pathway. Analysis of fibroblasts from one of the patients with homozygous C-terminal truncation of GLI1 demonstrated that the corresponding mutant GLI1 protein is fabricated by patient cells and becomes upregulated in response to Hh signaling. However, the transcriptional activity of the truncated GLI1 factor was found to be severely impaired by cell culture and in vivo assays, indicating that the balance between GLI repressors and activators is altered in affected subjects. Consistent with this, reduced expression of the GLI target PTCH1 was observed in patient fibroblasts after chemical induction of the Hh pathway. We conclude that GLI1 inactivation is associated with a phenotypic spectrum extending from isolated postaxial polydactyly to an EvC-like condition.

Introduction

Hedgehog (Hh) signaling is an evolutionary conserved pathway of intercellular communication which is one of the major regulators of vertebrate development. This pathway is also reactivated in the adult during tumorigenesis and is involved in stem cell maintenance and tissue regeneration (1–3).

Mammals express three tissue-specific Hh genes, Sonic-, Indian- and Desert-Hh, all of which encode signaling proteins that share the same intracellular transduction pathway (4). Hh proteins act in the receiving cells activating a dose-dependent transcriptional response, that results from Hh-induced modifications of the balance between repressor and activator functions of the GLI family of transcription factors: GLI1, GLI2 and GLI3 (5). GLI proteins interact with DNA through five C2H2 zinc (Zn) finger motifs and have a transactivation domain at the C-terminus. GLI3, and to a lesser extent GLI2, contain an additional N-terminal repressor domain, and in the absence of Hh are transformed into transcriptional repressors (GLIR) following proteolytic cleavage of their C-terminal region. In contrast, GLI1 lacks this repression domain and is considered to be a constitutive activator which functions as an amplifier of the pathway, given that its transcription is upregulated by Hh signaling (5–8). Genomic amplification of GLI1 has also been observed in several types of cancers where is thought to increase proliferation and self-renewal of cancer stem cells (5,9,10).

A distinctive characteristic of vertebrate Hh signaling is the requirement of primary cilia and the intraflagellar transport system (IFT) for signal transduction (11). Key components of the pathway including, PTCH1, SMO and the GLI factors move in and out the cilium in response to Hh ligands; and both, GLI3FL proteolytic cleavage and the disassembly of inhibitory SUFU-GLI complexes, which is a critical step for the generation of functionally competent GLI-activators (GLIA), depend on primary cilia (12–16).

In line with the central role of the GLI factors in Hh signaling, heterozygous mutations in GLI3 are accountable for several dominant traits of variable severity including Greig cephalopolysyndactyly syndrome (GCPS; MIM: 175700) (17,18), Pallister-Hall syndrome (MIM: 146510) (19), preaxial polydactyly type IV (MIM: 174700) and postaxial polydactyly types A1 and B (MIM: 174200) (20,21). These conditions mainly involve limb and craniofacial abnormalities and correlate with the location of the mutation (22). Similarly, GLI2 heterozygous mutations are associated with holoprosencephaly (HPE9; MIM: 610829) and Culler-Jones syndrome (MIM: 615849) (23,24), which are featured with hypopituitarism and craniofacial defects often in combination with polydactyly. Furthermore, pathogenic changes in Hh pathway components regulating the proportion between GLIR and GLIA also lead to birth defects that phenotypically overlap with GLI-disorders. The recessive entity Ellis–van Creveld syndrome (EvC; MIM: 225500) belongs to this group of GLI-related conditions. EvC patients have short limbs and ribs, postaxial polydactyly, teeth and nail defects and in 60% of cases abnormal cardiac septation. The large majority of EvC cases are described with two loss-of-function mutations in any of the two subunits of the EvC ciliary complex: EVC and EVC2 (25,26). Biallelic inactivation of either of these genes leads to increased levels of GLI3R after Hh pathway activation due to incomplete dissociation of SUFU-GLI3FL complexes, and reduced expression of Hh-target genes (27–29). A small fraction of EvC patients lack mutations in EVC or EVC2 (30,31) and three of these cases were identified having mutations in the retrograde IFT component WDR35. WDR35-EvC variants impair the ciliary localization of EVC, EVC2 and SMO and also lead to higher levels of GLI3R in the presence of an active pathway (32).

Remarkably, despite being a relevant component of Hh signaling, up to now no GLI1 mutations had been described in any malformation syndrome, thus maintaining enigmatic the phenotypic consequences of GLI1 inactivation in human development. In this article, we describe eight patients from three unrelated families carrying three different homozygous truncating variants in GLI1.

Results

Biallelic inactivation of GLI1 in patients with developmental phenotypes of variable severity

We commenced this study by analyzing a 27-days-old male patient, born to first-cousin Turkish parents, who was the only one affected among eight siblings. Clinically, this child was diagnosed as having features reminiscent of EvC including short stature (third percentile) with shortening of distal lower extremities, bilateral postaxial hexadactyly of hands and feet, mild nail dysplasia and atrial septal defect (ASD) (Fig. 1). To ascertain the genetic etiology underlying the condition of this infant, we sequenced all coding exons and flanking intronic regions of EVC and EVC2 in DNA of the patient but no pathogenic changes were identified. In addition, since cells carrying EVC, EVC2 or WDR35-EvC mutations, all exhibit increased levels of GLI3R following Hh signaling stimulation (27–29,32), we evaluated GLI3 processing in skin-derived fibroblasts from the proband to exclude pathogenic changes such as deep intronic variants. To assess this, we treated patient and control fibroblasts with the SMO agonist SAG or its vehicle (DMSO) and for each cell type we calculated the variation in the GLI3R/GLI3FL ratio (ΔGLI3R/GLI3FL) between the two conditions by subtracting the SAG value from the corresponding figure obtained in DMSO. This experiment revealed similar ΔGLI3R/GLI3FL between patient and normal control cells, thus suggesting no mutations in known EvC genes (Fig. 2a and b and Supplementary Material, Fig. S1a). Furthermore, western blot analysis and immunofluorescence in patient fibroblasts demonstrated normal protein levels and subcellular localization of EVC, and indirectly also of EVC2 given that the two proteins are codependent for localizing to cilia and for maintaining their levels of expression (Fig. 2a and c and Supplementary Material, Fig. S1b) (27). Apart from this, we evaluated the ciliation percentage, cilia length and ciliary translocation of SMO in patient and normal control fibroblasts and observed in all cases comparable numbers (Fig. 2d–g). These features are commonly altered in disorders associated with defects in IFT or transition zone proteins of the cilium, some of which lead to ciliopathy phenotypes partially overlapping with EvC (33). Altogether these results indicated that the causative variant in our proband was not in any of the two subunits of the EvC complex or WDR35, did not affect cilia morphology and should lie downstream of GLI3 processing and consequently downstream of SMO activation.

Figure 1.

 EvC-like phenotypic features of patient 1 with homozygous inactivation of GLI1 (GLI1-W780X). (a) Clinical image of patient 1 demonstrating shortening of lower limbs and postaxial polydactyly. (b) Bilateral postaxial polydactyly of hands, closer view. (c) X-Ray image showing shortening of long bones of the lower limbs, more pronounced in tibia. (d–g) Clinical and X-Ray images showing bilateral postaxial polydactyly of feet. There is no sixth metatarsal (f and g). In addition to these features, patient 1 had hydrocele and ASD. Reddish color of nails in (a) is due to henna application.

Figure 2.

 Evaluation of candidate genes in patient 1. (a) Top panel: GLI3 immunoblot from primary fibroblasts treated with SAG (+) or its vehicle (−). After SAG incubation, extinction of GLI3FL cleavage follows a normal pattern in patient 1 cells (P), whereas GLI3R is increased in SAG-treated EVC^−/− and EVC2^−/− cells (an independent replica of this experiment is in Supplementary Material, Fig. S1a). Normal control fibroblasts are labeled as C and tubulin (TUB) served as loading control. Lower panel: anti-EVC analysis of the same protein extracts demonstrating normal levels of EVC in patient 1. Densitometric quantification of EVC protein levels is in Supplementary Material, Fig. S1b. (b) Analysis of ΔGLI3R/GLI3FL = [(GLI3R/GLI3FL)_DMSO − (GLI3R/GLI3FL)_SAG] calculated by densitometry from the blot on the top. Data are normalized to ΔGLI3R/GLI3FL of normal control cells. (c) Sifted-overlay immunofluorescence images showing EVC at the base of cilia in normal control (C) and patient 1 (P) fibroblasts but not in EVC2^−/− or EVC^−/− cells. Red: EVC; green: acetylated-TUB + ɣ-TUB and blue: nuclei. Scale bars: 10 µm. (d, e) Ciliation percentage (d) and cilia length (e) in patient 1 and normal control cells. (f) Immunofluorescence (shifted-overlay) illustrating ciliary translocation of SMO following SAG incubation in both patient 1 (P) and normal control cells (C). Red: acetylated-TUB; green: SMO and blue: nuclei. Scale bars: 10 µm. (g) Quantification of ciliary presence of SMO showing normal trafficking of this protein in response to SAG in patient 1 fibroblasts. For (d) and (g), a minimum number of 100 cells or 100 cilia were analyzed, respectively; 50 cilia were measured in (e). DMSO cultures from (f) were used in (d, e) analyses. (h) Schematic representation of GLI1FL (full-length) and C-terminal variants identified in patients of this report. Number of residues, DNA binding motifs (Zn) and the critical transactivation domain (TAD) are indicated. (i) Representative GLI1 immunoblot of patient 1 (P) and normal control fibroblasts (C) treated with SAG (+) or DMSO (−) proving upregulation of GLI1-W780X in response to SAG. Densitometric quantification of GLI1 protein levels is in Supplementary Material, Fig. S1c. (j) qRT-PCR quantification of the increase in PTCH1 and GLI1 transcript levels induced by SAG in patient 1 and normal control cells using GUSB as the reference gene. Y axis represents 2^−ΔΔCt normalized to the control cells value. For each sample, the ΔCt of the corresponding DMSO condition was used as calibrator. Quantification of transcripts with respect to ACTB, GAPDH and HPRT1 yielded similar results (Supplementary Material, Fig. S1d).

As exclusion of known EvC genes and parental consanguinity were suggestive of a novel autosomal recessive gene, we performed SNP arrays in the proband and four of his normal siblings to determine the pool of patient-specific homozygous regions, and concurrently conducted whole-exome sequencing (WES) in the proband. As expected, the affected child was not homozygous for EVC-EVC2 or WDR35 in the arrays. Filtering of WES data was carried out assuming a recessive model of inheritance and that the mutation was homozygous by descent. This analysis identified 38 homozygous variants with MAF ≤ 0.01 within the entire collection of homozygosity blocks corresponding to the patient (Supplementary Material, Table S1a and b). Four variants were located inside patient-specific homozygous regions and among them, we selected a GLI1 nonsense change (NM_005269.2: c.2340G > A; p.Trp780*) predicting truncation of the protein at codon 780 (GLI1-W780X) as the best candidate due to the association of EvC with diminished Hh signaling (Fig. 2h). In addition, the GLI1-W780X variant had the highest predicted pathogenic impact, was absent from any of the public genetic variant databases including dbSNP, ExAC, EVS, gnomAD, 1000 Genomes and Kaviar (Supplementary Material) and was embedded into a 13 Mb stretch of homozygosity that was present exclusively in the proband. Using Sanger sequencing the GLI1-W780X variant was confirmed to be in the homozygous state in the patient, while analysis of four available unaffected siblings showed three brothers with wild-type sequence and a sister with the mutation in the heterozygous state (Fig. 3). Subsequently, we became aware that the proband (from now patient 1) belonged to a large kindred containing numerous individuals, many born to consanguineous couples, with upper limb postaxial polydactyly and hypothesized that they might carry the same homozygous mutation (Fig. 3). We had access to DNA from three additional patients (patients 2–4) and some of their relatives and we checked them for the presence of the GLI1-W780X variant. Patients 2 and 3 (8 years and 34 months of age, respectively) had bilateral postaxial polydactyly of hands and height values in the 10–25th percentile. Patient 4 was reported to have had either a unilateral or bilateral postaxial floating sixth finger since there was no information recorded after correction of polydactyly during infancy. Each of these patients was born to a different consanguineous couple. Following DNA sequencing, patient 2 was discovered to be homozygous for the GLI1-W780X variant and segregation analysis in his two unaffected siblings was consistent with GLI1-W780X being the causal change as in patient 1. However, patient 3 had the GLI1-W780X variant in the heterozygous state and patient 4 showed wild-type sequence (Fig. 3). Sequencing of all coding exons of GLI1 and adjacent intronic sequences in patients 3 and 4 revealed no additional significant changes. Hence, to exclude the possibility that a variant different to GLI1-W780X could be accountable for the phenotype of all four patients, we performed homozygosity mapping and WES in patients 2 and 3 and compared the results with those obtained in patient 1. Array hybridization showed almost no homozygous regions shared by the three patients except for a short fragment of 7.5 Kb on chromosome 8 (Supplementary Material, Table S2). In contrast, GLI1 was detected within a large 63.6 Mb region of homozygosity in patient 2, as expected for a disease variant inherited from a recent common ancestor to both parents (Fig. 3 and Supplementary Material, Fig. S2) (34). On the other hand, analysis of total WES data from patients 1–3 to reveal homozygous variants common to the three patients identified 11 homozygous changes with MAF ≤ 0.01, none of which was predicted as pathogenic (Supplementary Material, Table S3a and b). Further to this, we conducted SNP arrays in patient 4 and evaluated WES variants from patient 3 comprised within homozygous regions shared between patients 3 and 4. In the same manner, variants within homozygous blocks shared by patients 3 and 2 were also examined. In both cases, no damaging candidate changes were detected (Supplementary Material, Table S4a–d). Therefore, these results were consistent with GLI1-W780X being responsible for the phenotype of patients 1–2 and a different genetic variation or disease mechanism for the phenotypes of patients 3 and 4. These patients will be discussed later in the context of the remaining data. Finally, given that functional redundancy and compensation between GLI1-3 has been reported in mice (35), we investigated the presence of rare coding variants in GLI2 and GLI3 in patients 1–3 using WES data and Sanger sequencing. Of the three patients, only patient 1 was identified with two heterozygous missense changes in GLI2: M1352V and D1520N (Supplementary Material, Table S5). These two variants have relatively high ExAC allele frequencies ∼0.01 and are currently considered as variants with uncertain significance (36). Three unaffected siblings of patient 1 were also heterozygous for both M1352V and D1520N (Fig. 3).

Figure 3.

 Pedigree structure of the kindred corresponding to patient 1. The proband (patient 1, IV-8) who is indicated with an arrow, patient 2 (V-7), patient 3 (V-13) and patient 4 (II-10) are marked with asterisks. For II-10, it is uncertain whether she had unilateral or bilateral postaxial polydactyly due to lack of medical records. Chromatograms demonstrate the genotypes of the four patients with respect to the G > A transition (underlined nucleotides) causing GLI1-W780X. Red and blue + and − symbols underneath individuals indicate their corresponding allelic composition regarding GLI1-W780X (red) and GLI2: M1352V + D1520N (blue) variants. Wild-type alleles are designated with + whereas the variant alleles are indicated by −. Due to limited information on the age of each individual, the order of siblings in every generation is arbitrary.

We next sequenced GLI1 in 27 independent EvC or EvC-like cases previously excluded for EVC and EVC2 mutations and identified a second proband (patient 5) born to second-cousin parents with a different homozygous truncation variant in GLI1 (NM_005269.2: c.1930C > T; p.Gln644*; GLI1-Q644X) (Fig. 2h). The new change, like GLI1-W780X, was located in the C-terminal exon of the gene and was found to be absent from any of the public genomic databases (dbSNP, ExAC, EVS, gnomAD, 1000 Genomes and Kaviar). Patient 5 (II-5) was an 11-year-old Turkish boy featured with short stature (3rd–10th percentile), bilateral hexadactyly of hands and feet, genu valgum and high palate (Fig. 4a and b). He had a similarly affected older brother (patient 6; II-4) with short stature (third percentile) and bilateral postaxial polydactyly of hands and feet, and three more siblings of whom one was normal (II-1), and two (II-2, II-3) were characterized with high palate and intellectual disability. Sibling II-2 also had postminimi polydactyly of the left hand (Fig. 4b–e). SNP arrays had previously identified GLI1 inside of homozygosity blocks of 13.5 Mb and 11.3 Mb in patients 5 and 6, respectively, but not in the rest of brothers. Consistently, segregation analysis of GLI1-Q644X among the five siblings showed the two brothers with short stature and polydactyly (patients 5 and 6) being homozygous for the GLI1-Q644X mutation, whereas II-1 and II-2 had the mutation in the heterozygous state and II-3 presented wild-type sequence (Fig. 4b and e). To further confirm the causality of the GLI1 mutation in this family, we conducted WES in patient 5 and two of his siblings not affected with polydactyly (II-1 and II-3). During the filtering process, variants from patient 5 contained within chromosomal homozygous regions common to patients 5 and 6 were selected, and subsequently those present in the homozygous state in any of the non-polydactyly brothers were filtered out. As a result, we obtained a short list of variants in which GLI1-Q644X was identified as the most convincing pathogenic variant. Importantly, GLI1 was the only mutated gene in common with the filtered exome of the proband of family 1 (Supplementary Material, Tables S1b and S6a–c), thus providing additional support for the link between GLI1 and the polydactyly phenotypes of both families. Analysis of GLI2 by Sanger sequencing in patients 5 and 6 revealed three frequent polymorphisms only in patient 6 (rs3738880, rs12711538 and rs10167980), all detected in the homozygous state.

Figure 4.

 Phenotypic features of patients 5 and 6 with GLI1-Q644X homozygous truncation and their siblings. (a) Clinical image of patient 5 (II-5) at 11 years of age demonstrating short stature, genu valgum deformity (horizontal arrow) and scar from removal of the sixth toe (arrow). (b) Pedigree structure of the family of patients 5–6 (asterisks) showing inheritance of two different traits among siblings: (1) bilateral postaxial polydactyly of hands and feet + short stature and (2) mental retardation and high palate. Red + and − symbols refer to wild-type and Q644X alleles, respectively, and sequencing chromatograms illustrate wild-type sequence and the homozygous mutant C > T transition (underlined nucleotides) leading to GLI1-Q644X in patients 5–6. The order of siblings is arbitrary and the proband (patient 5) is designated with an arrow. (c) Hands of patient 6 (II-4) showing bilateral scars (arrows) from removal of the sixth finger in both hands. (d) Hands of individual II-2 showing rudimentary postminimi polydactyly of the left hand (arrow). (e) Summary of clinical findings and GLI1-Q644X genotypes of patients 5–6 and their siblings. In the GLI1 genotype column, wild-type and Q644X alleles are designated as + and −, respectively. P, percentile. Polydactyly of feet in patients 5–6 was reported to include an additional metatarsal but there was no extra metacarpal in hands.

Lastly, a third stop gain variant in GLI1 (NM_005269.2: c.337C > T; p.Arg113*; GLI1-R113X), truncating the protein before the Zn finger domains, was identified in a Pakistani family with four individuals (patients 7–10) from two different family branches affected with postaxial polydactyly type A. The patients of this family all had bilateral postaxial polydactyly of feet, except individual VI-2 who additionally presented postaxial polydactyly of the left hand. Stature was normal respecting the ethnically matched control population and no duplication of metacarpal/metatarsal was detected (Fig. 5a–d). Variant identification in this family was achieved by whole genome sequencing (WGS) in two affected subjects, VI-1 and VI-7, each from a different branch of the family (Fig. 5a). Based on the pedigree information and being the patients offspring of consanguineous parents the inheritance pattern of disease was predicted to be autosomal recessive. Consequently, annotated variants were filtered for homozygosity and a minor allele frequency < 0.01 and had to be shared by both affected individuals. As a result, the GLI1-R113X change was found in homozygosis in both patients, VI-1 and VI-7, and was considered the strongest candidate for the disease phenotype because of the implication of GLI1 in Hh signaling. Segregation of the c.337C > T variant within the pedigree was further confirmed by Sanger sequencing. All four affected subjects were found to be homozygous for the c.337C > T variant, while their parents and available unaffected siblings resulted to be either heterozygotes or had wild-type sequence (Fig. 5a). On database examination the GLI1-R113X variant was found to be present in a single carrier in gnomAD populations (allele count: 1; total allele number: 246 242), and to have a general allele frequency of 0.000004061; data that are comparable to many other recessive alleles associated with rare conditions.

Figure 5.

 Pedigree and phenotypic features of patients homozygous for the GLI1-R113X variant. (a) Pedigree structure corresponding to the family of patients 7–10 (VI-1, VI-2, VI-7 and VI-8). GLI1-R113X genotypes are indicated underneath the tested individuals. Normal and mutant alleles are referred with + and − red symbols, respectively. Chromatograms show normal, heterozygous and homozygous sequence for the c.337C > T variant (underlined nucleotides). (b and c) Clinical images of feet and left hand of patient VI-2 showing postaxial polydactyly type A. (d). X-Ray image of left hand of VI-2 proving duplication of phalanges without extra metacarpal.

GLI1-W780X transcriptional activity is severely impaired

As GLI1-W780X and GLI1-Q644X are both located in the last exon of the gene, the corresponding transcripts, despite containing a premature termination codon, are predicted to avoid non-sense mediated mRNA decay degradation (NMD) and be translated into mutant proteins (37). Indeed, similar to control cells, western blot analysis in fibroblasts from patient 1 demonstrated that the GLI1-W780X truncated protein is not only synthesized by patient cells but also upregulated in response to SAG (Fig. 2i and Supplementary Material, Fig. S1c). Testing of GLI1 in the rest of patients was not possible due to unavailability of biological material for cell culture analysis. Nevertheless, the GLI1-R113X variant is placed in the third coding exon of the gene and thus predicted to be subjected to NMD decay (37). Structurally, GLI1-Q644X and GLI1-W780X preserve the five C2H2 DNA binding motifs localized between residues 237 and 387 but lack the last 463 (GLI1-Q644X) or 327 (GLI1-W780X) C-terminal amino acids which contain the transactivation domain including the previously reported critical 1020–1091 residues (Fig. 2h) (38). To assess the pathogenic effect of these mutations, we treated patient 1 and control cells with SAG to activate the Hh pathway, and compared the extent of the increase in mRNA levels of the downstream targets PTCH1 and GLI1 in each culture by qRT-PCR using GUSB as reference gene. PTCH1 contains a well-characterized GLI-binding site in its promoter previously shown to be highly responsive to GLI1 and is considered to be a direct transcriptional target of GLI proteins (39). The result of this study showed the increase of PTCH1 expression in patient cells to be 44% lower than normal cells, while the upregulation of the GLI1 transcript resulted in non-significant differences between both cultures (Fig. 2j). This data indicated that Hh signaling was altered in the patient at a step downstream of the induction of GLI1, in keeping with GLI1-W780X pathogenicity. Validation of this expression assay was obtained using different reference genes (Supplementary Material, Fig. S1d). We then undertook a series of experiments to functionally characterize the more distal truncating variant GLI1-W780X. First, we tested the transcriptional activity of the mutant protein using luciferase reporter assays in HEK293T. Cotransfection of these cells with a DNA mixture comprising constructs expressing the entire coding region of GLI1 (GLI1FL) or the patient variant (GLI1-W780X), a GLI-responsive luciferase reporter harboring multiple GLI-binding sites and a third plasmid leading to constitutive expression of Renilla luciferase, showed that the transcriptional activity of GLI1-W780X was severely reduced (Fig. 6a and Supplementary Material, Fig. S3a). In addition, similar GLI1 luciferase reporter assays in which GLI1FL (human or mouse) or a murine activator form of GLI2 (GLI2ΔN) (40) were expressed alone or in combination with GLI1-W780X or its murine equivalent GLI1-W781X also proved reduced expression of the reporter in the presence of GLI1-W780X (or GLI1-W781X), indicating that the mutant variant could compete with other GLI factors for binding to promoters containing GLI-binding sites (Fig. 6a and Supplementary Material, Fig. S3b–d). Next, we analyzed the function of GLI1-W780X in the context of a Hh-driven physiological process in cultured cells. To do this, we retrovirally expressed human GLI1FL or GLI1-W780X in C3H10T1/2 mesenchymal cells, which undergo osteogenic differentiation in response to Hh signaling (41), and assessed expression of the osteoblast marker alkaline phosphatase (ALP) in the presence of SAG or its vehicle. Analysis of ALP enzymatic activity both in protein extracts and by cytochemical staining revealed clear differences between cells expressing GLI1-W780X or GLI1FL, with the mutant protein being practically unable to induce osteoblast differentiation in the absence of SAG (Fig. 6b and c and Supplementary Material, Fig. S3e). Likewise, retrotransduction of NIH3T3 cells with the same constructs demonstrated negligible induction of endogenous Ptch1 both at the transcript and protein levels in the GLI1-W780X cultures without SAG with respect to GLI1FL cells (Fig. 6d and Supplementary Material, Fig. S3f and g). Similar result was obtained for Gli1 following quantification of the endogenous Gli1 transcript by qRT-PCR in the NIH3T3 retrotransduced cells (Supplementary Material, Fig. S3h). Unlike in the luciferase assays, competition between GLI1-W780X and the endogenous GLI proteins was not observed in these experiments, since SAG-GLI1-W780X cells and SAG-empty vector control cultures showed analogous ALP activity and Ptch1 and Gli1 expression. This is very likely to be due to the GLI1-W780X competition effect being masked by the potent activation of the pathway at the level of SMO mediated by SAG which affects to all repressor and activator forms of GLIs.

Figure 6.

 Functional characterization of the GLI1-W780X variant. (a) Luciferase reporter assays in HEK293T transfected with HA-tagged GLI1 human constructs showing minimal reporter expression (12xGliBS-Luc) in cells transfected with GLI1-W780X. Note that cotransfection of GLI1FL and GLI1-W780X (1:1) significantly inhibits GLI1FL activity on reporter expression (n = 3). Values are normalized against the empty vector control (pcDNA3.1). RLU, relative luciferase units. Anti-HA immunoblot showing expression levels of the transfected GLI1 proteins corresponding to this experiment is underneath. TUB, tubulin. (b, c) C3H10T1/2 treated with retroviruses containing human GLI1 variants demonstrating that only GLI1FL is capable to induce osteogenic differentiation in the absence of SAG. ALP cytochemical staining of these cultures is in (b) and ALP catalytic activity assessed in protein extracts with CDP-star and normalized against protein concentration is in (c) (n = 3). AU, arbitrary units. A representative immunoblot demonstrating expression levels of the retrotransduced human GLI1 variants is in Supplementary Material, Fig. S3e. (d) Representative immunoblot showing expression of endogenous PTCH1 in NIH3T3 retrotransduced with human GLI1 variants using specific antibodies against the mouse protein. Densitometric quantification of endogenous PTCH1 protein levels with respect to TUB is in Supplementary Material, Fig. S3f. Expression of the retrotransduced human GLI1 variants monitored with a human-specific antibody is shown in the underneath panel. Quantification of Ptch1 and Gli1 endogenous transcripts in the same cultures by qRT-RTPCR is in Supplementary Material, Fig. S3g and h.

GLI1-W780X acts in vivo as a non-functional protein

Lastly, we checked the activity of the GLI1-W780X variant in vivo by analyzing the impact of GLI1-W780X expression during neural tube progenitor specification in chicken embryos. This process is essential for the correct differentiation of ventral interneurons and motoneurons in the spinal cord. To do this, the neural tube of chicken embryos of 42 h incubation (embryonic stage HH12–13) were electroporated with human GLI1FL or GLI1-W780X pCIG constructs, or with the empty pCIG vector as control (Fig. 7). Expression of the IRES-driven GFP cassette included in pCIG was used to monitor the electroporated cells. As expected, GLI1FL induced ventralization of the dorsal neural tube, promoting expression of direct targets such as NKX2.2 (Fig. 7b’) and OLIG2 (Fig. 7b’’), which depend on high and moderate levels of GLIA activity, respectively (42–44). Both, NKX2.2 and OLIG2 are responsible for the specification of progenitors for ventral interneurons subtype V3 and for motoneurons, respectively (45). Conversely, GLI1FL caused the loss of dorsal identity specifiers such as PAX6 (Fig. 7e). These results are similar to those previously reported upon SHH or GLIA enhancement (44,46). In contrast, the GLI1-W780X variant was not able to induce dorsal expansion of NKX2.2 expression (Fig. 7c’) and could only activate the expression of OLIG2, which does not require high levels of GLIA activity (Fig. 7c’’) (44). Additionally, GLI1-W780X showed no effect on PAX6 expression (Fig. 7f). Moreover, at more rostral levels of the spinal cord, when the ventral progenitor specification is more established and cells are less competent to change cell fate, GLI1-W780X was not able to induce ectopic OLIG2 expression and there was reduced number of NKX2.2 positive progenitor cells. This reduction in NKX2.2 expression might be due to GLI1-W780X interfering with the activity of endogenous GLI1 and/or GLI2 proteins (Supplementary Material, Fig. S4). These results demonstrated that the human GLI1-W780X variant is not capable of executing the complete differentiation program triggered by GLI1FL of the specification of ventral neural progenitors in the developing neural tube.

Figure 7.

 Immunofluorescence analysis showing that GLI1-W780X is not able to produce the ventralizing effects of GLI1FL in the embryonic central nervous system in vivo. Transversal sections of neural tubes from chicken embryos electroporated only on the right side with constructs expressing EGFP alone as control (pCIG, panels a, a’, a’’ and d) or in combination with human GLI1FL (panels b, b’, b’’ and e) or human GLI1-W780X (panels c, c’, c’’ and f). EGFP (green), NKX2.2 (red), OLIG2 (cyan) and PAX6 (red in d, e and f). Nuclei stained with bisbenzimide are in blue.

Discussion

The critical role of GLI(2–3) in development has been well-stablished and is manifested by the clinical features of patients with mutations in these genes. However, due to the absence of a condition linked to loss-of-function mutations in GLI1, the contribution of this gene to human embryonic development was unknown. Herein, we report eight patients with developmental defects having three different homozygous nonsense mutations in GLI1, hence indicating that all members of the GLI triad are required for human development. This is the first time that patients with biallelic germline inactivation of GLI1 are reported.

We identified GLI1-W780X as the causative variant in patients 1 and 2 using a combination of homozygosity mapping and WES. In the second family, homozygosity mapping followed by candidate gene screening led to the identification of GLI1-Q644X in the homozygous state exclusively in the two brothers featured with bilateral polydactyly of hands and feet and short stature. The causality of this mutation was further supported by WES. Taking together the offspring of the three couples from the two first families in which both, father and mother, were heterozygous for any of the GLI1 nonsense changes identified (that is the parents of patients 1, 2, 5 and 6), the corresponding mutations were demonstrated to segregate with the disease in a total number of 13 children, with the four homozygotes for any of the two GLI1 truncating variants presenting bilateral polydactyly and short stature. The third variant, GLI1-R113X was isolated by comparing WGS variants in two affected members of the same family. GLI1-R113X was proved to be in homozygosis in the four individuals of this family affected with polydactyly but not in the two available normal siblings tested. With regards to allele frequency, only GLI1-R113X was present in gnomAD populations but with an allele frequency (one single carrier) compatible with pathogenicity. Neither GLI1-W780X nor GLI1-Q644X were detected after interrogation of the current databases containing genomic data.

The pathogenic character of GLI1 mutations is further supported by functional studies. We provide cell culture and in vivo evidences demonstrating that GLI1-W780X acts as a loss-of function mutation having only very low residual activity. Consistent with this, patient 1 fibroblasts were characterized with reduced PTCH1 expression after SAG-activation regardless normal GLI3R: GLI3FL levels and normal cilia. Similar outcome is predicted for the GLI1-Q644X variant because it is also located in the final exon of GLI1 and leads to a shorter and structurally similar truncated protein. No protein is expected to be fabricated in patients homozygous for the GLI-R113X change (patients 7–10) due to NMD. Cotransfection assays and in vivo analysis in chicken embryos additionally revealed that GLI1-W780X can interfere with the function of other GLI factors. Thus, this aspect could also be contributing to the phenotype as the GLI1-W780X protein variant was found to be synthesized and upregulated by the Hh pathway in patient cells. These results are in line with previous studies by Ruiz I Altaba using theoretical GLI1 C-terminal truncations engineered to characterize the functional domains of the GLI1 protein (7). A schematic overview indicating the role of GLI proteins in Hh signaling and the pathogenic mechanism underlying GLI1 mutations is shown in Fig. 8.

Figure 8.

 Schematic diagram of vertebrate Hh pathway and proposed pathogenic mechanism associated with GLI1 mutations. (a) In the absence of Hh, the cytoplasmic protein SUFU binds to the full-length GLI(1–3)FL transcription factors stimulating the proteolytic processing of GLI3FL, and to a lesser degree of GLI2FL, into their repressor forms (GLIR). Subsequently, GLIR proteins act in the nucleus inhibiting the transcription of Hh-target genes including PTCH1 and GLI1. (b) Upon arrival of Hh ligands, the 7-pass transmembrane protein SMO translocates into the ciliary compartment where it interacts with the EVC-EVC2 complex at the base of cilia (27–29). SMO in combination with EVC-EVC2 transmit the signal to the SUFU-GLI(1–3)FL inhibitory complexes promoting their dissociation. As a result GLIR production is discontinued, and simultaneously, the GLI(1–3)FL factors are converted into functionally competent activators able to reach the nucleus and activate the transcription of target genes. Consequently, GLI1FL protein levels are upregulated and the pathway is amplified. (c) Similar to (b) but in a mutant GLI1 environment. In this scenario, the expression of Hh-target genes in response to Hh is reduced due to deficient GLI1FL activator function. In addition to this, subjects carrying GLI1 C-terminal truncation mutations synthetize GLI1-truncated proteins (GLI1-X) which can interfere with the action of GLI(2–3)FL by blocking shared promoter sites. GLI1-X protein levels are also increased following Hh stimulation. The repressor, Zn finger DNA binding and activator domains of GLI proteins are symbolized by R, Zn and A, respectively, and the strength of transcription activation in (b) and (c) is represented by a different number of vertical arrows in the nucleus.

Patient 3 deserves comment. This patient was heterozygous for the GLI1-W780X change which was inherited from his mother, while his father was confirmed to have normal sequence. As the family declined to donate tissue, GLI1 expression in this patient was not possible to be tested. Several biological explanations could account for the phenotype of patient 3. For example, there is the possibility that heterozygous mutations in GLI1 could act as a predisposition factor to postaxial polydactyly. Agreeing with this hypothesis, another heterozygote from family 2 (II-2) has postminimi polydactyly of the left hand. Hypomorphic variants in modifier genes, including exonic, intronic or promoter variants, that would normally cause no phenotype, or even environmental factors, could cooperate with the heterozygous GLI1 C-terminal mutations to result in polydactyly in some individuals. In this respect, since GLI1-W780X was shown to be expressed in patient cells and able to interfere with the function of the wild-type GLI1 factor, it might exert certain dominant negative effect in the heterozygotes, thus making them more susceptible to the effect of genetic modifiers or environmental factors. Alternatively, patient 3 could be explained by a second hit in a non-coding intronic or regulatory region of GLI1, by digenic inheritance or by a complete different homozygous genetic alteration not detectable by exome sequencing that could be shared by patient 4. Postaxial polydactyly has high prevalence: 1/3300 to 1/630 in Caucasians and from 1/300 to 1/100 in African-Americans (47), and there are many genes reported causing polydactyly (48). Hence, the possibility of another genetic component playing a role in the polydactyly of patients 3 and 4, who are part of an extraordinarily large and highly consanguineous kindred, is not surprising. It is well-known that the Hh pathway includes many different molecular actors and is influenced by modifier genes and external environmental factors, and that GLI1 is at the very end of this signaling cascade (49,50). Regarding family 2, two brothers had intellectual disability and high palate. However, since one of them was heterozygous (II-2) for the GLI1-Q644X mutation and the other had normal sequence (II-3), these features are not associated with GLI1 inactivation and may correlate with another distinct genetic etiology or result from alternative causes such as environmental or infectious factors. Consanguineous families are recognized to have increased probability of cooccurrence of two rare recessive disorders resulting from homozygous pathogenic variants in different genes (51).

All patients with homozygous inactivation of GLI1 had in common postaxial polydactyly but there is phenotypic variability with respect to the number of limbs affected and the presence of additional EvC-like features including, short stature, ASD, mild nail dysplasia or genu valgum. It is pretty likely that the position of the mutation may exert a role on the variability of the phenotype of GLI1 patients; however, more patients and mutations are required to be able to delineate a definitive theory. Variable expressivity, including inter and intrafamilial phenotypic variability, and even incomplete penetrance are frequently reported among GLI2- or SHH-disorders like holoprosencephaly (23,24,52–54), and phenotypic variability is also a common finding in GCPS (55) or GLI3-polydactyly (56). In the same manner, patients with EvC may not have the full spectrum of symptoms and can show differences for the presence of a heart defect or in how many extremities polydactyly is present even among patients with the same mutation (57,58). Phenotypical differences between subjects carrying identical mutation are normally attributed to different repertory of variants in modifier genes and/or environmental factors (59). In this context, the two GLI2 heterozygous variants (M1352V; D1520N) identified in patient 1, but not in patient 2, could account for the more severe phenotype of patient 1 since previous functional analysis indicated that the double M1352V; D1520N allele may act as a mild hypomorphic form of GLI2 (60). The clinical overlap between GLI1 patients and EvC is reasonable, because both conditions result from increased levels of GLIR versus GLIA, albeit due to a different mechanism. The levels of GLI3R are increased in patients with EvC, while the amount of GLIA is reduced in GLI1 patients.

Considering knockout mice, homozygous mutants for a Gli1 allele lacking the DNA binding domain were viable and phenotypically normal at birth (Gli1^zfd/zfd) (35). However, postnatally, Gli1^−/− mice were reported having lower survival rate during the first 10 days of life and a significantly lower growth rate and weight at 8 weeks of age, indicating that inactivation of Gli1 in mice is not innocuous (61). In keeping with this, the four GLI1 patients with homozygous C-terminal truncations of GLI1 have short stature or are in the lower part of the normal range (patient 1: 3rd percentile, patient 2: 10–25th percentile, patient 5: 3rd-10th percentile, patient 6: 3rd percentile). Heterozygous Gli1^+/− adult mice were also described with decreased osteoblast differentiation and decreased bone mass suggesting that GLI1 is especially required in the skeleton (61). Redundancy and functional compensation were proposed to explain the lack of developmental defects of Gli1^zfd/zfd newborn mice. This was based on the fact that Gli1^zfd/zfd; Gli2^+/zfd double mutants show a severe phenotype including a variable loss of ventral spinal cord cells but the double heterozygotes Gli1^+/zfd; Gli2^+/zfd were normal (35). In general, when human GLI-associated disorders are considered as a whole, the development of the human embryo seems more vulnerable to Hh signaling disturbances than embryonic mouse development. For example, GLI2 heterozygous mutations cause different birth defects but Gli2^zfd/+ mice are detected with no abnormalities (62), and similarly, although polydactyly is a constant finding in EvC, it is not present in Evc or Evc2 knockout mice, which also have no detectable heart defects (27,63,64).

In summary, we disclose the effect of GLI1 inactivation on human development and describe a phenotypic spectrum associated with a deficit of GLI1 activity extending from EvC-like features to isolated postaxial polydactyly, thus placing GLI1 among the genes responsible for these phenotypes.

Materials and Methods

SNP arrays

Analysis of regions of homozygosity was conducted by hybridizing 200 ng of genomic DNA to Illumina CytoSNP-850k BeadChip SNP-based arrays according to the recommendations of the manufacturer (Illumina). GenCall scores < 0.15 at any locus were considered ‘no calls’. GenomeStudio software (Illumina) was used to identify homozygous regions by evaluating B allele frequencies for all SNPs. Given the high degree of consanguinity between the parents of each patient, regions of homozygosity were defined using a minimum-length cut-off of 1.5 Mb of consecutive homozygous SNPs.

WES and Sanger sequencing

WES was provided by Sistemas Genomicos S.L. Targeted regions were captured using SureSelect Human All Exon Target Enrichment kit for 51 Mb (Agilent technologies) and sequenced on a HiSeq2000 sequencing platform (Illumina). Reads were aligned against the human genome assembly GRCh38/hg38. Initial filtering was conducted with Picard-tools (http://broadinstitute.github.io/picard/) and SAM-tools (65) and a combination of two different algorithms VarScan (66) and GATK (67) was used for variant calling. Identified variants were annotated using Ensembl database (www.ensembl.org). Total number of reads and quality metrics of all WES analysis which were performed either at 50× or 100× depth are in Supplementary Material, Table S7. Analysis of mutations by Sanger sequencing was carried out by direct sequencing of PCR products amplified from peripheral blood genomic DNA. PCR fragments were enzymatically treated with exonuclease I/shrimp ALP (GE Healthcare) and sequenced using a dye terminator cycle sequencing kit (Applied Biosystems).

Whole genome sequencing

DNAs of two affected individuals from separate branches of the pedigree (VI-1 and VI-7; Fig. 5a) were sent to Novogene (Sacramento, CA). The sequencing library was prepared using the Truseq Nano DNA HT Sample Preparation Kit (Illumina USA) following manufacturer’s instructions. Sequencing was performed on an Illumina HiSeq X system (HiSeq PE150, Q30 ≥ 80%) to an average depth of 9×. This depth was considered sufficient given the strong expectation of homozygosity. Genome sequencing data were aligned, annotated and filtered using the Genoox platform (https://genoox.com/).

Plasmids and constructs

The pLUT7 hemagglutinin (HA)-GLI1 plasmid (68) containing the full-length human GLI1 cDNA (GLI1FL) tagged at the N-terminal with the HA epitope was obtained from Addgene and used as template to generate a cDNA variant equal to the patient mutation (GLI1-W780X) by PCR. Amplification was performed with Pfx50 DNA polymerase (Invitrogene) using a forward primer placed upstream of the HA-tag and the reverse primer 5′-GATATCGGAACTCACCTCATGTTTC-3′. For HEK293T cotransfections, DNA fragments comprising HA-tagged GLI1FL (HA-GLI1FL) or the HA-tagged truncated variant (HA-GLI1-W780X) derived from pLUT7-HA-GLI1 were cloned in pcDNA3.1 (Invitrogen). For retroviral transduction and electroporation of chicken embryos, the entire coding region of GLI1FL and GLI1-W780X without HA tagging were amplified from pLUT7-HA-GLI1 and ligated to pBABE and pCIG, respectively. The murine Gli1 variant mimicking the human mutation (Gli1-W781X) was similarly generated by PCR using plasmid pcDNA3.1-HisB-Gli1 (Gli1FL) (8) as template, a forward primer corresponding to plasmid sequence and the reverse primer 5′-GCGGCCGCGGGAACTCACCTCAGTTTTC-3′. The resulting PCR product was cloned back in pcDNA3.1-HisB. Mouse Gli2ΔN was engineered by PCR amplification using pcDNA3.1-HisB-Gli2 as template (8) and it was designed so that Met 325 of Gli2 was the first ATG. The amplified DNA fragment was cloned in pcDNA3.1. All constructs were verified by DNA sequencing before use.

Cell culture and retroviral transduction

NIH3T3, C3H10T1/2, HEK293T (all from ATCC) and skin primary fibroblasts were cultured in growth medium [Dulbecco’s modified Eagle’s medium (DMEM) with 10% fetal bovine serum (FBS) and 1× antibiotic-antimycotic (Gibco)]. For SAG treatment, cells were seeded in growth medium and 24 h later changed to low serum medium (DMEM with 0.5% FBS and antibiotics) supplemented with 100 nM SAG (Calbiochem) or its vehicle (DMSO) and then maintained for another additional 24 h (NIH3T3 and primary cultures) or 48 h (C3H10T1/2) before being used for experimental procedures. Genotypes of EVC^−/− cells: NM_153717.2: c.940-150T > G/c.1098 + 1G > A and EVC2^−/−cells: NM_147127.4: c.1918delA/c.3660delC (69). Retroviral infections were performed as previously described (27). Briefly, retroviruses were generated by cotransfecting HEK293T by the calcium phosphate method with pBABE-GLI1-expressing constructs or the empty vector and the packaging plasmid pCL-Eco. Filtered supernatants of these cultures were applied to NIH3T3 or C3H10T1/2 cells three consecutive times, and 24 h after the last infection cells were subjected to puromycin selection (2 µg/ml).

Real-time RT-PCR

Total RNA was isolated with Tri Reagent Solution (Ambion) and used for cDNA synthesis with the High Capacity cDNA Reverse Transcription Kit (Applied Biosystems) and random primers. qRT-PCR was performed in a 7900HT Fast Real-Time PCR System (Applied Biosystems) using TaqMan real-time PCR gene expression assays from Applied Biosystems. A minimum of three independent experiments were performed for each analysis and samples were run in triplicates. Transcript levels were normalized against human or mouse GUSB expression. For qRT-PCR in human primary fibroblasts, three extra reference genes were used including ACTB (beta-Actin), GAPDH and HPRT1. Fold differences were calculated by the 2^−ΔΔCt method. Gene expression assays were as follows: Hs00181117_m1 (PTCH1), Hs01110766_m1 (GLI1), Hs99999908_m1 (GUSB), Hs99999903_m1 (ACTB), Hs99999905_m1 (GAPDH), Hs99999909_m1 (HPRT1), Mm00436026_m1 (Ptch1), Mm00494645_m1 (Gli1) and Mm00446953_m1 (GusB). Mm00494645_m1 spans between the first coding exon and the immediate upstream 5′-UTR exon of Gli1 which has no homology to the human GLI1 transcript.

Luciferase reporter assays

HEK293T were plated at a density of 3 × 10⁵/6-well and cotransfected by the calcium phosphate method using 6 µg of a plasmid DNA mixture containing: 12xGliBS-Luc (70) reporter (or 8 × 3′Gli-BS-Luc plasmid when using mouse Gli1 constructs) (71), pRL-TK Renilla luciferase (Promega) and the plasmid expressing the corresponding GLI1 variant at a 4.5:1:4.5 ratio. In competition experiments, two different GLI constructs (GLI1FL: GLI1-W780X; Gli1FL: Gli1-W781X or Gli2ΔN: Gli1-W781X) were added at 1:1 proportion in the transfection mix. Media was replaced 20 h after transfection and cells lysed 24 h later and processed with Dual luciferase reporter assay System (Promega) following manufacturer’s instructions. The reporter experiments were read on a Glomax 96 Microplate Luminometer 8 (Promega) in triplicates and the results were normalized for transfection efficiency using Renilla luciferase values. Three independent transfections were performed for each experiment.

ALP assay

ALP activity was assessed as described previously (72). Retrotransduced C3H10T1/2 cells were seeded at a density of 10⁴/24-well and 24 h later treated with SAG (100 nM) or DMSO for 48 h in low serum medium. Cells were lysed in alkaline buffer (100 mM Tris pH 9.5, 250 mM NaCl, 25 mM MgCl2, 1% Triton X-100) for 45 min and 10 µl of protein extracts were incubated in the dark with 50 µl of the ALP substrate CDP-Star (Roche) during 15 min at room temperature. Three independent retrotransductions were performed for each construct and the result from each retrotransduction was read in triplicate in a Glomax 96 Microplate Luminometer 8 (Promega). Protein concentration determined by the BCA Protein Assay (Pierce) was used to normalize ALP activity. Cytochemical ALP staining was performed using Leukocyte ALP Kit (Sigma) according to the manufacturer specifications.

Western blot

Cell lysis was performed in RIPA buffer (50 mM Tris–HCl pH 6.8, 150 mM NaCl, 2 mM EDTA, 1% NP40, 0.1% SDS, 1% Sodium deoxycholate) containing protease inhibitors (Sigma, P8340) and 1 mM PMSF. Proteins were quantified by the BCA Protein Assay (Pierce) and western blots conducted as indicated earlier (27). Primary antibodies for: human GLI1 (1:1000; Cell Signaling, C68H3), human GLI3 (0.4 µg/ml; R&D, AF3690), α-Tubulin (1:40000; Sigma, T9026), human EVC (1:1000; Sigma, HPA016046), HA (1:1000; Covance, 16B12) and mouse PTCH1 (1:2000) (27). HRP-conjugated secondary antibodies were from Jackson ImmunoResearch. Western blots were developed using ECL reagent (Amersham) and exposed to Agfa X-Ray films.

Immunofluorescence

Cells were plated onto cover slips, maintained 24 h in low serum medium to promote emission of cilia and then fixed on ice for 5 min in 4% PFA/PBS (for human EVC, SMO and acetylated tubulin). Detection of γ-tubulin included an extra incubation step in chilled methanol for 5 min. Cells were permeabilized in PBS/0.1% Triton X-100 for 15 min and blocked in 4% donkey serum/PBS/0.05% Triton X-100 for 30 min. Primary antibodies were incubated overnight at 4 °C and secondary antibodies were added for 1 h at room temperature. Images were acquired on a Nikon 90i microscope. Primary antibodies: acetylated tubulin (1:2000; Sigma, T7451), γ-tubulin (1:2000; Sigma, T5326), human EVC (1:200; Sigma, HPA016046) and SMO (1:1000; Abcam, ab 38686). Fluorescence-conjugated secondary antibodies (1:1000) were from Molecular Probes and nuclei were stained with DAPI (1:2000; Molecular Probes, D1306).

In ovo electroporation

Chick embryos were electroporated at stage 12–13 of Hamburguer and Hamilton classification (stage HH) with purified plasmid DNA (2 μg/μl) as previously described (73) and analyzed 24 h postelectroporation at stage HH18–20.

Immunohistochemistry

Embryos were fixed for 2–4 h at 4 °C with 4% PFA in PBS and then immersed in 30% sucrose solution, embedded in 7.5% gelatin/15% sucrose and sectioned on a cryostat (15 µm; CM1900; Leica). For immunohistochemistry, 15 μm cryostat sections were treated as previously described (73) and incubated overnight at 4 °C with the corresponding primary antibody. Polyclonal primary antibodies against OLIG2 (Millipore: AB9610) and PAX6 (Biolegend: 901301) were used and a monoclonal antibody against NKX2.2 (74.5A5) was obtained from the Developmental Studies Hybridoma Bank (developed under the auspices of the National Institute of Child Health and Development and maintained by the University of Iowa). Alexa Fluor 594- and Alexa Fluor 647-conjugated anti-mouse or -rabbit secondary antibodies (Invitrogen and Life Technologies) were used and nuclei were stained with bisbenzimide. Sections were mounted in Fluoromount (Southern Biotech) and photographed using a confocal scanner microscope (SP5; Leica). For each experimental condition, immunofluorescence analysis was conducted in three sections from at least four embryos.

Data analysis

For protein quantification, films were scanned and subjected to densitometry analysis using ImageJ. Cilia were measured with NIS-Elements software. Statistical analysis was conducted using unpaired Student’s t-test with the help of GraphPad software. P-values ≥ 0.05 were considered statistically non-significant (ns); ***P < 0.001; **P < 0.01; *P < 0.05. Bar graphs represent mean ± SD.

Human samples

This study was conducted in accordance with the declaration of Helsinki for Human Rights and was approved by Institutional Review Boards of Hospital La Paz, Consejo Superior de Investigaciones Científicas (CSIC), Quaid-i-Azam University Islamabad Pakistan and University of Michigan, Ann Arbor, MI, USA. Informed consent to conduct research was obtained from the patients or the patient’s guardians participating in the study.

Supplementary Material

Supplementary Material is available at HMG online.

Ackowledgements

We are grateful to patients and family members for their participation in this study. We thank Dr Karen E Heath, Dr Hiroshi Sasaki and Dr Rune Tofgård for helping with reagents and plasmids and Juan C Triviño for bioinformatics technical assistance. We are also highly grateful to Genoox (https://genoox.com/) and Dr Thomas Kubisiak for their technical support in WGS data analysis. This article resulted in part from a successful GeneMatcher match (74).

Conflict of Interest statement. None declared.

Funding

Spanish Ministry of Economy and Competitiveness (SAF2013-43365-R/SAF2016-75434-R to V.L.R.-P. and BFU2014-57494-R to A.V.M.); the Italian Ministry of Health (Ricerca Corrente 2016 and 2017 to A.D.L.); International Research Support Initiative Program (IRSIP), Higher Education Commission (HEC), Islamabad, Pakistan to A.U.