Analysis of Genes Associated with Both Neural Tube Defects and Neuroectodermal Tumors

Background

During embryonic development, a fertilized egg develops into an embryo by undergoing cell division and subsequent cell differentiation processes. Similarly, tumors develop by the abnormal proliferation and differentiation of cells due to the presence of mutations in key genes that regulate the cell cycle. The concept of the embryogenic origin of tumors, which was first suggested by the French biologists Lobstein and Recamier in 1892, states that tumor cells are very similar to embryonic cells and that tumor cells are derived from the continuous proliferation of embryonic cells [1]. In 1983, Pierce et al proposed that a tumor was a developmental biology problem and further highlighted the close relationship between the occurrence of tumors and embryo development [2].

Neural tube defects (NTDs) are a common congenital birth defect caused by the inadequate closure of the neural tube during embryo development. In most vertebrates, hollow neural tubes are formed as neuroectoderm cells proliferate, invaginate, and eventually migrate from the surface of the ectoderm [3–6]. Among all types of human tumor tissues, 6 classes derive from the neuroectoderm during early embryo development: primitive neuroectodermal tumor, glioma, glioblastoma (GBM), neuroblastoma, neurilemmoma, and medulloblastoma [7,8]. GBM is a brain tumor with high incidence and mortality rates, for which only a limited number of treatments are currently available.

The neural tube is the main primordium of the central nervous system. It gives rise to the glial cells, which retain the ability to proliferate throughout life. Most adult neurological tumors are of glial origin. These tumors are termed gliomas [9], and GBM is the most malignant form of glioma [10]. Therefore, from a developmental perspective, the embryonic origin of GBM tissue is the neural tube. It is perhaps not surprising that there are case reports of patients with both GBM and NTDs [11] and that individuals with congenital nervous system defects are at higher risk of developing neurological tumors [12]. One out of 10 babies with NTDs die before their first year [13]. GBM is also a disease with a high mortality rate, and the 5-year survival rate is less than 10% [14,15]. Although it is difficult to obtain samples of both diseases in the same patient, it does not prevent scientists from studying the shared genes and potential mechanisms of the 2 diseases. Previous studies have reported several genes associated with both NTDs and GBM. To the best of our knowledge, the p53 pathway contributes to NTDs and to GBM pathogenesis [16]; aberrant expression of the platelet-derived growth factor receptor-α (PDGFRA) gene has been associated with NTDs and GBM [17]; and insulin-like growth factor binding protein 2 (IGFBP2) mRNA levels remain high during the neural tube closure and are frequently overexpressed in GBMs [18]. Nonetheless, the underlying molecular mechanisms are still poorly understood, highlighting the need to better understand the mechanisms regulating the expression of key genes involved in cell proliferation during embryo development and tumorigenesis. Many studies have shown that microRNAs (miRNAs) such as miR-30, miR-9, and miR-375 play an important role in neural tube closure during the early embryonic development and in tumorigenesis [19–25].

Therefore, finding core genes and potential pathways associated with both NTDs and neuroectodermal tumors will contribute to understanding their pathogenesis as well as developing targeted treatments. Our study revealed the key genes associated with NTDs and neuroectodermal tumors, and their potential mechanism.

Material and Methods

ANIMALS:

ICR mice (9–11 weeks, 19–25 g) were provided by the Animal Center of Shanxi Medical University and raised in specific pathogen-free cages at the animal center with a 12-h light/dark cycle, keeping the temperature and humidity within the range of 20–26ºC and 40–70%, respectively. The NTDs mouse model was established by gavage of 28 mg/kg retinoic acid (RA) to pregnant mice at embryonic day 7.5 (E7.5). The control pregnant mice were treated with the same dose of sesame oil. The brain vesicles of mouse embryos were collected at E9.5 and E10.5. Total RNA was extracted using the TRIzol method (Ambion, USA), and cDNA was synthesized using a Revert Aid First Strand cDNA Synthesis Kit (Takara, Japan).

COLLECTION OF GBM SAMPLES:

All GBM samples were collected from the Neurosurgery Department of Shanxi Cancer Hospital. Normal adjacent tissue was obtained from the extended excision domain of GBM neurosurgical resections. Ethics approval was obtained from the Ethics Review Board of Shanxi Cancer Hospital.

GENE SET ANALYSIS:

The DisGeNET database was used to identify genes associated with NTDs and 6 kinds of neuroectodermal tumors. DisGeNET is a platform integrating standardized data on disease-associated genes from multiple sources. The current release covers more than 24 000 diseases and traits, 17 000 genes, and 117 000 genomic variants [26,27]. It provides a ‘score_gda’ that reflects how well established a particular association is based on current knowledge by giving the highest value to associations reported by several databases, particularly those reported by expert-curated resources, and with a large number of supporting publications. We selected different kinds of disease-associated genes with a score_gda higher than the mean score of all genes.

GENE EXPRESSION LEVEL ANALYSIS:

To confirm the gene expression levels of the identified disease-associated genes, we used the mRNA-sequencing data from mouse and human NTDs samples. Since the NTDs-susceptibility genes screened using animal models do not necessarily reflect the human pathological mechanisms, we considered only those genes that had similar expression trends both in mouse models and human NTDs samples as candidate genes.

mRNA sequencing data from the RA-induced NTDs mouse model were obtained by our research team, as reported in our previous study [28]. Transcriptome data were collected at E8.5, E9.5, and E10.5. Briefly, the samples were sequenced using the Illumina HiSeqTM 2000 sequencing platform. mRNA sequencing data of human NTDs samples were obtained by our research team in a previous study [29]. Briefly, in that previous study, embryos with a gestational age of approximately 17 weeks were diagnosed with spina bifida and age-matched normal embryos were obtained from abortions. Affymetrix HG-U133A 2.0 GeneChip arrays were used to examine the gene expression profile of these samples.

HUB GENES ANALYSIS:

The STRING database was used to construct a protein–protein interaction (PPI) network. This database mainly comprises information about interacting proteins collected from experimental data. The current version contains 3 123 056 667 interaction pairs, which originate from 24 584 628 proteins and 5090 organisms [30]. Cytoscape was used to visualize the network [31]. Hub genes were selected using the Maximal Clique Centrality algorithm [32].

PREDICTION AND CONFIRMATION OF THE POTENTIAL MIRNA-MRNA PAIRS:

StarBase was used to predict the miRNAs targeting the hub genes identified. StarBase is a database facilitating the comprehensive exploration of miRNA-target interaction maps from CLIP-Seq and Degradome-Seq data. The current version includes high-throughput sequencing data generated from 21 CLIP-Seq and 10 Degradome-Seq studies [33]. To confirm the potential miRNA-mRNA pairs, we used the miRNA-sequencing data from the mouse model of NTDs that was previously published by our research team [34]. Here, the NTDs mouse model was constructed identically according to our previous study. Using these data, we selected the miRNAs with opposite expression trends to the hub genes as candidate miRNAs.

QUANTITATIVE REVERSE TRANSCRIPTION PCR (RT-QPCR) ANALYSIS:

The mRNA expression level of hub genes in the RA-induced NTDs mice and human GBM samples was assessed using RT-qPCR. qPCR was performed using Maxima SYBR Green/ROX qPCR Master Mix (Takara, Japan). The data were analyzed with 2−ΔΔCt method. The mRNA expression level of the target genes was normalized to that of Gapdh. Primer sequences are shown in the Supplementary Tables 1 and 2.

DUAL-LUCIFERASE ANALYSIS:

A dual-luciferase report assay was used for the validation of the miRNA-gene pairs. The luciferase reporter plasmid was built using the pmirGLO vector, into which the wild-type and mutant candidate genes were cloned. For this purpose, an appropriate number of 293T cells were seeded into a 12-well plate and cultured at 37°C in an incubator overnight. Cells were co-transfected with the luciferase reporter plasmids and miRNA mimics/mimics control using Lipofectamine 2000 (Invitrogen), and the efficiency of transfection was assessed. After 24 h, Firefly and Renilla luciferase activity were detected with a Dual-Luciferase Reporter Assay kit (PR-E1910, Promega, Wisconsin, USA). Relative luciferase activity was defined as the ratio of Renilla luciferase activity to the Firefly luciferase activity with that of the control set as 1.0. Briefly, cells were washed with phosphate-buffered saline, and 400 μL of 1×Passive Lysis Buffer (Promega) was added to the cultured well and shaken gently for 15 min. The lysate was then transferred to a test tube. A total of 20 μL Passive Lysis Buffer lysis buffer was transferred into a white 96-well plate, then 100 μL of LARII (Promega) was added to detect Firefly luciferase in Multifunctional Enzyme Marker. A total of 100 μL of Stop&Glo (Promega) was used to detect Renilla luciferase activity.

STATISTICAL ANALYSIS:

Statistical analysis was performed using GraphPad Prism software, version 8.0. Data were expressed as the mean±S.D. The t test analysis was performed to compare the means between the 2 groups.

Results

INTERSECTING GENES BETWEEN NTDS AND NEUROECTODERMAL TUMORS:

A total of 304 NTDs-associated genes and 9572 neuroectodermal tumors-associated genes were downloaded from the DisGeNET database. Table 1 shows the number of genes associated with the 6 kinds of neuroectodermal tumors. Subsequently, we identified the genes that were associated with both NTDs and neuroectodermal tumors. A total of 104 overlapping genes were considered as candidate genes for further analyses (Figure 1).

BIOINFORMATICS ANALYSIS OF THE GENES:

Among the 104 overlapping genes, differentially expressed genes (DEGs) were selected using the threshold criteria of P value ≤0.05 and |Log2Ratio| ≥1 based on NTDs mouse transcriptome data. A total of 3 genes were differentially expressed at E8.5, 14 at E9.5, and 24 at E10.5, as shown in Figure 2.

A PPI interaction network was constructed and visualized for the 27 DEGs using the STRING database and Cytoscape software. The PPI revealed that 20 of these interacted very closely. The resulting 10 hub genes identified using the maximal clique centrality algorithm were: Pou5f1, Bmp4, Shh, Cdh1, Ncam1, Snai2, Igf2, Ascl1, Anxa5, and Ccl2 (Figure 3).

Four of the 10 hub genes had the same expression trend in human NTDs data as in mouse NTDs data. These genes were Ncam1 (up-regulated), Cdh1 (up-regulated), Shh (down-regulated), and Ascl1 (down-regulated). In summary, these 4 genes were closely associated with both NTDs and the 6 types of neuroectodermal tumors, and were also key genes in the PPI network, suggesting they are important candidate genes.

Hence, we further explored the potential mechanism whereby these genes are involved in the pathological mechanism of both NTDs and neuroectodermal tumors. As miRNAs are a type of small non-coding RNA that plays an important role in regulating gene expression [35], using StarBase V2.0, we predicted miRNAs regulate expression levels of the candidate genes. A total of 180 miRNAs with 245 targeted binding sites were predicted for Ncam1, 130 miRNAs with 170 targeted binding sites were predicted for Cdh1, 14 miRNAs with 14 sites were predicted for Shh, and 24 miRNAs with 34 sites were predicted for Ascl1.

We then used the RA-induced NTDs mouse model’s mRNA-sequencing and miRNA-sequencing data to filter the predicted results. For this purpose, only miRNA-mRNA pairs with opposite expression trends both at E9.5 and E10.5 were retained. A total of 38 miRNA binding sites for Ncam1, 1 miRNA binding site for Ascl1, and 32 miRNA binding sites for Cdh1 were obtained. miRNAs with the highest differential expression fold value were selected. Finally, we considered mmu-miR-30a-5p-Ncam1, mmu-miR-375-3p-Ascl1, and mmu-miR-9-5p-Cdh1 as candidate pairs for experimental verification, as shown in Figure 4.

VALIDATION OF THE EXPRESSION LEVEL OF THE TARGET GENES:

To validate the expression levels of the candidate genes, we assessed the mRNA levels of Ncam1, Cdh1, Shh, and Ascl1 in the brain vesicle tissue of E9.5 and E10.5 mouse embryos using RT-qPCR. Ncam1 and Cdh1 were upregulated in the NTDs mice, whereas Shh and Ascl1 were downregulated in NTDs mice. Additionally, Ncam1, Shh, and Ascl1 were significantly differentially expressed at E9.5 and E10.5; whereas Cdh1 was significantly differentially expressed at E10.5, as shown in Figures 5 and 6. Moreover, we validated expression levels of these 4 candidate genes in GBM samples using RT-qPCR. Shh and Ncam1 were significantly upregulated and Ascl1 was significantly downregulated in tumor samples, whereas Cdh1 was not differentially expressed in the tumor samples, as shown in Figure 7. In summary, among the 4 important susceptibility genes, 3 genes (Ncam1, Shh, and Ascl1) were differentially expressed both in RA-induced NTDs mice and human GBM samples.

VALIDATION OF THE INTERACTION BETWEEN THE MIRNAS AND THE TARGET GENES:

To validate the interaction between the selected miRNA-mRNA pairs, we performed a luciferase assay. For this purpose, we co-transfected 293T cells with the luciferase reporter plasmids and mimic control/miRNAs mimic inhibitor using Lipofectamine 2000. To ensure the transfection effect of the formal experiment, we conducted a preliminary experiment using the empty vector plasmid with a red fluorescent signal and acquired images using an ordinary fluorescence microscope to visually show the transfection effect. As shown in Figure 8, the transfection efficiency was high. The results of the dual-luciferase assay are shown in Figure 9. Transfection with the wild-type mmu-miR-30a-5p had an approximately 50% inhibitory effect on Ncam1, whereas the inhibitory effect was not seen with the mutated miRNA. Transfection with wild-type mmu-miR-375-3p inhibited approximately 30% of the expression level of Ascl1, and the inhibitory effect was not seen with the mutated miRNA. However, wild-type mmu-miR-9-5p had no significant effect on the expression of Cdh1, indicating that it did not have a significant specific target binding site on the Cdh1 transcript and that the results predicted by the software lacked sufficient accuracy. The predicted binding sites of the miRNAs and target genes, and the sequence alignments, are shown in Figure 10.

Discussion

With the continuous development of molecular biology tools and assays, a growing number of studies have revealed the similarities between early embryonic cells and tumor cells and the key genes, proteins, and metabolic pathways involved in both processes [36–38]. A landmark study published in the journal Science in 2019 revealed that the most common type of childhood kidney cancer is usually preceded by clonal expansions during embryonic development [39]. A follow-up study published in March 2020 showed extensive mutagenesis in the placenta and that the mutation sites were the same as in many childhood cancers [40]. In 2020, a research team from Singapore Genomics Institute and Singapore National Cancer Centre demonstrated that fetal liver and liver cancer tissues shared an immunosuppressive onco-fetal ecosystem, suggesting that atavistic remodeling might occur in the tumor microenvironment [41]. Based on a large-population cohort study, cancer risk in children with birth defects was precisely estimated in cancers, such as pyloric stenosis and medulloblastoma, several different cardiac phenotypes, and neuroblastoma [42].

Another large Nordic population-based case-control study indicated that cancer risk increased in individuals with birth defects and persisted into adulthood, both for non-chromosomal and chromosomal anomalies [43]. In our studies, through a series of bioinformatics analyses, we identified 4 important susceptibility genes associated with both NTDs and neuroectodermal tumors: Ncam1, Cdh1, Shh, and Ascl1. Three of these genes were differentially expressed both in RA-induced NTDs mice and human glioblastoma samples: Ncam1, Shh, and Ascl1. Ncam1 has been previously reported to be a candidate gene involved in NTDs [44], and the novel Ncam1 variant was reported to play an important role in cell signaling associated with tumor development in GBM [45]. Shh encodes a neural tube patterning-related protein [46], which is also a molecular marker for medulloblastoma [47]. Finally, the expression level of Ascl1 was also found to be altered in the neural tube of embryos from diabetic mice [48], and it played an important role in controlling the proliferation of GBM cells [49,50].

In conclusion, previous evidence supports a role for these 3 genes in both development and tumorigenesis. Future studies should evaluate the molecular features of patients with co-occurrence of NTDs and neuroectodermal tumors to further elucidate the mechanisms that lead to these complex outcomes. Upon further validation, our findings may help in the development of cancer surveillance protocols for early tumor detection in children with neural tube defects. However, if the models or samples we used have NTDs initially and eventually develop neuroectodermal tumors with the same genetic backgrounds, the DEGs we found will better prove the connection between NTDs and neuroectodermal tumors. Moreover, the RA-induced NTDs mouse model is not a representative NTDs mouse model; therefore, our results obtained from this model alone cannot fully explain the relationship between NTDs and neuroectodermal tumors. Additional models will be required to validate these findings. Furthermore, the susceptibility genes we detected have not been verified in human NTDs samples, and the pathophysiologic significance and clinical utility of our results remain to be determined.

miRNAs are a class of non-coding single-stranded RNA molecules about 22 nucleotides in length. miRNAs are known to participate in the regulation of gene expression at the post-transcriptional level by binding to the target transcript [51]. Through this mechanism, miRNAs drive embryonic development in a variety of species, such as zebrafish, Xenopus toads, and mice [52–54]. A comprehensive study by Guo et al supports a close relationship of miRNA expression between early human embryonic development and tumorigenesis [55]. Many studies have shown that miRNAs are involved in the occurrence of NTDs. miR-30 family members play an important role in regulating the WNT signaling pathway, TGF-β signaling pathway, and focal adhesion during neural differentiation [56]. miR-375 is upregulated in the NTDs mouse model and might be a potential pathway leading to spina bifida [24]. Similarly, miRNAs are also known to regulate gene expression during tumorigenesis, thereby regulating the process. miR-30 family members are downregulated in medulloblastoma, leading to enhanced cancer cell growth and other malignant behavior, and thus is a potential new therapeutic target [20]. In neuroblastoma, miR-375 regulates nerve cells differentiation, thereby contributing to tumorigenesis [57].

Our study verified 2 targeted binding pairs: mmu-miR-30a-5p-Ncam1 and mmu-miR-375-3p-Ascl1. These 2 pairs can regulate the development of neuroectodermal cells, and may also be related to neural tube development, making them potential biomarkers of NTDs and targets for GBM treatment. Further follow-up studies are needed to address the specific underlying mechanisms.

Conclusions

In summary, susceptibility genes related to both NTDs and 6 types of neuroectodermal tumors were obtained. Using bioinformatic analysis, we identified 4 potential candidate genes. Additionally, we demonstrated that 3 of these genes were differentially expressed both in NTDs mice and human GBM samples, and 2 were significantly regulated by miRNAs. This suggests that the dysregulation of their transcriptional control might lead to diseases such as NTDs or neuroectoderm tumors. Altogether, these results provide a theoretical framework supporting the relationship between early embryo development and tumorigenesis and provide the foundation for future studies aiming to identify novel therapies for NTDs or GBM.