SCExecute: custom cell barcode-stratified analyses of scRNA-seq data.

MOTIVATION: In single-cell RNA-sequencing (scRNA-seq) data, stratification of sequencing reads by cellular barcode is necessary to study cell-specific features. However, apart from gene expression, the analyses of cell-specific features are not sufficiently supported by available tools designed for high-throughput sequencing data. RESULTS: We introduce SCExecute, which executes a user-provided command on barcode-stratified, extracted on-the-fly, single-cell binary alignment map (scBAM) files. SCExecute extracts the alignments with each cell barcode from aligned, pooled single-cell sequencing data. Simple commands, monolithic programs, multi-command shell scripts or complex shell-based pipelines are then executed on each scBAM file. scBAM files can be restricted to specific barcodes and/or genomic regions of interest. We demonstrate SCExecute with two popular variant callers-GATK and Strelka2-executed in shell-scripts together with commands for BAM file manipulation and variant filtering, to detect single-cell-specific expressed single nucleotide variants from droplet scRNA-seq data (10X Genomics Chromium System).In conclusion, SCExecute facilitates custom cell-level analyses on barcoded scRNA-seq data using currently available tools and provides an effective solution for studying low (cellular) frequency transcriptome features. AVAILABILITY AND IMPLEMENTATION: SCExecute is implemented in Python3 using the Pysam package and distributed for Linux, MacOS and Python environments from https://horvathlab.github.io/NGS/SCExecute. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Transcriptional dysregulation in developing trigeminal sensory neurons in the LgDel mouse model of DiGeorge 22q11.2 deletion syndrome.

LgDel mice, which model the heterozygous deletion of genes at human chromosome 22q11.2 associated with DiGeorge/22q11.2 deletion syndrome (22q11DS), have cranial nerve and craniofacial dysfunction as well as disrupted suckling, feeding and swallowing, similar to key 22q11DS phenotypes. Divergent trigeminal nerve (CN V) differentiation and altered trigeminal ganglion (CNgV) cellular composition prefigure these disruptions in LgDel embryos. We therefore asked whether a distinct transcriptional state in a specific population of early differentiating LgDel cranial sensory neurons, those in CNgV, a major source of innervation for appropriate oropharyngeal function, underlies this departure from typical development. LgDel versus wild-type (WT) CNgV transcriptomes differ significantly at E10.5 just after the ganglion has coalesced. Some changes parallel altered proportions of cranial placode versus cranial neural crest-derived CNgV cells. Others are consistent with a shift in anterior-posterior patterning associated with divergent LgDel cranial nerve differentiation. The most robust quantitative distinction, however, is statistically verifiable increased variability of expression levels for most of the over 17 000 genes expressed in common in LgDel versus WT CNgV. Thus, quantitative expression changes of functionally relevant genes and increased stochastic variation across the entire CNgV transcriptome at the onset of CN V differentiation prefigure subsequent disruption of cranial nerve differentiation and oropharyngeal function in LgDel mice.

SCReadCounts: estimation of cell-level SNVs expression from scRNA-seq data.

BACKGROUND: Recent studies have demonstrated the utility of scRNA-seq SNVs to distinguish tumor from normal cells, characterize intra-tumoral heterogeneity, and define mutation-associated expression signatures. In addition to cancer studies, SNVs from single cells have been useful in studies of transcriptional burst kinetics, allelic expression, chromosome X inactivation, ploidy estimations, and haplotype inference. RESULTS: To aid these types of studies, we have developed a tool, SCReadCounts, for cell-level tabulation of the sequencing read counts bearing SNV reference and variant alleles from barcoded scRNA-seq alignments. Provided genomic loci and expected alleles, SCReadCounts generates cell-SNV matrices with the absolute variant- and reference-harboring read counts, as well as cell-SNV matrices of expressed Variant Allele Fraction (VAFRNA) suitable for a variety of downstream applications. We demonstrate three different SCReadCounts applications on 59,884 cells from seven neuroblastoma samples: (1) estimation of cell-level expression of known somatic mutations and RNA-editing sites, (2) estimation of cell- level allele expression of biallelic SNVs, and (3) a discovery mode assessment of the reference and each of the three alternative nucleotides at genomic positions of interest that does not require prior SNV information. For the later, we applied SCReadCounts on the coding regions of KRAS, where it identified known and novel somatic mutations in a low-to-moderate proportion of cells. The SCReadCounts read counts module is benchmarked against the analogous modules of GATK and Samtools. SCReadCounts is freely available ( https://github.com/HorvathLab/NGS ) as 64-bit self-contained binary distributions for Linux and MacOS, in addition to Python source. CONCLUSIONS: SCReadCounts supplies a fast and efficient solution for estimation of cell-level SNV expression from scRNA-seq data. SCReadCounts enables distinguishing cells with monoallelic reference expression from those with no gene expression and is applicable to assess SNVs present in only a small proportion of the cells, such as somatic mutations in cancer.

scReQTL: an approach to correlate SNVs to gene expression from individual scRNA-seq datasets.

BACKGROUND: Recently, pioneering expression quantitative trait loci (eQTL) studies on single cell RNA sequencing (scRNA-seq) data have revealed new and cell-specific regulatory single nucleotide variants (SNVs). Here, we present an alternative QTL-related approach applicable to transcribed SNV loci from scRNA-seq data: scReQTL. ScReQTL uses Variant Allele Fraction (VAFRNA) at expressed biallelic loci, and corelates it to gene expression from the corresponding cell. RESULTS: Our approach employs the advantage that, when estimated from multiple cells, VAFRNA can be used to assess effects of SNVs in a single sample or individual. In this setting scReQTL operates in the context of identical genotypes, where it is likely to capture RNA-mediated genetic interactions with cell-specific and transient effects. Applying scReQTL on scRNA-seq data generated on the 10 × Genomics Chromium platform using 26,640 mesenchymal cells derived from adipose tissue obtained from three healthy female donors, we identified 1272 unique scReQTLs. ScReQTLs common between individuals or cell types were consistent in terms of the directionality of the relationship and the effect size. Comparative assessment with eQTLs from bulk sequencing data showed that scReQTL analysis identifies a distinct set of SNV-gene correlations, that are substantially enriched in known gene-gene interactions and significant genome-wide association studies (GWAS) loci. CONCLUSION: ScReQTL is relevant to the rapidly growing source of scRNA-seq data and can be applied to outline SNVs potentially contributing to cell type-specific and/or dynamic genetic interactions from an individual scRNA-seq dataset. AVAILABILITY: https://github.com/HorvathLab/NGS/tree/master/scReQTL.

ReQTL: identifying correlations between expressed SNVs and gene expression using RNA-sequencing data.

MOTIVATION: By testing for associations between DNA genotypes and gene expression levels, expression quantitative trait locus (eQTL) analyses have been instrumental in understanding how thousands of single nucleotide variants (SNVs) may affect gene expression. As compared to DNA genotypes, RNA genetic variation represents a phenotypic trait that reflects the actual allele content of the studied system. RNA genetic variation at expressed SNV loci can be estimated using the proportion of alleles bearing the variant nucleotide (variant allele fraction, VAFRNA). VAFRNA is a continuous measure which allows for precise allele quantitation in loci where the RNA alleles do not scale with the genotype count. We describe a method to correlate VAFRNA with gene expression and assess its ability to identify genetically regulated expression solely from RNA-sequencing (RNA-seq) datasets. RESULTS: We introduce ReQTL, an eQTL modification which substitutes the DNA allele count for the variant allele fraction at expressed SNV loci in the transcriptome (VAFRNA). We exemplify the method on sets of RNA-seq data from human tissues obtained though the Genotype-Tissue Expression (GTEx) project and demonstrate that ReQTL analyses are computationally feasible and can identify a subset of expressed eQTL loci. AVAILABILITY AND IMPLEMENTATION: A toolkit to perform ReQTL analyses is available at https://github.com/HorvathLab/ReQTL. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Exosomes containing HIV protein Nef reorganize lipid rafts potentiating inflammatory response in bystander cells.

HIV infection has a profound effect on "bystander" cells causing metabolic co-morbidities. This may be mediated by exosomes secreted by HIV-infected cells and containing viral factors. Here we show that exosomes containing HIV-1 protein Nef (exNef) are rapidly taken up by macrophages releasing Nef into the cell interior. This caused down-regulation of ABCA1, reduction of cholesterol efflux and sharp elevation of the abundance of lipid rafts through reduced activation of small GTPase Cdc42 and decreased actin polymerization. Changes in rafts led to re-localization of TLR4 and TREM-1 to rafts, phosphorylation of ERK1/2, activation of NLRP3 inflammasome, and increased secretion of pro-inflammatory cytokines. The effects of exNef on lipid rafts and on inflammation were reversed by overexpression of a constitutively active mutant of Cdc42. Similar effects were observed in macrophages treated with exosomes produced by HIV-infected cells or isolated from plasma of HIV-infected subjects, but not with exosomes from cells and subjects infected with ΔNef-HIV or uninfected subjects. Mice injected with exNef exhibited monocytosis, reduced ABCA1 in macrophages, increased raft abundance in monocytes and augmented inflammation. Thus, Nef-containing exosomes potentiated pro-inflammatory response by inducing changes in cholesterol metabolism and reorganizing lipid rafts. These mechanisms may contribute to HIV-associated metabolic co-morbidities.

Analysis of Genomes and Transcriptomes of Hepatocellular Carcinomas Identifies Mutations and Gene Expression Changes in the Transforming Growth Factor-ß Pathway.

BACKGROUND & AIMS: Development of hepatocellular carcinoma (HCC) is associated with alterations in the transforming growth factor-beta (TGF-ß) signaling pathway, which regulates liver inflammation and can have tumor suppressor or promoter activities. Little is known about the roles of specific members of this pathway at specific of HCC development. We took an integrated approach to identify and validate the effects of changes in this pathway in HCC and identify therapeutic targets. METHODS: We performed transcriptome analyses for a total of 488 HCCs that include data from The Cancer Genome Atlas. We also screened 301 HCCs reported in the Catalogue of Somatic Mutations in Cancer and 202 from Cancer Genome Atlas for mutations in genome sequences. We expressed mutant forms of spectrin beta, non-erythrocytic 1 (SPTBN1) in HepG2, SNU398, and SNU475 cells and measured phosphorylation, nuclear translocation, and transcriptional activity of SMAD family member 3 (SMAD3). RESULTS: We found somatic mutations in at least 1 gene whose product is a member of TGF-ß signaling pathway in 38% of HCC samples. SPTBN1 was mutated in the largest proportion of samples (12 of 202, 6%). Unsupervised clustering of transcriptome data identified a group of HCCs with activation of the TGF-ß signaling pathway (increased transcription of genes in the pathway) and a group of HCCs with inactivation of TGF-ß signaling (reduced expression of genes in this pathway). Patients with tumors with inactivation of TGF-ß signaling had shorter survival times than patients with tumors with activation of TGF-ß signaling (P = .0129). Patterns of TGF-ß signaling correlated with activation of the DNA damage response and sirtuin signaling pathways. HepG2, SNU398, and SNU475 cells that expressed the D1089Y mutant or with knockdown of SPTBN1 had increased sensitivity to DNA crosslinking agents and reduced survival compared with cells that expressed normal SPTBN1 (controls). CONCLUSIONS: In genome and transcriptome analyses of HCC samples, we found mutations in genes in the TGF-ß signaling pathway in almost 40% of samples. These correlated with changes in expression of genes in the pathways; up-regulation of genes in this pathway would contribute to inflammation and fibrosis, whereas down-regulation would indicate loss of TGF-ß tumor suppressor activity. Our findings indicate that therapeutic agents for HCCs can be effective, based on genetic features of the TGF-ß pathway; agents that block TGF-ß should be used only in patients with specific types of HCCs.

RNA2DNAlign: nucleotide resolution allele asymmetries through quantitative assessment of RNA and DNA paired sequencing data.

We introduce RNA2DNAlign, a computational framework for quantitative assessment of allele counts across paired RNA and DNA sequencing datasets. RNA2DNAlign is based on quantitation of the relative abundance of variant and reference read counts, followed by binomial tests for genotype and allelic status at SNV positions between compatible sequences. RNA2DNAlign detects positions with differential allele distribution, suggesting asymmetries due to regulatory/structural events. Based on the type of asymmetry, RNA2DNAlign outlines positions likely to be implicated in RNA editing, allele-specific expression or loss, somatic mutagenesis or loss-of-heterozygosity (the first three also in a tumor-specific setting). We applied RNA2DNAlign on 360 matching normal and tumor exomes and transcriptomes from 90 breast cancer patients from TCGA. Under high-confidence settings, RNA2DNAlign identified 2038 distinct SNV sites associated with one of the aforementioned asymetries, the majority of which have not been linked to functionality before. The performance assessment shows very high specificity and sensitivity, due to the corroboration of signals across multiple matching datasets. RNA2DNAlign is freely available from http://github.com/HorvathLab/NGS as a self-contained binary package for 64-bit Linux systems.

A cellular and molecular mosaic establishes growth and differentiation states for cranial sensory neurons.

We compared apparent origins, cellular diversity and regulation of initial axon growth for differentiating cranial sensory neurons. We assessed the molecular and cellular composition of the developing olfactory and otic placodes, and cranial sensory ganglia to evaluate contributions of ectodermal placode versus neural crest at each site. Special sensory neuron populations-the olfactory and otic placodes, as well as those in vestibulo-acoustic ganglion- are entirely populated with cells expressing cranial placode-associated, rather than neural crest-associated markers. The remaining cranial sensory ganglia are a mosaic of cells that express placode-associated as well as neural crest-associated markers. We found two distinct populations of neural crest in the cranial ganglia: the first, as expected, is labeled by Wnt1:Cre mediated recombination. The second is not labeled by Wnt1:Cre recombination, and expresses both Sox10 and FoxD3. These populations-Wnt1:Cre recombined, and Sox10/Foxd3-expressing- are proliferatively distinct from one another. Together, the two neural crest-associated populations are substantially more proliferative than their placode-associated counterparts. Nevertheless, the apparently placode- and neural crest-associated populations are similarly sensitive to altered signaling that compromises cranial morphogenesis and differentiation. Acute disruption of either Fibroblast growth factor (Fgf) or Retinoic acid (RA) signaling alters axon growth and cell death, but does not preferentially target any of the three distinct populations. Apparently, mosaic derivation and diversity of precursors and early differentiating neurons, modulated uniformly by local signals, supports early cranial sensory neuron differentiation and growth.

Gigantism and acromegaly due to Xq26 microduplications and GPR101 mutation.

BACKGROUND: Increased secretion of growth hormone leads to gigantism in children and acromegaly in adults; the genetic causes of gigantism and acromegaly are poorly understood. METHODS: We performed clinical and genetic studies of samples obtained from 43 patients with gigantism and then sequenced an implicated gene in samples from 248 patients with acromegaly. RESULTS: We observed microduplication on chromosome Xq26.3 in samples from 13 patients with gigantism; of these samples, 4 were obtained from members of two unrelated kindreds, and 9 were from patients with sporadic cases. All the patients had disease onset during early childhood. Of the patients with gigantism who did not carry an Xq26.3 microduplication, none presented before the age of 5 years. Genomic characterization of the Xq26.3 region suggests that the microduplications are generated during chromosome replication and that they contain four protein-coding genes. Only one of these genes, GPR101, which encodes a G-protein-coupled receptor, was overexpressed in patients' pituitary lesions. We identified a recurrent GPR101 mutation (p.E308D) in 11 of 248 patients with acromegaly, with the mutation found mostly in tumors. When the mutation was transfected into rat GH3 cells, it led to increased release of growth hormone and proliferation of growth hormone-producing cells. CONCLUSIONS: We describe a pediatric disorder (which we have termed X-linked acrogigantism [X-LAG]) that is caused by an Xq26.3 genomic duplication and is characterized by early-onset gigantism resulting from an excess of growth hormone. Duplication of GPR101 probably causes X-LAG. We also found a recurrent mutation in GPR101 in some adults with acromegaly. (Funded by the Eunice Kennedy Shriver National Institute of Child Health and Human Development and others.).

SNPlice: variants that modulate Intron retention from RNA-sequencing data.

RATIONALE: The growing recognition of the importance of splicing, together with rapidly accumulating RNA-sequencing data, demand robust high-throughput approaches, which efficiently analyze experimentally derived whole-transcriptome splice profiles. RESULTS: We have developed a computational approach, called SNPlice, for identifying cis-acting, splice-modulating variants from RNA-seq datasets. SNPlice mines RNA-seq datasets to find reads that span single-nucleotide variant (SNV) loci and nearby splice junctions, assessing the co-occurrence of variants and molecules that remain unspliced at nearby exon-intron boundaries. Hence, SNPlice highlights variants preferentially occurring on intron-containing molecules, possibly resulting from altered splicing. To illustrate co-occurrence of variant nucleotide and exon-intron boundary, allele-specific sequencing was used. SNPlice results are generally consistent with splice-prediction tools, but also indicate splice-modulating elements missed by other algorithms. SNPlice can be applied to identify variants that correlate with unexpected splicing events, and to measure the splice-modulating potential of canonical splice-site SNVs. AVAILABILITY AND IMPLEMENTATION: SNPlice is freely available for download from https://code.google.com/p/snplice/ as a self-contained binary package for 64-bit Linux computers and as python source-code. CONTACT: pmudvari@gwu.edu or horvatha@gwu.edu SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

A genome-wide scan identifies mutations in the gene encoding phosphodiesterase 11A4 (PDE11A) in individuals with adrenocortical hyperplasia.

Phosphodiesterases (PDEs) regulate cyclic nucleotide levels. Increased cyclic AMP (cAMP) signaling has been associated with PRKAR1A or GNAS mutations and leads to adrenocortical tumors and Cushing syndrome. We investigated the genetic source of Cushing syndrome in individuals with adrenocortical hyperplasia that was not caused by known defects. We performed genome-wide SNP genotyping, including the adrenocortical tumor DNA. The region with the highest probability to harbor a susceptibility gene by loss of heterozygosity (LOH) and other analyses was 2q31-2q35. We identified mutations disrupting the expression of the PDE11A isoform-4 gene (PDE11A) in three kindreds. Tumor tissues showed 2q31-2q35 LOH, decreased protein expression and high cyclic nucleotide levels and cAMP-responsive element binding protein (CREB) phosphorylation. PDE11A codes for a dual-specificity PDE that is expressed in adrenal cortex and is partially inhibited by tadalafil and other PDE inhibitors; its germline inactivation is associated with adrenocortical hyperplasia, suggesting another means by which dysregulation of cAMP signaling causes endocrine tumors.

ScSNViz: a user-friendly toolset for visualization and analysis of Cell-Specific Expressed SNVs.

Motivation: Understanding genetic variation at the single-cell level is crucial for insights into cellular heterogeneity, clonal evolution, and gene expression regulation, but there is a scarcity of tools for visualizing and analyzing cell-level genetic variants. Results: We introduce scSNViz, a comprehensive R-based toolset for visualization and analysis of cell-specific expressed Single Nucleotide Variants (sceSNVs) within cell-barcoded single-cell RNA-sequencing (scRNA-seq) data. ScSNViz offers 3D sceSNV visualization capabilities for dimensionally reduced scRNA-seq gene expression data, compatibility with popular scRNA-seq processing tools like Seurat, cell-type classification tools such as SingleR and scType, and trajectory inference computation using Slingshot. Furthermore, scSNViz conducts estimation, summary, and graphical representation of statistical metrics pertaining to sceSNVs distribution and expression across individual cells. It also provides support for the analysis of individual sceSNVs as well as sets comprising multiple expressed sceSNVs of interest. Availability: ScSNViz is implemented as user-friendly R-scripts, freely available on https://horvathlab.github.io/NGS/scSNViz , supported by help utilities, and requiring no specialized bioinformatics skills for use.

TRIM28 promotes luminal cell plasticity in a mouse model of prostate cancer.

The Tripartite motif-containing 28 (TRIM28) transcriptional cofactor is significantly upregulated in high-grade and metastatic prostate cancers. To study the role of TRIM28 in prostate cancer progression in vivo, we generated a genetically-engineered mouse model, combining prostate-specific inactivation of Trp53, Pten and Trim28. Trim28 inactivated NPp53T mice developed an inflammatory response and necrosis in prostate lumens. By conducting single-cell RNA sequencing, we found that NPp53T prostates had fewer luminal cells resembling proximal luminal lineage cells, which are cells with progenitor activity enriched in proximal prostates and prostate invagination tips in wild-type mice with analogous populations in human prostates. However, despite increased apoptosis and reduction of cells expressing proximal luminal cell markers, we found that NPp53T mouse prostates evolved and progressed to invasive prostate carcinoma with a shortened overall survival. Altogether, our findings suggest that TRIM28 promotes expression of proximal luminal cell markers in prostate tumor cells and provides insights into TRIM28 function in prostate tumor plasticity.

BRCA1 deficiency in mature CD8+ T lymphocytes impairs antitumor immunity.

Women with BRCA1 germline mutations have approximately an 80% lifetime chance of developing breast cancer. While the tumor suppressor function of BRCA1 in breast epithelium has been studied extensively, it is not clear whether BRCA1 deficiency in non-breast somatic cells also contribute to tumorigenesis. Here, we report that mouse Brca1 knockout (KO) in mature T lymphocytes compromises host antitumor immune response to transplanted syngeneic mouse mammary tumors. T cell adoptive transfer further corroborates CD8+ T cell-intrinsic impact of Brca1 KO on antitumor adaptive immunity. T cell-specific Brca1 KO mice exhibit fewer total CD8+, more exhausted, reduced cytotoxic, and reduced memory tumor-infiltrating T cell populations. Consistent with the preclinical data, cancer-free BRCA1 mutation-carrying women display lower abundance of circulating CD8+ lymphocytes than the age-matched control group. Thus, our findings support the notion that BRCA1 deficiency in adaptive immunity could contribute to BRCA1-related tumorigenesis. We also suggest that prophylactic boosting of adaptive immunity may reduce cancer incidence among at-risk women.

Prostate cancer cell-platelet bidirectional signaling promotes calcium mobilization, invasion and apoptotic resistance via distinct receptor-ligand pairs.

Platelets play a crucial role in cancer and thrombosis. However, the receptor-ligand repertoire mediating prostate cancer (PCa) cell-platelet interactions and ensuing consequences have not been fully elucidated. Microvilli emanating from the plasma membrane of PCa cell lines (RC77 T/E, MDA PCa 2b) directly contacted individual platelets and platelet aggregates. PCa cell-platelet interactions were associated with calcium mobilization in platelets, and translocation of P-selectin and integrin αIIbß3 onto the platelet surface. PCa cell-platelet interactions reciprocally promoted PCa cell invasion and apoptotic resistance, and these events were insensitive to androgen receptor blockade by bicalutamide. PCa cells were exceedingly sensitive to activation by platelets in vitro, occurring at a PCa cell:platelet coculture ratio as low as 1:10 (whereas PCa patient blood contains 1:2,000,000 per ml). Conditioned medium from cocultures stimulated PCa cell invasion but not apoptotic resistance nor platelet aggregation. Candidate transmembrane signaling proteins responsible for PCa cell-platelet oncogenic events were identified by RNA-Seq and broadly divided into 4 major categories: (1) integrin-ligand, (2) EPH receptor-ephrin, (3) immune checkpoint receptor-ligand, and (4) miscellaneous receptor-ligand interactions. Based on antibody neutralization and small molecule inhibitor assays, PCa cell-stimulated calcium mobilization in platelets was found to be mediated by a fibronectin1 (FN1)-αIIbß3 signaling axis. Platelet-stimulated PCa cell invasion was facilitated by a CD55-adhesion G protein coupled receptor E5 (ADGRE5) axis, with contribution from platelet cytokines CCL3L1 and IL32. Platelet-stimulated PCa cell apoptotic resistance relied on ephrin-EPH receptor and lysophosphatidic acid (LPA)-LPA receptor (LPAR) signaling. Of participating signaling partners, FN1 and LPAR3 overexpression was observed in PCa specimens compared to normal prostate, while high expression of CCR1 (CCL3L1 receptor), EPHA1 and LPAR5 in PCa was associated with poor patient survival. These findings emphasize that non-overlapping receptor-ligand pairs participate in oncogenesis and thrombosis, highlighting the complexity of any contemplated clinical intervention strategy.

Embryonic stem cell factor FOXD3 (Genesis) defects in gastrointestinal stromal tumors.

Gastrointestinal stromal tumors (GISTs) are mesenchymal neoplasms, believed to originate from the interstitial cells of Cajal (ICC), often caused by overexpression of tyrosine kinase receptors (TKR) KIT or PDGFRA. Here, we present evidence that the embryonic stem cell factor FOXD3, first identified as 'Genesis' and involved in both gastrointestinal and neural crest cell development, is implicated in GIST pathogenesis; its involvement is investigated both in vitro and in zebrafish and a mouse model of FOXD3 deficiency. Samples from a total of 58 patients with wild-type GISTs were used for molecular analyses, including Sanger sequencing, comparative genomic hybridization, and methylation analysis. Immunohistochemistry and western blot evaluation were used to assess FOXD3 expression. Additionally, we conducted in vitro functional studies in tissue samples and in transfected cells to confirm the pathogenicity of the identified genetic variants. Germline partially inactivating FOXD3 sequence variants (p.R54H and p.Ala88_Gly91del) were found in patients with isolated GISTs. Chromosome 1p loss was the most frequent chromosomal abnormality identified in tumors. In vitro experiments demonstrate the impairment of FOXD3 in the presence of those variants. Animal studies showed disruption of the GI neural network and changes in the number and distribution in the ICC. FOXD3 suppresses KIT expression in human cells; its inactivation led to an increase in ICC in zebrafish, as well as mice, providing evidence for a functional link between FOXD3 defects and KIT overexpression leading to GIST formation.