Transcriptome profiling of megakaryocytes and platelets: Application to GP9- and IKZF5-related thrombocytopenia

Abstract

Platelets are anucleate cells produced in the bone marrow and derived from large progenitor cells called megakaryocytes (MKs). Platelets receive RNA transcripts from their progenitorial MKs during thrombopoiesis. However, the correspondence between platelet and MK transcriptomes is poorly understood, particularly in the context of germline mutations that cause platelet formation defects or thrombocytopenia. We have studied the effects of two such mutations on MK and platelet transcriptomes. We generated immortalized MK cell line (imMKCL)-based models of Bernard-Soulier syndrome and IKZF5-related thrombocytopenia. MKs derived from imMKCLs with either a homozygous deletion of GP9 (GP9 ^-/-) or a heterozygous Y121F variant in IKZF5 (IKZF5 ^WT/Y121F) exhibited reduced proplatelet formation (reductions of 96% and 57%, respectively). Platelets from patients with either GP9 ^-/- or IKZF5 ^WT/Y121F genotypes had broad transcriptomic dysregulation, suggesting that (pro)platelet formation defects due to mutations in glycoprotein receptor and transcription factor genes such as GP9 and IKZF5 already affect the MK transcriptome. RNA-seq data from MKs at four stages of differentiation revealed widespread but distinct changes in expression over time between the GP9 ^-/- and the IKZF5 ^WT/Y121F genotypes. Dysregulated genes in GP9 ^-/- MKs were enriched for RNA metabolism and actin/tubulin folding pathways, whereas those in IKZF5 ^WT/Y121F MKs were enriched for cell cycle pathways. Most of these genes were also dysregulated in the platelets of patients with the corresponding diseases. Our results suggest that patients with inherited forms of thrombocytopenia present with specific transcriptomic changes during platelet formation.

Figure 1

Schematic of the cell types used in the transcriptomic study of a glycoprotein receptor knockout (GP9 ^{− / −}) and a transcription factor defect (IKZF5 ^WT/Y121F). Wild‐type (WT) and mutant megakaryocytes (MKs) were harvested at Days 0, 1, 4, and 5 of differentiation and compared to platelets derived from peripheral blood from study participants. imMKCL, immortalized megakaryocyte cell line.

Figure 2

Functional characterization of GP9 ^{− / −} and wild‐type (WT) immortalized megakaryocyte cell lines (imMKCLs). (A) Western blot of WT, GP9KO1, and GP9KO2 imMKCL lysates on Day 0 and Day 4 of differentiation using antibodies for GPIbɑ (145 kDa), GPIX (22 kDa), and glyceraldehyde‐3‐phosphate‐dehydrogenase (GAPDH) (36 kDa) as a loading control. Dead or apoptotic cells were excluded using 4′,6‐diamidino‐2‐phenylindole (DAPI) staining. (B) Representative flow cytometry measurements of WT and GP9 ^{− /−} megakaryocytes (MKs) using markers for CD41a and CD42a. (C) Box plots showing the fraction of positive cells for markers CD41a and CD42a in the two GP9 ^{− / −} Day 4 differentiated MKs and in Day 4 differentiated WT MKs. (D) Histogram of fluorescent intensity measured by flow cytometry of propidium iodide (PI) stained Day 5 differentiated WT and GP9 ^{− / −} MKs and quantification of peaks corresponding to different levels of ploidy. Dead or apoptotic cells were excluded using DAPI staining. (E) Relative proportions of ploidy states in WT, GP9KO1, and GP9KO2 MKs on Day 5 of differentiation measured in duplicate. The vectors of proportions did not differ significantly between the cell lines (all P > 0.05, pairwise comparisons of Dirichlet regression models using likelihood ratio tests [LRTs]). (F) Representative images of proplatelet‐forming WT and GP9 ^{− / −} MKs on Day 6 of differentiation. (G) Fraction of proplatelet‐forming WT and GP9 ^{− / −} MKs obtained from three wells and three separate differentiations on Day 6 of differentiation.

Figure 3

Functional characterization of IKZF5 ^WT/Y121F and wild‐type (WT) immortalized megakaryocyte cell lines (imMKCLs). (A) Representative flow cytometry measurements of WT and Y121FHet megakaryocytes (MKs) using markers for CD41a and CD42a. Dead or apoptotic cells were excluded using 4′,6‐diamidino‐2‐phenylindole (DAPI) staining. (B) Box plots showing the fraction of positive cells for markers CD41a and CD42a in the two Y121FHet Day 4 differentiated MKs and in Day 4 differentiated WT MKs. The mean fractions did not differ by cell line (all P > 0.05, Wilcoxon rank sum test between pairs of vectors). (C) Histogram of fluorescent intensity measured by flow cytometry of propidium iodide (PI) stained Day 5 differentiated WT and IKZF5 ^WT/Y121F MKs and quantification of peaks corresponding to different levels of ploidy. Dead or apoptotic cells were excluded using DAPI staining. (D) Relative proportions of ploidy states in WT, Y121FHet1, and Y121FHet2 MKs on Day 5 of differentiation measured in duplicate. The vectors of proportions did not differ significantly between the cell lines (all P > 0.05, pairwise comparisons of Dirichlet regression models using likelihood ratio tests [LRTs]). (E) Representative images of proplatelet‐forming WT and Y121FHet MKs on Day 6 of differentiation. (F) Fraction of proplatelet‐forming WT and Y121FHet MKs obtained from three wells and three separate differentiations on Day 6 of differentiation.

Figure 4

Principal component analysis (PCA) of megakaryocyte (MK) and platelet transcriptomes. (A) Scatter plot of the first two principal components obtained by performing PCA on the 500 genes with the greatest variability of expression in wild‐type (WT) and mutant MKs. Control and patient platelet transcriptomes were projected onto the plane. (B) Gene set enrichment analysis of the genes contributing to the first principal component with a positive loading. The P‐values corresponding to the 20 most significantly enriched pathways are shown as bars. The number of genes overlapping each pathway is shown in brackets. ECM, extra cellular matrix; imMK, immortalized megakaryocyte; VWF, von Willebrand factor.

Figure 5

Polytomous model selection among five statistical models capturing different patterns of gene expression in GP9 ^{− / −} and wild‐type (WT) megakaryocytes (MKs) (as described in Supporting Information S1: Figure 3 ). For each gene, a model is selected based on posterior probability and estimated effect size (|log_e fold change|) thresholds, which were set to 0.7 for time and 0.5 for genotype. The first two columns depict weighted log₂ fold changes of effect time and genotype. The third column shows the number of genes attributed to this model. The fourth column shows the log expression parameter “mu” in the early and late stages of differentiation and the estimated regression line for WT and GP9 ^{− / −} MKs.

Figure 6

Polytomous model selection among five statistical models capturing different patterns of gene expression in IKZF5 ^WT/Y121F and wild‐type (WT) megakaryocytes (MKs) (as described in Supporting Information S1: Figure 3 ). For each gene, a model is selected based on posterior probability and estimated effect size (|log_e fold change|) thresholds, which were set to 0.7 for time and 0.5 for genotype. The first two columns depict weighted log₂ fold changes of effect time and genotype. The third column shows the number of genes assigned to this model. The fourth column shows the log expression parameter “mu” in the early and late stages of differentiation and the estimated regression line for WT and IKZF5 ^WT/Y121F MKs.

Figure 7

(A) Gene set enrichment analysis for downregulated genes in early stage maturation of GP9^{− / −} megakaryocytes (MKs). Pathways with false discovery rate (FDR) < 0.05 and a minimum of three genes overlapping with those in the pathways (indicated between brackets) are shown. Scatter plot on the right panel shows corresponding log₂ fold changes in early‐stage MKs (mutant vs. wild‐type) compared with log₂ fold changes in platelets (patients vs. controls) of genes involved in “metabolism of RNA” pathways (including ribosomal RNA [rRNA] processing), marked with an asterisk. (B) Gene set enrichment analysis for downregulated genes in late‐stage maturation of GP9 ^{− / −} MKs. Pathways with FDR < 0.05 and a minimum of three genes overlapping with those in the pathways (indicated between brackets) are shown. Scatter plot on the right panel shows corresponding log₂ fold changes in late‐stage MKs (mutant vs. wild‐type) compared with log₂ fold changes in platelets (patients vs. controls) of genes involved in “metabolism of RNA” and actin/tubulin folding pathways, marked with an asterisk. (C) Gene set enrichment analysis for genes downregulated in time and by genotype in IKZF5 ^WT/Y121F MKs. The top 20 pathways with FDR < 0.05 and a minimum of three genes overlapping with those in the pathways (indicated between brackets) are shown. Scatter plot on the right panel shows the corresponding log₂ fold changes in MKs (mutant vs. wild‐type) compared with log₂ fold changes in platelets (patients vs. controls) of genes involved in cell cycle pathways, marked with an asterisk. (D) Histogram of fluorescent intensity measured by flow cytometry of propidium iodide (PI) stained Day 1 differentiated WT and IKZF5 ^WT/Y121F MKs and quantification of peaks corresponding to different levels of ploidy. Dead or apoptotic cells were excluded using 4′,6‐diamidino‐2‐phenylindole (DAPI) staining. (E) The percentage of MKs in G1, S, and G2 cell cycle stages quantified through flow cytometry and PI staining on Day 1 differentiated MKs (pre‐endoreplication stage). This experiment was performed in three independent differentiation replicates. P‐values are shown for pairwise t‐tests between WT and mutant conditions. iMK, immortalized megakaryocyte; tRNA, transfer RNA.