Genetic variants in canonical Wnt signaling pathway associated with pediatric immune thrombocytopenia

Pediatric immune thrombocytopenia (ITP) is the most common acquired bleeding disorder in childhood, with ∼3000 to 4000 new cases diagnosed per year.¹ Despite the first classical case of ITP being described in 1735 by Paul Gottleib Werlof, substantial gaps in our understanding of its pathogenesis remain.² Pediatric ITP is thought to develop from an interplay of a triggering event and subsequent immune activation in the context of predisposing genetic factors.³ Particularly in chronic, refractory, or multilineage immune cytopenias, causative inborn errors of immunity are increasingly being identified.^1,4-8 

The exact immune mechanism driving disease is identified for few cases of pediatric ITP. Children with symptomatic or prolonged courses of ITP are often treated as a monolith and broad immune suppressive medications are applied with variable efficacy and disease control. Fortunately, the landscape of ITP care is evolving. Genetic testing has led to the identification of new disorders that predispose to ITP. Targeted therapies are now available for a subset of patients with germ line causes of immune cytopenias, for example, the use of abatacept for patients with CTLA-4 haploinsufficiency or Lipopolysaccharide-responsive and beige-like anchor (LRBA) deficiency.^4,9,10 With increased accessibility and decreased cost of genetic sequencing, there is substantial opportunity to improve ITP care by tailoring treatments to individual patients.

Beyond the potential to identify genetic etiologies of ITP and refine therapeutic choices, genetic testing holds the promise of improving global understanding of the ITP pathophysiology. The scope of ITP pathophysiology extends well beyond early descriptions of autoantibody–mediated platelet destruction, with pathologic mechanisms now understood to involve a multitude of immune abnormalities including defects in the megakaryocyte compartment^11-13; shifts toward proinflammatory T-cell responses^14-16; loss of regulatory immune populations^17-21; direct cytotoxic T-cell destruction^22-24; loss of sialic acid residues^13,25,26; defects in innate immunity^27-30; and, recently, complement-mediated platelet destruction.^31-33 With the goal of better understanding the genetic predisposition for ITP and to shed light on pathophysiology, we performed genotyping on DNA samples obtained from a large multi-institutional cohort of pediatric patients with ITP.

Pediatric patients with ITP diagnosed at large tertiary referral care centers across the United States were enrolled in the Pediatric ITP Consortium of North America’s biorepository from August 2015 to March 2021. Informed consent was obtained from all patients and/or their parents or legal guardians, and assent was obtained when applicable. In addition to collecting demographic, clinical, and laboratory data, biologic specimens were obtained from enrolled patients. DNA was extracted from peripheral blood samples and preserved in the central Pediatric ITP Consortium of North America’s biorepository. Patients with secondary ITP, Evans syndrome and ITP of all disease phases (newly diagnosed, persistent, and chronic) were included in the biorepository and this analysis.³⁴ Eight centers, representing diverse populations across the United States (Boston Children’s Hospital, Children’s Hospital of Philadelphia, Cincinnati Children’s Hospital, Children’s Hospital Colorado, Columbia University, Nationwide Children’s Hospital, Texas Children’s Hospital, and University of California San Francisco), contributed samples used in this study.

Reporting race and ethnicity in this study was mandated by the US National Institutes of Health, consistent with the Inclusion of Women, Minorities, and Children policy. Individuals participating in the study were categorized as Asian, Black or African American, Hispanic or Latino, other/unknown, or White based on the US National Institutes of Health Policy on Reporting Race and Ethnicity Data. Data were based on the parents’ report. Sex was reported as biologic (chromosomal) sex.

Patients with chronic ITP (cITP) were defined as having ITP for >1 year. Spontaneously resolving ITP (srITP) was defined as disease that self-resolved within 1 year. Disease resolution was defined as having 2 consecutive platelet counts of >150 × 10⁹/L at least 1 month apart, with a mean platelet volume (MPV) of ≤11 fL, or a single platelet count of >200× 10⁹/L with a normal MPV or ≥3 consecutive platelet counts of >150 × 10⁹/L with an MPV of ≤11.5 fL. For all the prior categories, the referenced platelet counts were only considered if they were obtained absent of ITP-directed therapy, which was defined as: 6 weeks since IV immunoglobulin or Rh immune globulin, 3 weeks since completion of corticosteroid treatment, off immunomodulators and thrombopoietin receptor agonists for at least 1 week, or after achieving B-cell repopulation for those who received rituximab therapy.

Samples were genotyped using the Illumina Global Screening Array BeadChip version 1.0. The global screening array includes >650 000 variants and is designed to tag common and low-frequency variants (minor allele frequency [MAF] of >1%) in global populations. Genotype calls were generated using Illumina’s GenomeStudio software. Genotype and limited demographic data on otherwise healthy pediatric controls (n = 8,372) were obtained from the National Longitudinal Study of Adolescent to Adult Health (Add Health) cohort,³⁵ which used 2 Illumina platforms for genotyping: the Illumina Human Omni1-Quad BeadChip for the majority of samples, and the Illumina Omni-2.5 Quad BeadChip. Quality control checks of the genotypes (for both ITP cases and Add Health controls) were performed in Plink 1.9.³⁶ Single nucleotide variants (SNVs) with a call rate of <95%, a Hardy-Weinberg equilibrium P value of <5 × 10⁻⁵, or MAF of <1% were excluded. Individuals with sample call rates of <95%, heterozygosity rates greater than ±2.5 standard deviations from the array median, or evidence of cryptic relatedness (identity by decent, P > .2) were excluded from the analysis.³⁷ 

Whole-genome level SNV information was obtained through imputation of disomic autosomal chromosomes through the Michigan imputation server.³⁸ Phasing of genotyped SNVs, which passed quality control, were completed with Eagle version 2.4 before imputation with Minimac3 using the Haplotype Reference Consortium r1.1 as reference panel.³⁹ After imputation, genetic variants with a MAF of ≥0.01, estimated imputation quality score (R2) of >0.30, and Hardy-Weinberg equilibrium P > 1.0 × 10⁻⁴ were retained for association analyses. After quality control, 476,159 disomic SNVs were kept in the analysis.

Logistic regression models were used to estimate adjusted odds ratios (aORs) comparing (1) ITP and control cohorts, and (2) cITP and srITP cohorts. To account for population structure, principal component analysis was performed on the imputed data using Plink 1.9. Models for autosomal variants assumed additive allelic effects and accounted for the top 3 principal components, based on optimizing λ from Q-Q plots. Models for variants on the X chromosome were evaluated in 3 ways using standard analysis options in Plink: (1) including sex as a covariate and coding male X chromosome genotypes 0 of 2 instead of 0 of 1 (model option 2 in Plink); (2) including sex as a covariate and also including a dosage–sex interaction term (model option 3 in Plink); and (3) conduct stratified analysis based on sex and report results separately for males and females. Sensitivity analyses were performed to assess the robustness of the association results and included adjustment for additional principal components and covariates. Genome-wide significance for all analyses was defined as P < 5 × 10⁻⁸ and suggestive genome-wide significance, P < 1 × 10⁻⁵. All statistical analyses were conducted using Plink 1.9, VCFTools 0.1.15, and R Studio 4.2.2.

Potential functional impacts of top variants were evaluated using 3 primary tools: the combined annotation dependent depletion tool,⁴⁰ ENSEMBL variant effect predictor,⁴¹ and the get_qtl function in the “Qtlizer” R package.⁴² Colocalization analyses were conducted for top suggestive expression quantitative trait loci (eQTL) using ezQTL (https://analysistools.cancer.gov/ezqtl/).

Overall, 610 unique samples from patients with ITP were genotyped; 19 (3.1%) did not pass quality control, leaving 591 for analysis. This cohort included 307 patients with cITP (52.0%), 209 patients with srITP (35.4%), and 75 (12.7%) individuals without complete follow-up data to allow for classification of cITP or srITP. Demographic data for the srITP and cITP cohort are outlined in Tables 1 and 2. The ITP cohort was 50% male compared with 47% male among the controls (Table 1). The ITP cohort was self-reported 50.9% non-Hispanic White, 29.4% Hispanic, 6.3% Black, 3.0% Asian, and 4.9% additional racial/ethnic groups. The Add Health controls were self-reported 61.5% non-Hispanic White, 14.8% Hispanic, 22.5% Black, 0.8% Asian, and 1.12% additional racial/ethnic groups. Sex, race, and ethnicity were similar in the srITP and cITP cohorts (P = .602, P = .254, and P = .054, respectively; Table 2).

Compared with Add Health controls, 6 loci were associated with ITP at genome-wide significance (Figure 1: Manhattan plot; Figure 2: forest plot case-control analysis). No variants on the X chromosome reached significance or suggestive significance by any of the 3 aforementioned approaches (supplemental Figure 1: consolidated Manhattan plots for X chromosome). An intronic variant, rs2165798, in NAV2 was inversely associated with ITP (aOR, 0.52; P = 3.2 × 10⁻¹¹). Immediately upstream of NKD1, rs4785216 was also inversely associated with ITP (aOR, 0.43; P = 8.86 × 10⁻¹⁵). Another intergenic variant, rs4785426, in close proximity to NKD1, was also inversely associated with ITP (aOR, 0.48; P = 1.84 × 10⁻¹⁶). Additionally, rs75650997, located within the intronic region of SLC24A4 was associated with an increased risk of ITP compared with controls (aOR, 2.38; P = 3.9 × 10⁻⁸). Two additional variants reached genome-wide significance; 1 in the noncoding RNA LOC105375957 on chromosome 9 (rs10967515; aOR, 0.40; P = 2.99 × 10⁻¹⁵) and 1 in an intergenic region on chromosome 16 (rs8054390; aOR, 3.98; P = 1.25 × 10⁻⁸). An additional 20 variants met suggestive genome-wide significance (Tables 3 and 4).

Pairwise linkage disequilibrium (LD) was evaluated among all SNVs that met genome-wide significance and suggestive level of significance. On chromosome 16, there was a weak association between rs4785426 and rs4785216 (r² = 0.11; D’ = 0.358), however, there were no other SNVs of interest in high LD. Further analysis of the LD for these SNVs across the 32 1000 Genomes Project populations available in LD Link, showed a range of r² values from 0.0067 to 0.3645, suggesting low levels of LD for these variants across populations. The highest r² values were among various European populations (ie, Toscani in Italia [TSI], British in England and Scotland [GBR], and Iberian population in Spain [IBS]).

All of the top variants are fairly conserved according to both the PhyloP and Genomic Evolutionalry Rating Profiling (GERP) scores, except rs75650997, which shows evidence of potentially being fast evolving. Two of the variants were predicted to be in regulatory regions; rs75650997 in the promoter region for SLC24A4I, and rs4785426 in an enhancer region for NKD1. The combined annotation-dependent depletion Phred scores for 3 of the variants indicate that they are likely to be deleterious: rs2165798 in NAV2 (Phred = 16.09), rs4785426 in NKD1 (PHRED = 17.46), and rs478516 in NKD1 (Phred = 21.1). This indicates that rs2165798 and rs4785426 are predicted to be among the 10% most deleterious substitutions that could be made to the human genome (ie, a score of ≥10) and that rs478516 is among the 1% most deleterious (ie, a score of ≥20). Additionally, 3 of the top 6 variants show some evidence of being in relevant eQTL regions: rs4785426 for ADCY7 (in the peripheral blood), rs75650997 for SLC24A4 (in the whole blood), and rs8054390 for C16orf95 (in the brain frontal cortex). Full output for the eQTL analysis of the top hits is provided in supplemental Table 1. Suggestive eQTLs on chromosome 14 and chromosome 16 were further analyzed using the colocalization analysis in ezQTL. The colocalization plots for each are provided in supplemental Figure 2.

No SNV reached genome-wide significance, but 1 variant, rs117609649, in ABCC4, met suggestive genome-wide significance (MAF srITP vs cITP: 10.5% vs 2.9%; aOR, 0.26; P = 4.49 × 10⁻⁷).

This analysis is, to our knowledge, the largest, unbiased genetic study of pediatric patients with ITP. The analysis, comparing pediatric ITP vs healthy controls, identified 6 SNVs meeting genome-wide significance and 20 meeting suggestive genome-wide significance. No genetic variants met genome-wide significance in a comparison of patients with srITP vs cITP. Collectively, these findings suggest there may be genetic variants implicated in ITP risk, but those genetic features are not specific to srITP or cITP.

The SNVs identified in this study, whose frequencies were significantly different in ITP vs control populations, were common with MAF of >1%. Importantly, the etiology of ITP is complex. Although there are known monogenic immunodeficiency syndromes with features of immune cytopenias, it is probable that some ITP cases develop in the context of a triggering event and individual immune milieu, in concert with genetic predisposition. These variants may increase ITP risk but may not be the sole cause of an individual patient’s disease.

The significant SNVs identified in comparison of ITP and controls may indicate a novel molecular pathway contributing to ITP pathogenesis. Three of the SNVs meeting genome-wide significance and 2 SNVs suggestive of genome-wide significance identified in this study are in genes that affect Wnt/β-catenin signaling or the canonical Wnt signaling pathway. This pathway is involved in cell proliferation, differentiation, and migration, which are mediated through T-cell factor/lymphoid enhancer-binding factor transcription factors, that activate target genes.⁴³ 

NAV2 may play a role in cellular growth and migration.^44,NAV2 variants have been implicated in autoimmune conditions including Celiac disease and rheumatoid arthritis.^45,46 Furthermore, 1 study has demonstrated that NAV2 activates Wnt/β-catenin signaling in the setting of immune-mediated arthritis.^44,NKD1 passively inhibits the canonical Wnt/β-catenin signaling pathway.^47,48 

Of the SNVs that met the threshold for suggestive genome-wide significance, CCDC85A, whose presence conferred a nearly 4-times increased risk for ITP, is also linked to the Wnt/β-catenin pathway. One study postulated that CCDC85A may play a role in cell-cell adhesion through interaction with β-catenin proteins.⁴⁹ An intronic variant in GRIP1 was present >2.5 times more often in patients with ITP vs healthy controls. GRIP1, as a component of a coactivator molecule, binds β-catenin, which then binds T-cell–specific transcription factor/lymphoid enhancer-binding factor and activates target gene transcription.^50,51 

Identified in the case-with-case comparison, ABCC4, also known as MRP4, is highly expressed by platelets and regulates intracellular cyclic adenosine monophosphate and guanosine 3′,5′-cyclic monophosphate.^52,53 It contributes to platelet activation and function, with variants in ABCC4 being implicated in δ-granule disorders.⁵³ It is unclear how this variant could be affecting srITP or cITP but may reflect an increased risk for aberrant platelet function in ITP. Interestingly, several studies in prostate, endometrial, and breast cancer models suggest that ATP binding cassette subfamily C member 4)/Multidrug resistance protein 4 (ABCC4/MRP4) activity correlates β-catenin levels and signaling.^54-56 

Wnt/β-catenin signaling has been implicated in platelet function, as well as the development, proliferation, and maturation of megakaryocytes.^57-59 Electron microscopy studies have long since established that megakaryocyte maturation and morphology are abnormal in ITP.¹¹ Potentially Wnt/β-catenin function is a contributor to abnormalities in the megakaryocyte compartment in ITP. Wnt-signaling affects the development and regulation of multiple immune cell populations.⁶⁰ Wnt/β-catenin upregulation may induce dendritic cell activity that then induces T regulatory cells, promote anti-inflammatory cytokine signaling, and inhibit T helper type 1/T helper type 17 cells and cytotoxic CD8⁺ T cells.^61,62 One study suggested that Wnt/β-catenin may directly regulate T regulatory cell function through modulating the activity of Foxp3, a lineage defining transcription factor.⁶² 

Based on our genomic findings, we postulate alterations in Wnt/β-catenin signaling could affect T regulatory cell populations and may result in loss of immune tolerance. However, the exact impact of these SNVs on Wnt/β-catenin signaling is unknown and should be a focus of future ITP study. Potentially, the findings of this study could lead to novel therapies for ITP; multiple agents targeting Wnt/β-catenin signaling in cancers are in development.⁶³ 

First, despite being, to our knowledge, the largest cohort of pediatric patients with ITP to undergo broad genotyping analysis, because pediatric ITP is a rare disease, the cohort remains small for a genome-wide analysis study. As such, variants with relevance in ITP pathogenesis could be missed, especially when comparing cITP with srITP. Additionally, the DNA microarray used is designed to capture common and low frequency SNVs but not rare variants (MAF of <1%). Variants that are present at <1% frequency or private to a single individual may be drivers of ITP pathogenesis but would be missed by this microarray. This could explain the lack of variants identified in genes implicated in inborn errors of immunity, which have been linked previously to development of ITP and Evans syndrome.^4-8,64 Next-generation sequencing methods could better identify rare variant associations with ITP. The DNA microarray does cover HLA alleles and killer cell immunoglobulin-like receptors (KIR) genes. However, the most rigorous analysis of this region should use different analytic approaches because of the presence of short, repetitive sequence in this region.^65-68 The HLA region and KIR genes have long been a genomic area of interest in autoimmune disease, including ITP. KIR expression patterns may also be relevant to ITP but would not be captured in a genotyping assay.²⁹ The results of this study does not rule out their role in ITP risk.^29,69-72 

Overall, the few SNVs strongly linked to cITP or ITP emphasize the need to study other areas of the genome that have not been well interrogated by this DNA microarray analysis.

Candidate variants should be validated in an independent patient cohort. If validated, evaluation of Wnt/β-catenin signaling activity in patients with ITP should be performed to understand the functional relevance of the variants identified in this study. We advise further sequencing of HLA and KIR genes in addition to next-generation sequencing given the potential limitations of the DNA microarray.

The clinical course of pediatric ITP is highly heterogenous, and its molecular pathophysiology is likely equally complex. This is, to our knowledge, the largest pediatric ITP cohort to undergo genotyping and the first to compare children with cITP with those with self-resolving disease. We identified SNVs meeting genome-wide significance potentially associated with the development of ITP but not specific to srITP or cITP. Variants in genes that affect the canonical Wnt-signaling pathway may play a role in ITP risk. This represents a novel pathway that should be further evaluated in patients with ITP.

The authors acknowledge the hard work of Pediatric ITP Consortium of North America’s (ICON) coinvestigators and research teams at contributing ICON centers, as well as the contribution of our patients and their families.

Funding for the Pediatric ITP Consortium biorepository was provided by Baylor College of Medicine. This research uses data from Add Health, funded by grant P01 HD31921 (Harris) from the Eunice Kennedy Shriver National Institute of Child Health and Human Development, with cooperative funding from 23 other federal agencies and foundations. Add Health is currently directed by Robert A. Hummer and funded by the National Institute on Aging cooperative agreements U01 AG071448 (Hummer) and U01AG071450 (Aiello and Hummer) at the University of North Carolina at Chapel Hill. Add Health was designed by J. Richard Udry, Peter S. Bearman, and Kathleen Mullan Harris at the University of North Carolina at Chapel Hill.³⁵ 

Contribution: T.O.K. designed the research, enrolled patients, prepared samples, interpreted results, and cowrote the manuscript; J.M.G. planned and executed statistical analysis, analyzed, interpreted results, and cowrote the manuscript; J.M.F. designed the research, and prepared and processed samples; C.O. collected data and cowrote the manuscript; R.F.G., M.P.L., M.J.R., K.A.S., O.N., C.N., T.A.N., and S.E.K. enrolled patients, collected samples, and reviewed and edited the manuscript; D.M. planned the research, interpreted results, and cowrote the manuscript; B.D. enrolled patients and prepared samples; J.M.D., E.J.N., M.P.L., and R.F.G. designed the research; A.B.G. enrolled patients, interpreted results, and cowrote the manuscript; and M.E.S. interpreted results and cowrote the manuscript.

Conflict-of-interest disclosure: R.F.G. reports research funding from Agios, Novartis, and Sobi; and reports consultancy role with Agios and Sanofi. M.P.L. reports advisory board membership with Octapharma, Dova, Principia, Rigel, Argenx, Platelet Disorder Support Association (PDSA), 22qSociety, and CdLS Foundation; serves as a consultant for Novartis, Dova, Principia, Argenx, Rigel, Sobi, Sanofi, and Janssen; and reports research funding from Foundation for Women & Girls with Blood Disorders, PDSA, the US National Institutes of Health, Sysmex, Novartis, Principia, Argenx, Dova, Octapharma, and Sanofi. K.A.S. reports research funding from Novartis, Pfizer, Sobi, and Sanofi. C.N. reports research funding from Novartis; and received an honorarium from Sobi, Novartis, and Sanofi. S.E.K. serves on the advisory board of Sanofi; and reports honoraria from X4 Pharmaceuticals and BioMarin. D.M. reports salary support in part of Division of Clinical Research at National Institute of Allergy and Infectious Disease/the US National Institutes of Health. E.J.N. serves on data and safety monitoring board for Merck, LFB, Sobi, and Agios; serves on the advisory boards of Genentech, Takeda, Saliogen, and Pfizer; and reports stock options in Saliogen. A.B.G. reports research funding from Novartis. The remaining authors declare no competing financial interests.

Correspondence: Taylor Olmsted Kim, Children’s Hospital Los Angeles, University of Southern California Keck School of Medicine, Department of Pediatrics, Cancer and Blood Disease Institute, 4650 Sunset Blvd, Los Angeles, CA 90027; email: tkim@chla.usc.edu.