Inflammatory recruitment of healthy hematopoietic stem and progenitor cells in the acute myeloid leukemia niche

Inflammation in the bone marrow (BM) microenvironment is a constitutive component of leukemogenesis in acute myeloid leukemia (AML). Current evidence suggests that both leukemic blasts and stroma secrete proinflammatory factors that actively suppress the function of healthy hematopoietic stem and progenitor cells (HSPCs). HSPCs are also cellular components of the innate immune system, and we reasoned that they may actively propagate the inflammation in the leukemic niche. In two separate congenic models of AML we confirm by evaluation of the BM plasma secretome and HSPC-selective single-cell RNA sequencing (scRNA-Seq) that multipotent progenitors and long-lived stem cells adopt inflammatory gene expression programs, even at low leukemic infiltration of the BM. In particular, we observe interferon gamma (IFN-γ) pathway activation, along with secretion of its chemokine target, CXCL10. We show that AML-derived nanometer-sized extracellular vesicles (EVAML) are sufficient to trigger this inflammatory HSPC response, both in vitro and in vivo. Altogether, our studies indicate that HSPCs are an unrecognized component of the inflammatory adaptation of the BM by leukemic cells. The pro-inflammatory conversion and long-lived presence of HSPCs in the BM along with their regenerative re-expansion during remission may impact clonal selection and disease evolution.


Transcriptional gene expression assays
RNA from FACS-sorted HSPCs were extracted using RNEasy Micro Plus Kit (Qiagen, Hilden, Germany) according to manufacturer's protocol.qRT-PCR: Reverse transcription were performed using SuperScript IV RT (Thermo Fisher Scientific, Waltham, MA, USA) and qRT-PCR using appropriate primers (Supplementary Table 1).RNA-Seq: Libraries were generated using NEX Rapid Dir RNA-Seq 2.0 Kit (Perkin Elmer) and sequenced using NovaSeq (Illumina, San Diego, CA, USA).RNA-seq data were preprocessed using the nf-core/RNAseq pipeline version 3.9.In brief, FastQC version 0.11.9 was used to assess the quality of 100 bp paired-end reads.Trim Galore! version 0.6.7 was then used to trim adapters and low quality bases from raw reads.
Following preprocessing, data were filtered using the R package WGCNA version 1.69 to remove genes with many missing entries or zero variance.Differential expression between pairs of groups was estimated using DESeq2 version 1.38.3.Genes were considered to have significant differential expression if there was an adjusted p-value of less than 0.05.Pathway enrichment analysis was performed using genes with significant differential expression between each group pairing ranked by their log2foldchange using the R package fGSEA version 1.24.0.scRNA-Seq: Gel-in-beads and libraries were generated using Chromium Next GEM Single Cell 3' Reagent Kits v3.1 (10x Genomics) and sequenced using NovaSeq (Illumina).Cell Ranger 7.0.1.was used to map the sequencing data.The raw sequencing data (FASTQ files) was mapped to the mouse reference genome 2020-A (GRCm38, GENCODE vM23/Ensembl 98).Data from all samples were loaded in R (version 4.2.2) and processed using the Seurat package (version 4.3.0-Hao et al. Cell 2021).Cells with low UMI counts were filtered out.The filtering thresholds were determined independently for each library to exclude outliers (A1: 1800 UMIs/cell; A2: 3000 UMIs/cell; A3: 1800 UMIs/cell; B1: 2000 UMIs/cell; B3: 2000 UMIs/cell).Cells with more than 10% mitochondrial genes were also removed.Cells from library B2 were excluded due to the low quality of the library.
The dataset was normalized for sequencing depth per cell and log-transformed using a scaling factor of 10,000.Principal component analysis (PCA) was performed on the 2,000 most variable features.The first 29 PCs were used to generate the uniform manifold approximation and projection (UMAP).Cell cycle was scored using the CellCycleScoring function of Seurat.
Inspection of the UMAP revealed differential distribution assigned to cell cycle.Cell cycle effects were then regressed out using the vars.to.regress parameter in Seurat's ScaleData function.After regression, principal component analysis (PCA) was performed again on the 2,000 most variable features, and the first 27 PCs were used for dimensionality reduction and clustering.Batch effects between libraries were removed using the Harmony package (version 0.1.1).Clustering was performed using the Louvain algorithm, based on a shared nearest neighbour graph, as implemented in the FindClusters function in Seurat with parameters k.param = 5 and resolution = 0.3.The differentially expressed genes per cluster were calculated using the Wilcoxon rank sum test as implemented in the FindAllMarkers function in Seurat, with parameters logfc.threshold= 0.1, return.thresh= 0.01, only.pos= TRUE, min.pct= 0.1.Differentially expressed genes between conditions (AML vs PBS) were calculated using the Wilcoxon rank sum test as implemented in the FindAllMarkers function in Seurat with parameters logfc.threshold= 0.1, return.thresh= 0.01 and min.pct= 0.1.We then performed an over-representation analysis of the significantly differentially expressed genes using the enricher function in the clusterProfiler package (version 4.6.2) with the collections "MH: hallmark gene sets", "M2: curated gene sets", "M5: ontology gene sets" and "M8: cell type signature gene sets" from the Mouse Molecular Signatures Database.Cells were annotated using gene signatures for HSPCs subpopulations (HSCs, MPP-2, MPP-3, and MPP-4) reported by Pietras et al (26).The gene signature for each population was defined by the top 100 genes ranked by the Significance Analysis of Microarrays (SAM) scores.We calculated a module score for each cell based on the average expression levels of every gene signature, subtracted by the aggregated expression of randomly selected control genes.Each cell was assigned to the population with the highest positive score.Cells with only negative scores were annotated as "Unclassified".Differentially expressed genes between conditions (AML vs PBS) per cell type were calculated using the Wilcoxon rank sum test.Over-representation analysis of the differentially expressed genes was performed using the enricher function in the clusterProfiler package (version 4.6.2) with the collections "MH: hallmark gene sets", "M2: curated gene sets", "M5: ontology gene sets" and "M8: cell type signature gene sets" from the Mouse Molecular Signatures Database.Gene expression profiling assay: Extracted RNA were subjected to gene expression analysis using RT2 Profiler PCR Array Mouse Inflammatory Response & Autoimmunity Panel according to manufacture's recommended protocol.miRNA extraction: miRNA were extracted from either cells or EVs using miRNeasy Mini Kit (Qiagen) according manufacture's protocol.For selected miRNA targeted analysis, miRNA were first subjected to RT using TaqMan MicroRNA Reverse Transcription Kit (Thermo), followed by real-time PCR using TaqMan Advanced miRNA Assays (Thermo).For miRNA panel profiling, extracted miRNA were assayed using "Mouse Inflammatory Response & Autoimmunity Focus V2, miRCURY LNA miRNA Focus PCR Panel" (Qiagen) according to manufacture's protocol.

FACS cell sorting or analysis
Cells were sorted and analyzed using either FACSAria Fusion cell sorter or FACS Canto II (Becton Dickinson, Franklin Lakes, NJ, USA).Lineage-negative cells were enriched using EasySep Mouse Hematopoietic Progenitor Cell Isolation Kit (Stem Cell Technologies) according to manufacturer's protocol prior to sorting.Cells were stained with appropriate mouse antibodies (Supplementary Table 2) for LT-HSC (Lin -cKit + Sca-1 + CD150 + CD48 -), ST-HSC (Lin -cKit + Sca-

Extracellular vesicles isolation and characterization
For EV AML generation, EV-depleted FBS were used instead of regular FBS to make complete DMEM and RPMI medium.EV-depleted FBS were generated by subjecting 0.1µm-filtered FBS (GeminiBio) to ultracentrifugation at 100,000 x g for 20 hours at 4°C, where the supernatant fraction served as EV-depleted FBS.AML cells (7E6) were seeded and cultured for 72 hours, followed by EV AML harvests from conditioned medium using differential centrifugation methods involving: 400 x g for 10 mins, 2,000 x g for 20 mins, 10,000 x g for 20 mins, and 100,000 x g for 2 hours.Ultracentrifugation were performed using either Optima L-60 or L-80 instrument with either SW-28, SW-40, or SW-55-Ti rotors (Beckman Coulter, Brea, CA, USA).For BM-EVs and PB-EV: BM plasma was obtained by flushing hind legs (tibia and femur) using HBSS buffer supplemented with Pen/Strep (100U/mL), followed by 400 x g centrifugation for 10 minutes at 4°C.PB plasma was obtained by subjecting PB to 400 x g centrifugation for 10 minutes at 4°C.EVs from both BM and PB plasma were pelleted using differential centrifugation methods described above.Pelleted EVs were resuspended in 0.1µm-filtered sterile PBS supplemented with 1% P/S.EVs were quantified utilizing tunable resistance pulse sensing method using qNano Gold (Izon, Christchurch, New Zealand), with either NP150 or NP200 nanopores.

BM plasma protein assays
Murine BM plasma fluids were collected by flushing 2 hind legs (femur and tibia) with Hanks' Balanced Salt Solution (HBSS) buffer in fixed volumes (10mL).Cells from flushed BM were pelleted using centrifugation, and the supernatant fraction were collected and considered as "BM plasma".BM plasma were directly subjected to protein quantification using Luminex-based multiplex assay mouse cytokine/chemokine 44-Plex Discovery Assay, which was performed by Eve Technologies (Calgary, Alberta, Canada).Human BM plasma samples from AML patients and healthy donors were obtained from Dr. Martin Carroll and Dr. Nicolas Skuli from the Stem Cell and Xenograft Core (Perelman School of Medicine, University of Pennsylvania).Human BM plasma samples were collected from BM aspirates following cell separation using Ficoll.The aqueous plasma layer on the top were subjected to protein quantification using enzyme-linked immunosorbent assay (ELISA) (Thermo Fisher Scientific).

Statistical analysis
Statistical analyses for data generated by qRT-PCR and flow cytometry were performed using GraphPad Prism software (GraphPad Software Inc, La Jolla, CA, USA).Statistical analysis methods used in RNA-Seq and scRNA-Seq analyses is reported in the "Transcriptional gene expression assay" section above.Statistical analysis was calculated by either one-way ANOVA with Tukey test, 2-way ANOVA with Dunnett correction, or two-tailed Student's t-tests where appropriate.P-values of less than 0.05 were considered statistical significance and denoted as * p < 0.05;** p < 0.01; *** p < 0.001; **** p < 0.0001.No tests were used in sample size estimation, and methods of randomization and blinding were utilized.

Graphical illustrations
Illustrations were generated utilizing a licensed application BioRender with rights to use in publications.