A Neanderthal/Denisovan GLI3 variant contributes to anatomical variations in mice

Introduction

Understanding the genetic mechanisms underlying phenotypic variation is a fundamental challenge in developmental and evolutionary biology. Several empirical studies suggest that adaptive mutations occurring in cis-regulatory regions are a major force driving morphological evolution (Carroll, 2000; Chan et al., 2010; Werner et al., 2010; Kvon et al., 2016). In contrast, mutations in protein-coding regions have been thought to play little role in anatomical variations, since structural mutations often disrupt protein function and cause severe malformations or pathological conditions in organisms. However, recent studies have revealed that changes in protein structures result in novel pleiotropic functions (Ronshaugen et al., 2002; Protas et al., 2006; Hoekstra and Coyne, 2007). Thus, nucleotide variants that alter amino acid sequences contribute significantly to phenotypic diversity.

Adaptive mutations in coding regions are behind the phenotypic diversifications in human evolution as well (Cargill et al., 1999; Chasman and Adams, 2001; Enard et al., 2002; Asgari et al., 2020). Accumulating data have shown that several admixtures resulted in gene flow between modern human ancestors and archaic hominins such as Neanderthals and Denisovans, who diverged from modern humans approximately 520–630 thousand years ago (Prufer et al., 2014; Prufer et al., 2017). Fossil evidence has indicated distinct anatomical characteristics in Neanderthals compared to modern humans, including elongated and low crania, larger brow ridges, and broader rib cages (Reilly et al., 2022). Putative morphological profiles inferred from DNA methylation patterns predict that Denisovans may have shared some Neanderthal traits (Gokhman et al., 2019). Some of the genetic diversity of extant humans is derived from the extinct hominins (Green et al., 2010; Meyer et al., 2012; Prufer et al., 2014; Sankararaman et al., 2016; Skov et al., 2020). Recent studies have unveiled some of the functional consequences of coding variations specific to archaic or modern humans (Zeberg et al., 2020; Greer et al., 2021; Maricic et al., 2021; Trujillo et al., 2021; Mora-Bermudez et al., 2022; Pinson et al., 2022). However, the impact of archaic genomic variants in anatomical structures has not been directly demonstrated.

Previous studies have revealed that relative to modern humans, all sequenced Neanderthals and Denisovans carry an amino acid substitution in the C-terminal region of GLI3 (R1537C) (Castellano et al., 2014; Prufer et al., 2014). GLI3 is a GLI-Krüppel family transcription factor that plays essential roles in embryonic development by mediating Hedgehog signaling (Matissek and Elsawa, 2020). In the absence of Hedgehog, GLI3 is cleaved to form the N-terminal repressor form, while the presence of Hedgehog prevents protein processing, which triggers intracellular accumulation of full-length GLI3. Functional disruptions of GLI3 normally result in severe abnormalities in various anatomical structures in mice and humans (Vortkamp et al., 1991; Kang et al., 1997; Naruse et al., 2010). R1537C variant of GLI3 corresponds to a single nucleotide polymorphism (SNP) that exists in modern human populations (rs35364414, 4609G>A results in arginine to cysteine substitution). However, the consequences of the amino acid substitution of the GLI3 protein and phenotypic outcomes have not been clarified.

Here, we investigated the biological effects of the Neanderthal/Denisovan GLI3 variant (R1537C). Our study shows that R1537C substitution affects the regulation of downstream target genes associated with developmental processes without interfering Hedgehog-dependent transcriptional activities. Genome-edited mice carrying the Neanderthal/Denisovan GLI3 mutation exhibited various alterations in skeletal morphology depending on genetic backgrounds. The data suggest that the extinct hominin-type GLI3 variant affects species-specific anatomical variations, which were accommodated by concurrent changes in developmental constraints during human evolution.

Results

The R1537C variant might be accommodated in archaic and modern humans

R1537C is localized at C-terminal transactivation domain 1 (TA1) of human GLI3 (Figure 1A). At this position, arginine (R) is extensively conserved among mammals and a bird (Figure 1B), suggesting that the corresponding cysteine residue is a derived mutation in Neanderthal/Denisovan lineages. In silico analysis of protein functions predicted that R1537C is deleterious to human GLI3 (Figure 1B). Of note, the SIFT program suggested that substitution to histidine was tolerable at this position, which corresponds to the altered amino acid in Koala Gli3 (Figure 1B). These data suggest that relaxed constraint on the residue occurred in the evolution of archaic human lineages. Next, we investigated the R1537C variant exists in modern human population. Present-day human genomes outside Africa contain 1%–4% DNA derived from Neanderthals, while approximately 6% of the Melanesian genomes are derived from Denisovans (Green et al., 2010; Reich et al., 2010). We checked geographic distribution of the SNP corresponding to the R1537C variant in the modern human population [rs35364414 (G/A)]. The allele frequency of the R1537C variant (A allele) ranges from highest in the European population (3.7%–7.7%) to lowest in the African population (0.8%–2.1%) (The 1000G genome dataset, Figure 1C). A previous study suggested that the chromosomal fragment containing the SNP was inferred to be derived from archaic hominins (Skov et al., 2020). These lines of evidence suggest that the R1537C variant is accommodated in genetic background of archaic and present-day humans.

FIGURE 1

FIGURE 1. GLI3^R1537C is a derived mutation in archaic and present-day humans (A) Protein structures of the modern human GLI3 (GLI3^MH) and a Neanderthal and Denisovan GLI3 variant (GLI3^R1537C). RD, repressor domain; ZnF, Zinc finger type DNA binding domain; TA1, transactivation domain 1; TA2, transactivation domain 2. (B) Upper panel: amino acid sequences of GLI3/Gli3 C-terminal regions in mammals and a bird (chick). A red arrowhead indicates the corresponding arginine position. Lower panel: PROVEAN score and SIFT prediction of GLI3^R1537C. (C) Allele frequencies of R1537C (rs35364414; A allele) in 1000G populations. ACB, African Caribbean in Barbados; ASW, African ancestry in Southwest US; CLM, Colombian in Medellin, Columbia; CEU, Utah residents with Northern and Western ancestry; FIN, Finnish in Finland; GBR, British in England and Scotland; IBS, Iberian population in Spain; TSI, Toscani in Italy; BEB, Bengali in Bangladesh; GIH, Gujarati Indian in Houston; ITU, Indian Telugu in the UK; PJL, Punjabi in Lahore, Pakistan.

The R1537C variant does not compromise protein stability and Hedgehog-dependent GLI3 functions

Next, we investigated whether R1537C substitution affects GLI3 protein expression. To test this, we overexpressed modern human-type GLI3 (GLI3^MH) or a Neanderthal/Denisovan GLI3 (GLI3^R1537C) in human embryonic kidney cells (HEK293T cells) and performed Western blotting. The amounts of full-length and truncated forms of GLI3 (GLI3 full and GLI3R, respectively) were comparable between GLI3^MH and GLI3^R1537C (Figures 2A, B). Specificity of the GLI3 antibody was confirmed by Western blotting with GLI3 KO cell lines (Supplementary Figures S1A–C). To further examine the stability of GLI3^MH and GLI3^R1537C, we overexpressed luciferase tagged GLI3^MH and GLI3^R1537C and measured luciferase activity after halting protein expression. We confirmed that the degradation rate was not different between of GLI3^MH and GLI3^R1537C (Figure 2C). Thus, in contrast to in silico prediction, R1537C substitution does not disturb the protein stability of GLI3.

FIGURE 2

FIGURE 2. R1537C does not compromise GLI3 protein stability (A) Western blotting of GLI3^MH and GLI3^R1537C. (B) Quantification of GLI3FL (full length, upper panel) and GLI3FL/GLI3R ratio (lower panel). (C) Comparison of the degradation speed of nanoluc (NlucP)-tagged GLI3^MH and GLI3^R1537C after halting protein synthesis by doxycycline administration.

Next, we test whether R1537C substitution affect GLI-mediated Hedgehog signaling. Previous studies reported that full-length GLI3 predominantly acts as an inhibitor of Hedgehog signaling, as GLI3 negatively regulates GLI1-dependent transcriptional activity (Wang et al., 2000). Thus, we examined the effect of R1537C on the suppressive potential of GLI3 in the presence of GLI1. To test this, we overexpressed a GLI reporter vector together with GLI1 and GLI3^MH or GLI3^R1537C in HEK293T cells (Figure 3A). Compared to the control, GLI1 strongly enhanced reporter activity, while co-transfection of GLI3^MH or GLI3^R1537C significantly decreased GLI1-dependent reporter activity, indicating that the repressive activity of GLI3 is not curtailed by the amino acid substitution.

FIGURE 3

FIGURE 3. R1537C does not interfere repressor or activator-dependent transcriptional activities (A) GLI-reporter activities in HEK293T cells overexpressing GLI3^MH or GLI3^R1537C together with GLI1. (B) Protein structures of full-length and repressor domain-truncated GLI3 (GLI3^MHΔRD or GLI3^R1537CΔRD). (C) GLI-reporter activities in HEK293T cells overexpressing GLI3^MHΔRD or GLI3^R1537CΔRD. (D) Induction of Olig2 expression in the developing chick neural tube by overexpression of GLI3^MHΔRD or GLI3^R1537CΔRD. Scale bar: 50 µm (D).

Full-length GLI3 contains both transcriptional repressor and activator domains. To directly assess the effect of R1537C on the transcriptional activator domain, we removed the N-terminal repressor domain and examined luciferase activity in response to GLI3 variants (GLI3^MHΔRD and GLI3^R1537CΔRD, Figure 3B). Compared to the control, both GLI3^MHΔRD and GLI3^R1537CΔRD significantly increased reporter activity at comparable levels (Figure 3C). Furthermore, we corroborated that overexpression of GLI3^MHΔRD and GLI3^R1537CΔRD induced the expression of Olig2, a direct downstream target of GLI proteins, in the developing chick neural tube (Figure 3D). These data indicate that R1537C does not affect the function of transcriptional activator domain of GLI3.

R1537C affects unique downstream gene expression related to anatomical structures

Our luciferase assay suggest that R1537C preserves transcriptional activity of GLI3. To elucidate the impact of R1537C on GLI3-dependent unique target genes, we overexpressed full-length GLI3^MH or GLI3^R1537C in HEK293T cells and performed RNA sequencing (RNA-seq) analysis. Because GLI3 is weakly expressed in HEK293T cells, we also performed RNA-seq with control HEK293T cells transfected with empty vector to capture background level of gene expression. First, we checked the expression of direct downstream genes of Hedgehog signaling, such as NKX2.2, NKX6.1, PAX6, IRX3, and OLIG2. Compared to the control, these genes were slightly up- or downregulated by full-length GLI3 overexpression. Importantly, we could not detect differences in Fold change (Fc) of these genes between GLI3^MH and GLI3^R1537C overexpression (Figure 4A), corroborating that R1537C does not interfere with Hedgehog-dependent signaling. Next, we focused on differentially expressed genes [DEGs; |Fc|≥1.5, false discovery rate (FDR) < 0.05] between the control and GLI3^MH or GLI3^R1537C overexpression. Compared to the control HEK293T cells, a total of 846 and 671 genes were differentially expressed by the overexpression of GLI3^MH or GLI3^R1537C, respectively (Supplementary Table S1). We further conducted Gene Ontology (GO) enrichment analysis to identify the functional characteristics of DEGs. Compared to control, GLI3^MH overexpression enriched GO terms related to various biological processes, consistent with the pleiotropic functions of GLI3. Specifically, GLI3^MH overexpression enriches GO terms associated with developmental process, such as ossification, renal system development, and extracellular matrix organization (Figure 4B). In contrast, GLI3^R1537C did not enriched these GO terms, but increased terms related to chromatin and nucleosome assembly (Figure 4B). To further investigate altered gene regulations by R1537C, we directly compared transcriptomes with GLI3^MH and GLI3^R1537C overexpression and identified 91 DEGs that were significantly up- or downregulated by GLI3^R1537C (DEGs in GLI3^R1537C vs. GLI3^MH; |Fc| ≥ 1.5, FDR < 0.05, Supplementary Table S2). These DEGs included unique transcripts, such as LINC00294, a human specific noncoding RNA that is involved in cell proliferation (Qiu et al., 2020; Zho et al., 2020; Dong et al., 2021), and SLC9A3 (NHE3), a member of the sodium-hydrogen exchanger family that is responsible for the maintenance of sodium balance (Orlowski and Grinstein, 1997). Among them, 23 genes were associated with developmental process (Figure 4C). Notably, H4 clustered histone 3 (H4C3; Fc = 3.7935, p = 0.0465) plays an essential role in chromatin organization and disruption of this gene results in growth retardation and skeletal abnormalities in humans (Tessadori et al., 2022). Stanniocalcin 1 (STC1; Fc = −2.93183, p = 7.15551E-05) is a glycoprotein involved in calcium and phosphate metabolism, and STC1 transgenic mice exhibited dwarf phenotype due to altered osteogenesis (Filvaroff et al., 2002). These data suggest that the missense mutation of GLI3 alters transcriptional regulations related to specific organ development during embryogenesis.

FIGURE 4

FIGURE 4. R1537C alters the regulation of GLI3-dependent gene expression (A) Fold changes of direct downstream genes of Hedgehog signaling in HEK293T cells overexpressing GLI3^MH or GLI3^R1537C compared with control samples. (B) GO enrichment analysis in comparisons of HEK293T cells overexpressing GLI3^MH or GLI3^R1537C versus control samples. Spots represent the top-ranking terms with gene ratios over 0.02. (C) Genes associated with developmental process that were down- or upregulated by GLI3^R1537C compared with GLI3^MH overexpression.

Mice with a Neanderthal/Denisovan-type Gli3 variant exhibit altered skeletal morphology

To further characterize the effect of R1537C on embryonic development, we created knock-in mice with the Neanderthal/Denisovan GLI3 variant by CRISPR-mediated genome engineering technology. In mice, the amino acid that corresponds to the missense variant is naturally arginine (1540R) (Figures 1A, 5A). Thus, we introduced the point mutation c4618t by Cas12 (Cpf1)-dependent homologous recombination, resulting in an R1540C substitution in the mouse GLI3 protein (Figure 5B, Supplementary Figures S2A–C). Western blotting in HEK293T cells confirmed that the amino acid substitution did not compromise protein expression compared to wild type GLI3 (1540R) (Supplementary Figure S2D). On the C57BL6 background, mice carrying 1540C were born at the expected Mendelian ratio. Both full-length and repressor form of Gli3 were detected in homozygous 1540C/1540C mice as well as heterozygous (1540R/1540C) and wild-type mice (Figure 5C), corroborating that protein processing is not impaired in the Gli3 variant. Homozygous 1540C mice had similar body sizes to wild type or heterozygous (1540R/1540C) littermates (Supplementary Figure S2E). However, we observed altered head morphology in 1540C homozygous mice compared with the wild type or heterozygous littermates. Whole mount skeletal analysis indicated enlarged crania with frontal and parietal bossing in the 1540C mice (Figures 5D, G and Supplementary Table S3). Previous studies reported that Gli3 null mice (Gli3^xt/xt) exhibited craniosynostosis, a premature closure of cranial sutures, which frequently results in abnormal shape of skulls (Veistinen et al., 2012). Consistently, deformed skulls in 1540C mice were frequently associated with premature ossification at the coronal suture (Figures 5E, F).

FIGURE 5

FIGURE 5. Altered skeletal morphology of the mice carrying a Neanderthal/Denisovan type GLI3 variant (A) Structures of modern human or extinct hominins GLI3, wild type mouse Gli3 (Gli3^1540R), and mouse Gli3 with a Neanderthal/Denisovans variant (Gli3^1540C). (B) Single nucleotide substitution in mouse Gli3 by CRISPR-mediated genome editing. (C) Western blotting of wild-type (WT), heterozygous (1540R/1540C) and homozygous (1540C/1540C) mouse brains (E17.5 cerebellar primordium) with anti-GLI3 antibody. (D–F) The skull of heterozygous (1540R/1540C) and homozygous (1540C/1540C) mice on C57BL6 background. An arrowhead indicates an enlarged cranium in 1540C mouse. (E,F) Coronal sutures (arrowheads) of heterozygous (E) and homozygous (F) mice. (G) Quantification of skull height and length. (H) Comparison of rib cages. Arrowheads indicate abnormal rib torsion in 1540C/1540C mice. (I) Reduced number of lumber vertebrae in 1540C/1540C mouse. (J–M) Skeletal morphology of mice on CD-1 background. The 1540C/1540C mice exhibit a twist of the posterior angle of the rib [arrowheads in (J,K)], hypoplastic ossification in the 2nd cervical vertebrae [arrowheads in (L)], extra thoracic ribs [T14, arrowheads in (M)] and reduced number of lumber vertebrae (M). Scale bars: 1 mm.

Some homozygous 1540C mice exhibited asymmetric shapes of rib cages associated with scoliosis (Figure 5H and Supplementary Table S3). Furthermore, compared to wild type mice, the number of lumbar vertebrae was reduced in heterozygous and homozygous mice carrying the 1540C allele (Figure 5I and Supplementary Table S3). Polydactyly was shown to be a typical phenotype caused by Gli3/GLI3 null conditions in both mice and humans (Naruse et al., 2010). In contrast, the heterozygous and homozygous 1540C mice did not exhibit any abnormalities in limb and digit formation (Supplementary Figure S2E). These results indicated that a Neanderthal/Denisovan Gli3 variant resulted in altered skeletal morphology, which differed from the phenotypes associated with the loss of Gli3 functions.

It has been shown that calvarial ossification is mediated by Indian hedgehog (Ihh)-dependent signaling (Veistinen et al., 2017). To ask whether ossification phenotypes in 1540C mice is due to altered Ihh signaling, we examined the expression of Ptch1, a Hedgehog responsive gene, in the developing cartilage. We confirmed that Ptch1 expression in the developing calvarial or thoracic mesenchymal tissues was comparable between wild type and 1540C homozygous mice (Supplementary Figures S3A–E), suggesting that altered ossification in 1540C mice was not due to dysregulation of Ihh signaling.

Several lines of studies on archaic human genomes indicated that Neanderthals and Denisovans had distinct genetic backgrounds compared to modern humans, although these extinct hominins and ancestral modern humans had frequently intermixed with each other (Green et al., 2010; Meyer et al., 2012; Prufer et al., 2014; Reilly et al., 2022). Accordingly, the phenotypic consequences caused by the missense variant might be variable depending on specific genetic backgrounds. To examine this possibility, we introduced the 1540C substitution into mice with a CD-1 (ICR) background, which is an outbred strain with heterogeneity of individual genomic sequences compared to inbred strains. Homozygous 1540C mice were viable and their body sizes were similar to those of wild-type or heterozygous knock-in mice (Supplementary Figures S2F–H). As in the case of the knock-in mice on the C57BL6 background, CD-1 mice carrying the 1540C allele exhibited abnormal rib torsion (Figures 5J, K). In particular, the posterior angle of the ribs in 1540C mice was abnormally kinked, which resulted in the crank-shaft sternum (Figure 5K, Supplementary Figure S2I and Supplementary Table S3). Furthermore, we found additional phenotypes in CD-1 background mice, such as hypoplasia of the second cervical vertebrae and extra ribs at the 14th thoracic vertebrae (Figures 5L, M). The number of lumber vertebrae was reduced in the mice with extra ribs (Figure 5M). These skeletal abnormalities were also observed in a lower proportion of the heterozygous mice, indicating a dose-dependency of the phenotypes correlated with the presence of the 1540C allele (Supplementary Table S3). Curiously, extra 14th ribs occurred in a small proportion of the wild-type CD-1 mice (Supplementary Table S3), suggesting that the missense Gli3 variant increased phenotypic variations sustained in the background population. In contrast to C57BL6 background, cranial abnormalities were not detected in the homozygous mice on CD-1 background, suggesting that these phenotypes are highly affected by genetic background (Figure 6A).

FIGURE 6

FIGURE 6. The arginine to cysteine substation in GLI3/Gli3 affects common anatomical structures in human and mice (A) Variation of skeletal phenotypes in Gli3^R1540C mice on different genetic background. The proportion represents phenotypes in homozygous mice. (B) HPO terms shared with Neanderthal (Gokhman et al., 2020) and Gli3^R1540C mouse anatomical characteristics. GLI3 related clinical features in modern humans are also represented. Images of Neanderthal, mouse and human are from PHYLOPIC.

Similarities in human and mouse transcriptomes depending on the GLI3/Gli3 variants

To identify downstream genes affected by the missense Gli3 variant in mice, we isolated mesenchymal tissues of wild-type and 1540C mouse embryos on the C57BL6 background and performed RNA-seq analysis. A total of 453 and 238 genes were up- or downregulated in 1540C compared to wild-type mice [|Fc| ≥ 1.5, false discovery rate (FDR) < 0.05, Supplementary Table S4]. We found that DEGs included Msx1, Dlx1 and Dlx2, which were known to be essential for skeletal morphology (Alappat et al., 2003; Levi et al., 2022). Notably, several histone H4 cluster genes, such as H4c2, H4c6 and H4c14 were upregulated in 1540C mice, suggesting that the arginine to cysteine substitution of Gli3 commonly affects specific histone components in mice and humans (Supplementary Table S4). To further address the commonalities of downstream target genes affected by the GLI3^R1537C and Gli3^R1540C, we compared the GO terms of HEK293T cells transfected with human GLI3 variants (Figure 4) and Gli3 knock-in mice. Although these data sets were derived from distinct species and cell types, the terms for cellular components showed remarkable similarity to the top-ranking enriched terms, including CENP-A-containing chromatin and CENP-A-containing nucleosomes (Supplementary Figures S4A, B). These results suggest that the missense variant of GLI3/Gli3 disturbed common downstream targets associated with nucleosome and chromatin assembly between human and mice.

The arginine to cysteine substation in GLI3/Gli3 affects common anatomical structures in human and mice

We demonstrate that mice with a Neanderthal/Denisovan GLI3 variant exhibit altered skeletal structures, such as enlarged cranium, altered shapes of vertebrae, and rib malformations (Figure 6A). As shown in Figure 6A, enlarged skull phenotype showed higher penetrance in C57BL6 mice, while such phenotype was not observed in CD-1 background mice. In contrast, hypoplastic vertebrae and extra ribs frequently occurred in CD-1 background mice, despite these phenotypes represented lower penetrance.

A recent study reported the association of unique anatomical phenotypes in Neanderthal with Human Phenotypic Ontology (HPO), a database of disease-associated human phenotypes (Gokhman et al., 2020). Notably, some of neanderthal-specific HPO phenotypes, such as macrocephaly (HP: 0000256), accelerated skeletal maturation (HP: 0005616) and kyphoscoliosis (HP: 0002751) are linked with the phenotypes of R1540C mice (Figure 6B and Supplementary Table S5). We find that the mice carrying the Neanderthal/Denisovan GLI3 variant also showed altered shapes of rib cages due to abnormal rib torsion (Figures 5H, J, K, 6A). Fossil evidence indicated a stronger torsion of central ribs (6th–8th ribs) in Neanderthal infants, which contributes to distinct rib cages (shorter and deeper forms) compared to those of modern humans (Garcia-Martinez et al., 2020). Even if not all skeletal phenotypes of R1540C mice phenocopied anatomical characteristics of extinct humans, phenotypic alterations in the knock-in mice suggest that the corresponding missense mutation of GLI3 might affect conserved developmental processes of homologous skeletal structures between mice and humans.

Discussion

Here, we show that an extinct hominin-specific amino acid alters GLI3 functions by conferring differential transcriptional regulation of downstream genes and results in altered skeletal morphology in mice. Several studies have shown that the ratio of GLI3 activator versus repressor is critical for mediating Hedgehog-dependent signaling during embryogenesis (Wang et al., 2007; Wen et al., 2010). Despite in silico prediction, we demonstrate that R1537C substitution does not compromise GLI3 protein stability or activator-dependent transcriptional activity. Furthermore, R1537C does not alter the regulations of direct target genes of Sonic hedgehog signaling. Thus, GLI3^R1537C does not interfere with the potential to mediate Hedgehog signaling during embryogenesis, thereby detouring deleterious phenotypes. In contrast, disruption of the GLI3 C-terminal region causes the pathogenic phenotypes in human (Demurger et al., 2015), suggesting that the C-terminal region is a prerequisite for GLI3 functions, possibly through interaction with various cofactors for transcriptional regulation (Zhou et al., 2006).

We also demonstrate that GLI3^R1537C affects the transcription of various nucleosome components, providing first evidence on the association of C-terminus region of GLI3 on chromatin assembly. Furthermore, transcriptome analysis of GLI3^R1540C knock-in mice revealed evolutionary conservation of these regulatory pathways between human and mice. It has been shown that de novo missense variants of histone H4 clusters results in various neurodevelopmental and skeletal abnormalities (Tessadori et al., 2022). Notably, a part of these clinical features, such as craniosynostosis, hypertelorism and syndactyly, are also reported in patients with GLI3 variants (Vortkamp et al., 1991; Kang et al., 1997; Demurger et al., 2015). Although empirical data explaining the relationship between these pathways and developmental phenotypes are still limited, these lines of evidence suggest that phenotypic alterations by the GLI3 variant is mediated by global consequences of chromatin organization.

The corresponding amino acid substitution is highly unique in archaic hominins and modern humans. Various skeletal abnormalities induced by archaic hominin-type substitution might not be adaptive in mice. These lines of evidence suggest that concurrent changes in developmental programs occurring in archaic hominins accommodated R1537C-dependent phenotypic abrogation. Alternatively, the variant-dependent functional alterations might be tolerated by additional compensatory mutations that occurred in archaic and modern human lineages. R1537C-dependent phenotypes could be attributed to the interaction with other variants reside in archaic and modern human-specific genetic backgrounds. We confirmed that differential genetic background significantly affects Gli3 variant-dependent phenotypic variations in mice. Phenotypic differences might be due to genetic and/or epigenetic differences in other genomic regions, which affect Gli3-dependent developmental processes. Notably, the R1540C variant frequently contributed to 14th rib formation in CD-1 background, which also occurs occasionally in wild type population. It has been reported that variable number of vertebrae depends on genetic backgrounds as well as maternal influences (Mc and Michie, 1958; Rengasamy and Padmanabhan, 2004). A recent study reported that the R1537C variant is associated with neural tube defects in European population (Renard et al., 2019). Thus, it is possible that the GLI3^R1537C/Gli3^R1540C increase susceptibility to variable morphological characteristics or pathogenic phenotypes in response to genetic and environmental perturbations.

It remains unclear if R1537C variant was fixed in extinct hominin groups by positive selection. Recent studies have suggested that the small population size of Neanderthals was not effective for natural selection, by which weakly deleterious mutations might be accumulated by genetic drift (Prufer et al., 2014; Prufer et al., 2017). Notably, datasets of associations between various human traits and variants based on UK Biobank cohorts (http://geneatlas.roslin.ed.ac.uk, https://pheweb.org/UKB-SAIGE/) revealed that R1537C (A allele) is associated with specific anatomical traits, such as anisotropy of the superior longitudinal fasciculus (p = 0.00038) waist-hip ratio (p = 0.00057), and other unspecified back disorders (p = 0.0000045), while the modern human variant (G allele) increases the risk of lethargy (p = 0.00058) and osteoarthritis (p = 0.00067). Although most of these associations did not pass statistical correction for multiple comparisons, these traits are linked with the predicted lifestyles of Neanderthals, suggesting that R1537C provided beneficial traits for extinct hominins. To distinguish these possibilities, further analyses are required to identify the biological effects of the variant in human development, as well as genetic and epigenetic modifiers that affect variant-dependent phenotypes.

Materials and methods

Animals

Male and female mice (CD-1 and C57BL6/J backgrounds) that were originally obtained from Japan SLC and Charles River Laboratories Japan were maintained in a 12 h dark/light cycle at the Experimental Animal Facility of Kyoto Prefectural University of Medicine. Noon of the day the vaginal plug was identified was designated as E0.5. Fertilized chicken eggs (Gallus gallus) were obtained from a local farm (Yamagishi, Japan) and incubated at 37.0°C ± 0.2°C. Embryonic stages were determined according to the Hamburger-Hamilton stages (Hamburger and Hamilton, 1951). All animal experiments were approved by the experimental animal committee of Kyoto Prefectural University of Medicine and were performed in accordance with the relevant guidelines (M2021–233, M2022-209, M2021-511, M2022-210, M2021-217, M2022-201).

Plasmids

Full-length human GLI3 and mouse Gli3 cDNAs were obtained from pEGFPC3-hGli3 [a gift from Aimin Liu, Addgene plasmid (Zeng et al., 2010)] and pCMV-Gli3-Myc-DDK (ORIGENE, MR226591), and then subcloned into the pCAG-RB vector by using an In-Fusion HD cloning kit. For construction of pCAG-GLI3^R1537C and pCAG-mGli3^R1540C, the 3’ regions of GLI3 or mGli3 were replaced with PCR-amplified fragments containing c4609t (GLI3^R1537C) or c4618t (mGli3^R1540C), respectively. pCAG-GLI3^MHΔRD and pCAG-GLI3^R1537CΔRD were generated by subcloning the C-terminal region of GLI3^MH or GLI3^R1537C into the pCAGRB vector, according to the method of a previous report (Stamataki et al., 2005). For construction of pTRE-GLI3^MH-NLucP and pTRE-GLI3^R1537C-NLucP, coding sequences of GLI3^MH or GLI3^R1537C were subcloned into pTRE3G-NlucP [a gift from Masaharu Somiya, addgene plasmid (Somiya and Kuroda, 2021)]. All vectors were verified by sequencing.

Luciferase reporter assay

Gli-reporter NIH3T3 cells were transfected with pCAG-GLI3^MH or pCAG-GLI3^R1537C together with pTK-RL by using Lipofectamine 2000 (Thermo Fisher Scientific). To monitor Gli-dependent luciferase activity in HEK293T cells, p8xGli-BS-delta51-LucII (Sasaki et al., 1997) was transfected with other expression vectors. As a control condition, the pCAGRB empty vector was cotransfected with reporter vectors. For analysis of GLI1-dependent reporter activities, pcDNA3.1-His-Gli1 (Sasaki et al., 1999) was transfected with GLI3 expression vectors. Transfected cells were cultured in Dulbecco’s modified Eagle’s medium (high glucose, Nacalai Tesque) containing 10% fetal bovine serum and 1% penicillin and streptomycin (Fujifilm Wako Pure Chemical Corporation) for 24 h. For analysis of the degradation rates of GLI3^MH and GLI3^R1537C, pTRE-GLI3^MH-NLucP or pTRE-GLI3^R1537C-NLucP was transfected into the HEK293 Tet-Off Advanced cell line, and doxycycline was then added to the culture medium 12 h after transfection. Luciferase activity was examined with the Dual-Luciferase Reporter Assay System and analyzed with a luminometer (GENE LIGHT GL210A, Microtec, Inc.). All firefly luciferase values were normalized to Renilla luciferase activities. At least three biologically independent samples were analyzed in each experimental condition, and all experiments were confirmed by at least three technical replicates.

GLI3 KO cell line

A GLI3 knockout cell line (HZGHC007467c021, Horizon Genomics) was generated by manufacture’s protocols. Briefly, 166 bp of bovine minisatellite genomic sequence was inserted into GLI3 genomic locus of HAP1 cells by CRISPR/Cas9 system. A frame-shift mutation with ectopic stop codons in exon 5 was confirmed by sanger sequencing. Parental cell line (WT HAP1) carries modern human type GLI3 variant (1540R). Sequence information of the cell line is available on Mendeley data.

Western blotting

HEK293T cells transfected with pCAG-GLI3^MH, pCAG-GLI3^R1537C, pCAG-mGli3, or pCAG-mGli3^R1540C were lysed in RIPA buffer (20 mM Tris, 150 mM NaCl, 1 mM EDTA, 1% Nonidet P40, 1% SDS, 0.1% deoxycholate, 1 mM NaF and protein inhibitor). GLI3 KO or WT HAP1 cells were directly lysed in RIPA buffer. and After electrophoresis, proteins were transferred to Poly Vinylidene Di-Fluoride (PVDF) membranes, blocked with Bullet Blocking One for Western blotting (Nacalai Tesque), and then incubated with anti-human GLI3 or anti-α-tubulin antibody. The membranes were further incubated with biotinylated anti-rabbit IgG or biotinylated anti-rat IgG antibody, followed by a Vectastain Elite ABC kit and developed with Chemi-Lumi One Super and analyzed with a luminescent image analyzer (LAS-2000, Fujifilm). For Western blotting with mouse brains, the cerebellar primordium of E17.5 mice were manually dissected and lysed in RIPA buffer.

In ovo electroporation

In ovo electroporation of developing chick embryos was performed according to a previous study (Stamataki et al., 2005). Briefly, a window was opened in the shell of an egg, and then a small amount of DNA solution (less than 0.05 µL) was injected into the neural tube of E2 (HH stage 13–15) embryos with a fine glass needle. Needle-type electrodes (CUY200S, BEX) were placed on the neural tube and square electric pulses (28 V for 50 ms, 2 or 3 times) were applied to the thoracic neural tube with a pulse generator (CUY21 EDITII, BEX). After electroporation, the extraembryonic cavity was filled with sterilized Hank’s balanced salt solution (HBSS, Nacalai Tesque) containing antibiotics (penicillin/streptomycin and gentamycin), and the window was sealed with parafilm. Electroporated embryos were incubated in an incubator at 37°C for 24–48 h.

Immunohistochemistry

Electroporated embryos were fixed with 4% paraformaldehyde dissolved in phosphate-buffered saline (PBS) at 4°C overnight. After PBS washes, the embryos were cryoprotected with a 20% sucrose solution and immersed in Tissue-Tek. The frozen samples were sectioned at a thickness of 20 µm using a cryostat (Leica CM 1850, Germany), and incubated with primary antibodies, including anti-Olig2 and anti-GFP antibodies. After washing, the sections were incubated with secondary antibodies, including Alexa-Fluor 488 or 594 conjugated anti-rat or anti-rabbit antibodies. Fluorescence images were captured with fluorescence microscopes (BX51, Olympus) equipped with a cooled CCD camera (DP80, Olympus). All captured images were processed with cell Sense standard (v.1.17, Olympus), ImageJ and Adobe Photoshop.

RNA-seq analyses

Total RNA was prepared from HEK293T cells transfected with pCAG-RB, pCAG-GLI3^MH or pCAG-GLI3^R1537C by using an RNeasy kit. To isolate total RNA from mouse tissue, thoracic mesenchymal tissue of E17.5 wild-type and 1540C mice was dissected and preserved in RNAlater until RNA extraction. Residual DNA was eliminated by DNase treatment. The quality and quantity of RNAs were assessed by using an Agilent Technologies 2100 Bioanalyzer or 2200 TapeStation (Agilent). The cDNA library was constructed using a TrueSeq Standard mRNA LT Sample Prep Kit according to the manufacturer’s protocol (15031047 Rev. E) and sequenced on an Illumina platform and 101 bp paired-end reads were generated. The sequence data were mapped to a reference genome sequence (Homo sapiens GRCh38; NCBI_109.20200522 or Mus musculus GRCm38, GCA_000001635.2) with a splice-aware aligner (HISAT2 v 2.1.0). The transcripts were assembled by StringTie (v 2.1.3b) with aligned reads. GO enrichment analysis was performed based on Gene Ontology. A significant gene list was constructed by the g:Profiler tool.

SNP population genetics

The genetic population of rs35364414 (GLI3 4609G/A) was analyzed by Ensemble genome browser 108 (https://www.ensembl.org/index.html).

Mouse genome engineering

For generation of a point mutation in mouse Gli3 that corresponds to a Neanderthal GLI3 variant, genome edited mouse lines were established by CRISPR/Cas12a (Cpf1)-mediated homology directed repair (HDR). A crRNA sequence that targets mouse Gli3 4618c was designed by CHOPCHOP. We chose a target sequence with no off targets based on mouse genome information (mm10/GRCm38). As an HDR template for homologous recombination, single strand DNA (ssDNA) that carries a nucleotide substation (4618t) was designed by using the mouse Gli3 sequence.

Generation of genome edited mice

The cumulus complexes were isolated from superovulated female C57BL6/J mice treated with pregnant mare serum gonadotropin (PMSG) and human chronic gonadotropin (hCG), and then mixed with sperm taken from the epididymis of male C57BL6/J mice. After in vitro fertilization, the CRISPR/Cas12a solution containing Cpf1 protein (Integrated DNA Technology, Inc.), crRNA, and ssDNA were introduced into zygotes by electroporation (Genome Editor Plus. BEX, Tokyo, Japan). After electroporation, the zygotes were transferred into the oviduct of pseudopregnant CD-1 mice. Anesthesia was performed by intraperitoneal injection of medetomidine, midazolam and butorphanol. Genome-edited mice on the CD-1 background were generated by the i-GONAD method (Gurumurthy et al., 2019). Briefly, pregnant female mice (E0.75) were anesthetized with 2% isoflurane and the CRISPR/Cas12a solution was injected into the oviduct lumen. Then, square electric pulses were applied to the oviduct by using a pulse generator (super electroporator NEPA21, Nepa gene).

Genotyping of mice

Genomic DNA from embryonic or adult mouse tissue (skin or tail fragments) was prepared using a DNeasy Blood and Tissue kit. The DNA fragments containing the CRISPR/Cas12a target region were amplified by polymerase chain reaction (PCR). PCR fragments were digested by BsmBI for 1 h, and the restriction fragment length was determined by electrophoresis with 3% agarose. For determination of the edited sequences of candidate mouse lines, PCR amplicons were subcloned into the pBluescriptSK vector, and six to ten randomly selected clones were examined by Sanger sequencing. F0 mice heterozygous for the 4618c > t substitution (Gli3^R1537C) were propagated and used to generate homozygous mice.

Skeletal analyses

Embryonic and postnatal mice were collected (embryonic Day 16.5, postnatal Day 1, or Day 7) and dissected in phosphate buffered saline. The samples were fixed with 100% ethanol for 3 days and degreased by acetone overnight. Skeletal staining was performed by Alcian blue and Alizarin red according to a standard protocol (Rigueur and Lyons, 2014). Sample images were captured by using the microscope (SZX7, Olympus) equipped with a cooled CCD camera (DP80, Olympus).

Image processing and analyses

Images of Western blotting were examined by using ImageJ. All captured images were processed with Adobe Photoshop 2021. Original Western blotting images were deposited in Mendeley Data.

Statistical analyses

For statistical analysis, at least three independent samples from each experimental group were compared. Comparisons between experimental groups were performed using Microsoft Excel and Prism 9. All data are presented as the mean ± SE. Statistical significance was determined using the two-tailed unpaired Student’s t-test, ordinary one-way ANOVA with Tukey’s multiple comparisons test, or two-way ANOVA with Sidak’s multiple comparison test. Statistical analysis of differential gene expression conducted by RNAseq data was performed using fold change and the nbinomWald Test using the DESeq2 R package (Bioconductor) per comparison pair, and significant results were selected based on |fc| ≥ 1.5 and nbinomWald Test with a raw p-value < 0.05.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ddbj.nig.ac.jp/, DRA015189, https://www.ddbj.nig.ac.jp/, DRA015774, https://data.mendeley.com/datasets/mwdwkddfh8, 10.17632/mwdwkddfh8.2. Key resources used in this study is listed in Table 1.

TABLE 1

TABLE 1. Key resources used in this research.

Ethics statement

The animal study was approved by the Kyoto Prefectural University of Medicine Ethics Committee for Animal Experiments. The study was cunducted in accordance with the local legistration and institutional requirements.

Author contributions

AA and TN contributed to the experimental design, data collection, interpretation, and manuscript preparation. SO contributed to the generation and maintenance of mice carrying GLI3^R1540C. RN contributed to the collection and analysis of RNA-seq samples, HG and KO contributed to the acquisition of the reagents and interpretation. All authors contributed to the article and approved the submitted version.

Funding

This work was supported by Japanese Grants-In-Aid for Scientific Research (KAKENHI, 21H02591 to TN), the Koyanagi Foundation (to TN), the Ohsumi Frontier Science Foundation (to TN), and the Japanese Association of University Woman (JAUW) to AA.

Acknowledgments

We thank Drs. Ikuo K. Suzuki and Takuma Kumamoto for critically reading the manuscript and Ms. Kawami Misato and Ms. Mariko Yazaki for technical assistance.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fcell.2023.1247361/full#supplementary-material

References

Alappat, S., Zhang, Z. Y., and Chen, Y. P. (2003). Msx homeobox gene family and craniofacial development. Cell Res. 13, 429–442. doi:10.1038/sj.cr.7290185

PubMed Abstract | CrossRef Full Text | Google Scholar

Asgari, S., Luo, Y., Akbari, A., Belbin, G. M., Li, X., Harris, D. N., et al. (2020). A positively selected FBN1 missense variant reduces height in Peruvian individuals. Nature 582, 234–239. doi:10.1038/s41586-020-2302-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Cargill, M., Altshuler, D., Ireland, J., Sklar, P., Ardlie, K., Patil, N., et al. (1999). Characterization of single-nucleotide polymorphisms in coding regions of human genes. Nat. Genet. 22, 231–238. doi:10.1038/10290

PubMed Abstract | CrossRef Full Text | Google Scholar

Carroll, S. B. (2000). Endless forms: the evolution of gene regulation and morphological diversity. Cell 101, 577–580. doi:10.1016/s0092-8674(00)80868-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Castellano, S., Parra, G., Sánchez-Quinto, F. A., Racimo, F., Kuhlwilm, M., Kircher, M., et al. (2014). Patterns of coding variation in the complete exomes of three Neandertals. Proc. Natl. Acad. Sci. U. S. A. 111, 6666–6671. doi:10.1073/pnas.1405138111

PubMed Abstract | CrossRef Full Text | Google Scholar

Chan, Y. F., Marks, M. E., Jones, F. C., Villarreal, G., Shapiro, M. D., Brady, S. D., et al. (2010). Adaptive evolution of pelvic reduction in sticklebacks by recurrent deletion of a Pitx1 enhancer. Science 327, 302–305. doi:10.1126/science.1182213

PubMed Abstract | CrossRef Full Text | Google Scholar

Chasman, D., and Adams, R. M. (2001). Predicting the functional consequences of non-synonymous single nucleotide polymorphisms: structure-based assessment of amino acid variation. J. Mol. Biol. 307, 683–706. doi:10.1006/jmbi.2001.4510

PubMed Abstract | CrossRef Full Text | Google Scholar

Demurger, F., Ichkou, A., Mougou-Zerelli, S., Le Merrer, M., Goudefroye, G., Delezoide, A. L., et al. (2015). New insights into genotype-phenotype correlation for GLI3 mutations. Eur. J. Hum. Genet. 23, 92–102. doi:10.1038/ejhg.2014.62

PubMed Abstract | CrossRef Full Text | Google Scholar

Dong, X., Pi, Q., Yuemaierabola, A., Guo, W., and Tian, H. (2021). Silencing LINC00294 restores mitochondrial function and inhibits apoptosis of glioma cells under hypoxia via the miR-21-5p/CASKIN1/cAMP Axis. Oxid. Med. Cell Longev. 2021, 8240015. doi:10.1155/2021/8240015

PubMed Abstract | CrossRef Full Text | Google Scholar

Enard, W., Przeworski, M., Fisher, S. E., Lai, C. S. L., Wiebe, V., Kitano, T., et al. (2002). Molecular evolution of FOXP2, a gene involved in speech and language. Nature 418, 869–872. doi:10.1038/nature01025

PubMed Abstract | CrossRef Full Text | Google Scholar

Filvaroff, E. H., Guillet, S., Zlot, C., Bao, M., Ingle, G., Steinmetz, H., et al. (2002). Stanniocalcin 1 alters muscle and bone structure and function in transgenic mice. Endocrinology 143, 3681–3690. doi:10.1210/en.2001-211424

PubMed Abstract | CrossRef Full Text | Google Scholar

Garcia-Martinez, D., Bastir, M., Gómez-Olivencia, A., Maureille, B., Golovanova, L., Doronichev, V., et al. (2020). Early development of the Neanderthal ribcage reveals a different body shape at birth compared to modern humans. Sci. Adv. 6, eabb4377. doi:10.1126/sciadv.abb4377

PubMed Abstract | CrossRef Full Text | Google Scholar

Gokhman, D., Mishol, N., de Manuel, M., de Juan, D., Shuqrun, J., Meshorer, E., et al. (2020). Reconstructing denisovan anatomy using DNA methylation maps. Cell 180, 601. doi:10.1016/j.cell.2020.01.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Gokhman, D., Mishol, N., de Manuel, M., de Juan, D., Shuqrun, J., Meshorer, E., et al. (2019). Reconstructing denisovan anatomy using DNA methylation maps. Cell 179, 180–192. doi:10.1016/j.cell.2019.08.035

PubMed Abstract | CrossRef Full Text | Google Scholar

Green, R. E., Krause, J., Briggs, A. W., Maricic, T., Stenzel, U., Kircher, M., et al. (2010). A draft sequence of the Neandertal genome. Science 328, 710–722. doi:10.1126/science.1188021

PubMed Abstract | CrossRef Full Text | Google Scholar

Greer, C., Bhakta, H., Ghanem, L., Refai, F., Linn, E., and Avella, M. (2021). Deleterious variants in genes regulating mammalian reproduction in Neanderthals, Denisovans and extant humans. Hum. Reprod. 36, 734–755. doi:10.1093/humrep/deaa347

PubMed Abstract | CrossRef Full Text | Google Scholar

Gurumurthy, C. B., Sato, M., Nakamura, A., Inui, M., Kawano, N., Islam, M. A., et al. (2019). Creation of CRISPR-based germline-genome-engineered mice without ex vivo handling of zygotes by i-GONAD. Nat. Protoc. 14, 2452–2482. doi:10.1038/s41596-019-0187-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Hamburger, V., and Hamilton, H. L. (1951). A series of normal stages in the development of the chick embryo. J. Morphol. 88, 49–92. doi:10.1002/jmor.1050880104

PubMed Abstract | CrossRef Full Text | Google Scholar

Hoekstra, H. E., and Coyne, J. A. (2007). The locus of evolution: evo devo and the genetics of adaptation. Evolution 61, 995–1016. doi:10.1111/j.1558-5646.2007.00105.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Kang, S., Graham, J. M., Olney, A. H., and Biesecker, L. G. (1997). GLI3 frameshift mutations cause autosomal dominant Pallister-Hall syndrome. Nat. Genet. 15, 266–268. doi:10.1038/ng0397-266

PubMed Abstract | CrossRef Full Text | Google Scholar

Kohler, S., Gargano, M., Matentzoglu, N., Carmody, L. C., Lewis-Smith, D., Vasilevsky, N. A., et al. (2021). The human phenotype Ontology in 2021. Nucleic Acids Res. 49, D1207–D1217. doi:10.1093/nar/gkaa1043

PubMed Abstract | CrossRef Full Text | Google Scholar

Kvon, E. Z., Kamneva, O. K., Melo, U. S., Barozzi, I., Osterwalder, M., Mannion, B. J., et al. (2016). Progressive loss of function in a limb enhancer during snake evolution. Cell 167, 633–642. doi:10.1016/j.cell.2016.09.028

PubMed Abstract | CrossRef Full Text | Google Scholar

Labun, K., Montague, T. G., Krause, M., Torres Cleuren, Y. N., Tjeldnes, H., and Valen, E. (2019). CHOPCHOP v3: expanding the CRISPR web toolbox beyond genome editing. Nucleic Acids Res. 47, W171–W174. doi:10.1093/nar/gkz365

PubMed Abstract | CrossRef Full Text | Google Scholar

Levi, G., Narboux-Neme, N., and Cohen-Solal, M. (2022). DLX genes in the development and maintenance of the vertebrate skeleton: implications for human pathologies. Cells 11, 3277. doi:10.3390/cells11203277

PubMed Abstract | CrossRef Full Text | Google Scholar

Maricic, T., Helmbrecht, N., Riesenberg, S., Macak, D., Kanis, P., Lackner, M., et al. (2021). Comment on "Reintroduction of the archaic variant of NOVA1 in cortical organoids alters neurodevelopment. Science 374, eabi6060. doi:10.1126/science.abi6060

PubMed Abstract | CrossRef Full Text | Google Scholar

Matissek, S. J., and Elsawa, S. F. (2020). GLI3: a mediator of genetic diseases, development and cancer. Cell Commun. Signal 18, 54. doi:10.1186/s12964-020-00540-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Mc, L. A., and Michie, D. (1958). Factors affecting vertebral variation in mice. 4. Experimental proof of the uterine basis of a maternal effect. J. Embryol. Exp. Morphol. 6, 645–659.

PubMed Abstract | Google Scholar

Meyer, M., Kircher, M., Gansauge, M. T., Li, H., Racimo, F., Mallick, S., et al. (2012). A high-coverage genome sequence from an archaic Denisovan individual. Science 338, 222–226. doi:10.1126/science.1224344

PubMed Abstract | CrossRef Full Text | Google Scholar

Mora-Bermudez, F., Kanis, P., Macak, D., Peters, J., Naumann, R., Xing, L., et al. (2022). Longer metaphase and fewer chromosome segregation errors in modern human than Neanderthal brain development. Sci. Adv. 8, eabn7702. doi:10.1126/sciadv.abn7702

PubMed Abstract | CrossRef Full Text | Google Scholar

Naruse, I., Ueta, E., Sumino, Y., Ogawa, M., and Ishikiriyama, S. (2010). Birth defects caused by mutations in human GLI3 and mouse Gli3 genes. Congenit. Anom. (Kyoto) 50, 1–7. doi:10.1111/j.1741-4520.2009.00266.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Orlowski, J., and Grinstein, S. (1997). Na+/H+ exchangers of mammalian cells. J. Biol. Chem. 272, 22373–22376. doi:10.1074/jbc.272.36.22373

PubMed Abstract | CrossRef Full Text | Google Scholar

Pinson, A., Xing, L., Namba, T., Kalebic, N., Peters, J., Oegema, C. E., et al. (2022). Human TKTL1 implies greater neurogenesis in frontal neocortex of modern humans than Neanderthals. Science 377, eabl6422. doi:10.1126/science.abl6422

PubMed Abstract | CrossRef Full Text | Google Scholar

Protas, M. E., Hersey, C., Kochanek, D., Zhou, Y., Wilkens, H., Jeffery, W. R., et al. (2006). Genetic analysis of cavefish reveals molecular convergence in the evolution of albinism. Nat. Genet. 38, 107–111. doi:10.1038/ng1700

PubMed Abstract | CrossRef Full Text | Google Scholar

Prufer, K., de Filippo, C., Grote, S., Mafessoni, F., Korlević, P., Hajdinjak, M., et al. (2017). A high-coverage neandertal genome from vindija cave in Croatia. Science 358, 655–658. doi:10.1126/science.aao1887

PubMed Abstract | CrossRef Full Text | Google Scholar

Prufer, K., Racimo, F., Patterson, N., Jay, F., Sankararaman, S., Sawyer, S., et al. (2014). The complete genome sequence of a Neanderthal from the Altai Mountains. Nature 505, 43–49. doi:10.1038/nature12886

PubMed Abstract | CrossRef Full Text | Google Scholar

Qiu, J., Zhou, S., Cheng, W., and Luo, C. (2020). LINC00294 induced by GRP78 promotes cervical cancer development by promoting cell cycle transition. Oncol. Lett. 20, 262. doi:10.3892/ol.2020.12125

PubMed Abstract | CrossRef Full Text | Google Scholar

Raudvere, U., Kolberg, L., Kuzmin, I., Arak, T., Adler, P., Peterson, H., et al. (2019). g:Profiler: a web server for functional enrichment analysis and conversions of gene lists (2019 update). Nucleic Acids Res. 47, W191–W198. doi:10.1093/nar/gkz369

PubMed Abstract | CrossRef Full Text | Google Scholar

Reich, D., Green, R. E., Kircher, M., Krause, J., Patterson, N., Durand, E. Y., et al. (2010). Genetic history of an archaic hominin group from Denisova Cave in Siberia. Nature 468, 1053–1060. doi:10.1038/nature09710

PubMed Abstract | CrossRef Full Text | Google Scholar

Reilly, P. F., Tjahjadi, A., Miller, S. L., Akey, J. M., and Tucci, S. (2022). The contribution of Neanderthal introgression to modern human traits. Curr. Biol. 32, R970–R983. doi:10.1016/j.cub.2022.08.027

PubMed Abstract | CrossRef Full Text | Google Scholar

Renard, E., Chéry, C., Oussalah, A., Josse, T., Perrin, P., Tramoy, D., et al. (2019). Exome sequencing of cases with neural tube defects identifies candidate genes involved in one-carbon/vitamin B12 metabolisms and Sonic Hedgehog pathway. Hum. Genet. 138, 703–713. doi:10.1007/s00439-019-02015-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Rengasamy, P., and Padmanabhan, R. R. (2004). Experimental studies on cervical and lumbar ribs in mouse embryos. Congenit. Anom. (Kyoto) 44, 156–171. doi:10.1111/j.1741-4520.2004.00029.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Rigueur, D., and Lyons, K. M. (2014). Whole-mount skeletal staining. Methods Mol. Biol. 1130, 113–121. doi:10.1007/978-1-62703-989-5_9

PubMed Abstract | CrossRef Full Text | Google Scholar

Ronshaugen, M., McGinnis, N., and McGinnis, W. (2002). Hox protein mutation and macroevolution of the insect body plan. Nature 415, 914–917. doi:10.1038/nature716

PubMed Abstract | CrossRef Full Text | Google Scholar

Sankararaman, S., Mallick, S., Patterson, N., and Reich, D. (2016). The combined landscape of denisovan and neanderthal ancestry in present-day humans. Curr. Biol. 26, 1241–1247. doi:10.1016/j.cub.2016.03.037

PubMed Abstract | CrossRef Full Text | Google Scholar

Sasaki, H., Hui, C., Nakafuku, M., and Kondoh, H. (1997). A binding site for Gli proteins is essential for HNF-3beta floor plate enhancer activity in transgenics and can respond to Shh in vitro. Development 124, 1313–1322. doi:10.1242/dev.124.7.1313

PubMed Abstract | CrossRef Full Text | Google Scholar

Sasaki, H., Nishizaki, Y., Hui, C., Nakafuku, M., and Kondoh, H. (1999). Regulation of Gli2 and Gli3 activities by an amino-terminal repression domain: implication of Gli2 and Gli3 as primary mediators of Shh signaling. Development 126, 3915–3924. doi:10.1242/dev.126.17.3915

PubMed Abstract | CrossRef Full Text | Google Scholar

Shumate, A., Wong, B., Pertea, G., and Pertea, M. (2022). Improved transcriptome assembly using a hybrid of long and short reads with StringTie. PLoS Comput. Biol. 18, e1009730. doi:10.1371/journal.pcbi.1009730

PubMed Abstract | CrossRef Full Text | Google Scholar

Skov, L., Coll Macià, M., Sveinbjörnsson, G., Mafessoni, F., Lucotte, E. A., Einarsdóttir, M. S., et al. (2020). The nature of Neanderthal introgression revealed by 27,566 Icelandic genomes. Nature 582, 78–83. doi:10.1038/s41586-020-2225-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Somiya, M., and Kuroda, S. (2021). Reporter gene assay for membrane fusion of extracellular vesicles. J. Extracell. Vesicles 10, e12171. doi:10.1002/jev2.12171

PubMed Abstract | CrossRef Full Text | Google Scholar

Stamataki, D., Ulloa, F., Tsoni, S. V., Mynett, A., and Briscoe, J. (2005). A gradient of Gli activity mediates graded Sonic Hedgehog signaling in the neural tube. Genes Dev. 19, 626–641. doi:10.1101/gad.325905

PubMed Abstract | CrossRef Full Text | Google Scholar

Tessadori, F., Duran, K., Knapp, K., Fellner, M., and Smithson, S. (2022). Recurrent de novo missense variants across multiple histone H4 genes underlie a neurodevelopmental syndrome. Am. J. Hum. Genet. 109, 750–758. doi:10.1016/j.ajhg.2022.02.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Trujillo, C. A., Rice, E. S., Schaefer, N. K., Chaim, I. A., Wheeler, E. C., Madrigal, A. A., et al. (2021). Reintroduction of the archaic variant of NOVA1 in cortical organoids alters neurodevelopment. Science 371, eaax2537. doi:10.1126/science.aax2537

PubMed Abstract | CrossRef Full Text | Google Scholar

Veistinen, L. K., Mustonen, T., Hasan, M. R., Takatalo, M., Kobayashi, Y., Kesper, D. A., et al. (2017). Regulation of calvarial osteogenesis by concomitant de-repression of GLI3 and activation of IHH targets. Front. Physiol. 8, 1036. doi:10.3389/fphys.2017.01036

PubMed Abstract | CrossRef Full Text | Google Scholar

Veistinen, L., Takatalo, M., Tanimoto, Y., Kesper, D. A., Vortkamp, A., and Rice, D. P. C. (2012). Loss-of-Function of Gli3 in mice causes abnormal frontal bone morphology and premature synostosis of the interfrontal suture. Front. Physiol. 3, 121. doi:10.3389/fphys.2012.00121

PubMed Abstract | CrossRef Full Text | Google Scholar

Vortkamp, A., Gessler, M., and Grzeschik, K. H. (1991). GLI3 zinc-finger gene interrupted by translocations in Greig syndrome families. Nature 352, 539–540. doi:10.1038/352539a0

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, B., Fallon, J. F., and Beachy, P. A. (2000). Hedgehog-regulated processing of Gli3 produces an anterior/posterior repressor gradient in the developing vertebrate limb. Cell 100, 423–434. doi:10.1016/s0092-8674(00)80678-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, C., Ruther, U., and Wang, B. (2007). The Shh-independent activator function of the full-length Gli3 protein and its role in vertebrate limb digit patterning. Dev. Biol. 305, 460–469. doi:10.1016/j.ydbio.2007.02.029

PubMed Abstract | CrossRef Full Text | Google Scholar

Watanabe, K., Takebayashi, H., Bepari, A. K., Esumi, S., Yanagawa, Y., and Tamamaki, N. (2011). Dpy19l1, a multi-transmembrane protein, regulates the radial migration of glutamatergic neurons in the developing cerebral cortex. Development 138, 4979–4990. doi:10.1242/dev.068155

PubMed Abstract | CrossRef Full Text | Google Scholar

Wen, X., Lai, C. K., Evangelista, M., Hongo, J. A., de Sauvage, F. J., and Scales, S. J. (2010). Kinetics of hedgehog-dependent full-length Gli3 accumulation in primary cilia and subsequent degradation. Mol. Cell Biol. 30, 1910–1922. doi:10.1128/MCB.01089-09

PubMed Abstract | CrossRef Full Text | Google Scholar

Werner, T., Koshikawa, S., Williams, T. M., and Carroll, S. B. (2010). Generation of a novel wing colour pattern by the Wingless morphogen. Nature 464, 1143–1148. doi:10.1038/nature08896

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeberg, H., Dannemann, M., Sahlholm, K., Tsuo, K., Maricic, T., Wiebe, V., et al. (2020). A neanderthal sodium channel increases pain sensitivity in present-day humans. Curr. Biol. 30, 3465–3469. doi:10.1016/j.cub.2020.06.045

PubMed Abstract | CrossRef Full Text | Google Scholar

Zeng, H., Jia, J., and Liu, A. (2010). Coordinated translocation of mammalian Gli proteins and suppressor of fused to the primary cilium. PLoS One 5, e15900. doi:10.1371/journal.pone.0015900

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Y., Park, C., Bennett, C., Thornton, M., and Kim, D. (2021). Rapid and accurate alignment of nucleotide conversion sequencing reads with HISAT-3N. Genome Res. 31, 1290–1295. doi:10.1101/gr.275193.120

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, X., Lv, L., Zhang, Z., Wei, S., and Zheng, T. (2020). LINC00294 negatively modulates cell proliferation in glioma through a neurofilament medium-mediated pathway via interacting with miR-1278. J. Gene Med. 22, e3235. doi:10.1002/jgm.3235

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, H., Kim, S., Ishii, S., and Boyer, T. G. (2006). Mediator modulates Gli3-dependent Sonic hedgehog signaling. Mol. Cell Biol. 26, 8667–8682. doi:10.1128/MCB.00443-06

PubMed Abstract | CrossRef Full Text | Google Scholar