Isabelle Allaeys, Guillaume Lemaire, Mickaël Leclercq, et. al. Blood Adv (2024) 8 (11): 2777–2789.

• In a mouse model of SARS-CoV-2 infection, cytokines in the bone marrow are different from those in blood and lungs.
• SARS-CoV-2–mediated inflammation promotes changes to the megakaryocyte transcriptome in the bone marrow.

Megakaryocytes (MKs), integral to platelet production, predominantly reside in the bone marrow (BM) and undergo regulated fragmentation within sinusoid vessels to release platelets into the bloodstream. Inflammatory states and infections influence MK transcription, potentially affecting platelet functionality. Notably, COVID-19 has been associated with altered platelet transcriptomes. In this study, we investigated the hypothesis that severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) infection could affect the transcriptome of BM MKs. Using spatial transcriptomics to discriminate subpopulations of MKs based on proximity to BM sinusoids, we identified ∼19 000 genes in MKs. Machine learning techniques revealed that the transcriptome of healthy murine BM MKs exhibited minimal differences based on proximity to sinusoid vessels. Furthermore, at peak SARS-CoV-2 viremia, when the disease primarily affected the lungs, MKs were not significantly different from those from healthy mice. Conversely, a significant divergence in the MK transcriptome was observed during systemic inflammation, although SARS-CoV-2 RNA was never detected in the BM, and it was no longer detectable in the lungs. Under these conditions, the MK transcriptional landscape was enriched in pathways associated with histone modifications, MK differentiation, NETosis, and autoimmunity, which could not be explained by cell proximity to sinusoid vessels. Notably, the type I interferon signature and calprotectin (S100A8/A9) were not induced in MKs under any condition. However, inflammatory cytokines induced in the blood and lungs of COVID-19 mice were different from those found in the BM, suggesting a discriminating impact of inflammation on this specific subset of cells. Collectively, our data indicate that a new population of BM MKs may emerge through COVID-19–related pathogenesis.

Introduction

Megakaryocytes (MKs) reside primarily in the bone marrow (BM) in which they derive from the myeloid branch of hematopoiesis. They are rare (<0.05% of the BM cellular population) and recognizable due to their large size (up to 100 μm), polyploid nucleus, and expression of surface glycoproteins (eg, IIb and IIIa of the IIb/IIIa complex). Through several rounds of DNA replication without cellular division and an intricate maturation process, MKs generate proplatelets, which protrude into blood sinusoids and fragment into anucleate platelets.^1-3 Thus, MKs are the cellular precursors of platelets, whose RNA content (coding and noncoding) mirrors that of the mother MKs.

MKs help maintain hemostasis through their prominent role in platelet production; however, evidence suggests that they also regulate other functions. They support the hematopoietic stem cell niche by producing cytokines and chemokines (platelet factor 4 [PF4], transforming growth factor β [TGF-β], fibroblast growth factor 1 [FGF1], and insulin-like growth factor [IGF1]).^2,4-6 Their role in immunity is suggested by the expression of immune receptors, including Toll-like receptors that recognize pathogen-associated molecular patterns, receptors capable of direct recognition of pathogens (e.g. C-type lectin-like type II receptor [CLEC-2] and glycoprotein VI [GPVI)), and Fc receptors that can bind antibody-coated antigens, and they can present antigens and promote T-cell proliferation.^7-12 Furthermore, MKs express receptors for cytokines, such as interferon (IFN), interleukin-1 (IL-1), IL6, TNF, and RANTES.^7,8,13,14 

The study of the platelet transcriptome confirms that gene expression is regulated in MKs during inflammation or infection.¹⁵ During sepsis as well as in the murine cecal ligation and puncture model, human and murine platelets, respectively, show hundreds of altered genes, including enhanced expression of ITGA2B (glycoprotein IIb encoding gene).¹⁶ In dengue and influenza virus infection, platelets are enriched in IFN-stimulated genes such as IFTIM3,¹⁷ whereas platelets in patients with systemic lupus erythematosus (SLE) and in murine models of SLE also show enhanced expression of IFN-stimulated genes,^18-20 thus pointing to a common role of IFN in MK functions in both viral infection and autoimmune disease.

MKs are a heterogenous population of cells that vary notably in terms of ploidy and size.² Recent single-cell transcriptomic analyses revealed that distinct subpopulations of MKs exist in the BM, each with potentially different specialized functions. Clustering of MKs on the basis of transcript enrichment identified 4 main subpopulations of MKs related to: (1) “cell cycling,” (2) “platelet generation,” (3) “the hematopoietic stem cell niche,” and the more rare (4) “inflammation.”^2,21,22 Although the “platelet generation”–related MKs are suggested to be found near BM sinusoids, the injection of lipopolysaccharide in mice led to the stimulation and expansion of “inflammation”–related MKs located farther from blood vessels.²¹ Cellular heterogeneity can also result from organ distribution. MKs are found in lungs,^23-25 in which they also generate platelets, and may have a privileged access to allergens, microbiota, and infectious agents. Single-cell transcriptomic analyses shed light on distinct populations of MKs, such as those with an immune phenotype in lungs as well as in the spleen under conditions of sepsis.^24-26 

The severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) is a highly contagious virus that causes COVID-19. Although manifestations and severity of the disease vary between infected individuals, millions of lives were lost due to COVID-19. The virus mainly affects the lungs, which can cause severe respiratory complications and hypoxia. Systemic inflammation can occur even after the virus is largely eliminated from the lungs due to the cytokine storm that sometimes accompanies the disease. Enhanced coagulopathy and thrombosis are observed in patients with COVID-19,^27,28 and studies performed at the onset of the pandemic rapidly reported the increased platelet activation and presence of IFN-stimulated gene signature in platelets.^29-31 Furthermore, single-cell RNA sequencing analyses of cells in the blood circulation of patients with COVID-19 identified a subset of MKs with altered gene transcription, such as enhanced IFN-stimulated genes, major histocompatibility complex II, and calprotectin (S100A8/A9),^32,33 pointing to the effects of SARS-CoV-2 or the dysregulated inflammatory response on MK transcription.

Localized near the sinusoids in which they are ready to release their platelets, MKs may be more likely to sense inflammatory mediators from the blood circulation. We hypothesized that the systemic inflammation prevalent in SARS-CoV-2 infection would affect the transcriptome of the BM MKs, principally those with direct access to blood. To verify this, we used a spatial transcriptomic approach to examine MKs at different stages of pathogenesis.

Methods

Additional details can be found as supplemental Material.

Mice

All animal experiments were conducted under biosafety level 3 confinement and approved by the Université Laval animal committee (approval number: 2021-90). B6.Cg-Tg(K18-hACE2)2Prlmn/J (stock number 3034860) (K18-ACE2 mice) mice were purchased from the Jackson Laboratories (Bar Harbor, ME). Nine-week- to 10-week-old male and female mice were infected with 25 μL of medium containing 250 or 500 tissue culture infectious dose 50 (TCID₅₀) of SARS-CoV-2 or 25 μL of medium for mock-infected mice. Mouse body temperature and weight were recorded every day for 10 days. Mice (19.55 ± 2.95 g at day 0) were euthanized when a weight-loss threshold of 20% was reached. On days 3 and 7, mice (n = 5 per group) were euthanized, and lungs were collected for assessment of viral loads, RNA extraction, and tissue homogenization for cytokine analysis and histological studies. BM cells and BM plasma were obtained by flushing 1 femur in 300 μL of phosphate-buffered saline containing 2 mM EDTA, and cells were spun down for 5 minutes at 300g at room temperature. Supernatant was considered as BM plasma, and cell pellet was used to extract RNA.

Immunofluorescence

Femurs and lungs were harvested, fixed for 24 hours in 10% buffered formalin at 4°C, and femur were decalcified for 14 days in 14% EDTA at 4°C. Lungs or decalcified bones were saturated 24 hours in phosphate-buffered saline containing 30% sucrose at 4°C and cryopreserved in M-1 embedding matrix (Epredia). Staining was performed on 6 μm cryosections using 1 μg/mL of rabbit anti-CD41 (EPR17876, Abcam) or 25 μg/mL of rabbit anti-N (Rockland Immunochemicals Inc) and counterstain with Hoechst 33342 (10 μg/mL; Thermofisher). AF488-goat anti-rabbit (2 μg/mL; Thermofisher) was used as secondary antibody. Tissue autofluorescence was quench for 3 minutes using TrueVIEW Autofluorecence Quenching kit (Vector Laboratories). Sections were mounted using VECTASHIELD Vibrance Antifade Mounting Medium (Vector Laboratories). For quantification, images were acquired at 20× with a slide scanner (Zeiss Axio Scan.Z1) connected to a Colibri 2 camera (Zeiss, Germany). CD41⁺ cells (megakaryocytes) and SARS-CoV2–covering surfaces were quantified using a total surface of at least 1.5 mm², randomly assigned and using Zen 3.3 software and imageJ. MK mean size and mean area were quantified using 70 cells randomly assigned, and frequency distribution was calculated. Representative images were acquired using a Z2 confocal microscope with LSM 800 scanning system (Zeiss, Germany). Images were captured with a 20× objective. The mean size of all MKs located in the region of interest (ROI) and captured for the transcriptomic analyzes was evaluated using Fiji.

Immunohistology

Lungs were fixed 24 hours in 10% neutral buffered formaldehyde solution, embedded in paraffin, and sectioned. The 5-μm thick sections were stained with a modified Carstairs’ method to distinguish platelets and fibrin in histological sections: platelet thrombi (gray), fibrin (orange-red), red blood cells (orange-yellow), and collagen (bright blue).^34,35 A histologic scoring index was developed to evaluate lung damage using the presence of platelet thrombi (0-3) and cell infiltration (0-3).

GeoMx digital spatial profiling for whole transcriptome analysis

Femurs from female K18-hACE2 mice infected with SARS-CoV2 (500 TCID₅₀) for 3 (COVID_3d) or 7 (COVID_7d) days or with conditioned media (Mock) for 7 days were fixed in formalin 24 hours, decalcified 14 days with EDTA, and embedded in paraffin. Ten μm sections were prepared and mounted onto Leica Bond Plus microscope slides. Digital spatial profiling was performed by NanoString Technologies using the GeoMX platform with the Nanostring Mouse Whole Transcriptome Atlas. Immunofluorescence staining was performed using anti-CD41 (EPR17876, Abcam) and anti-endomucin (clone V.7C7, SantaCruz Biotechnology) as morphological markers along with Syto13 DNA nucleic acid staining, to design ROI. ROI were selected based on the fluorescence of interest (CD41) and whether MKs were adjacent to sinusoid vessels (ASV) or nonadjacent to sinusoid vessels (non-ASV). Oligos from the hybridized probes were released by photo-cleaving UV light (385 nm). Sequencing libraries were constructed, and subsequent sequencing and counting were performed on an Illumina NovaSeq6000, and fastq files were processed by the NanoString DND pipeline, resulting in count data for each target probe in each ROI. Limit of quantitation was set to detect gene outliers. The limit was defined by the geometric mean of the negative probes multiplied by the square of geometric standard deviations (SDs) of the negative probes, as generally recommended by the manufacturer. The data for each ROI were adjusted to achieve the same 75th percentile (Q3) of RNA expression signal using a normalization strategy generally used for whole transcriptome data analysis.³⁶ A total of 19 072 genes normalized by third quartile (Q3) were expressed above the limit of quantitation in at least 10% of the ROI. Differentially expressed genes between ROIs and conditions were selected on the basis of statistically significant differences (false discovery rate (FDR); pvalue-adjusted P < .05; and an absolute log₂ fold-change > 1).

Bioinformatic analyses

Gene expression data normalized by third quartile were retrieved from the Nanostring platform. We then identified gene signatures associated with various comparisons, including COVID-19 at 3 days vs control and COVID-19 at 7 days vs control. In these conditions, we also compared cells that were ASV and non-ASV. Signatures were identified with BioDiscML,³⁷ a machine learning sequential minimal optimization algorithm designed to identify predictive features within a given condition (supplemental Table 1). This tool trains thousands of models using many machine learning classifiers (eg, Naïve Bayes, Random Forest, and support vector machines (SVMs). All models are evaluated by crossvalidation procedures, (eg, k-fold, Bootstrapping, and repeated holdout) to prevent overfitting, on which SD was calculated. It also performs feature selection using stepwise learning and top-k features based on information gain ranking. Furthermore, because stepwise feature selection tends to remove all correlated features, for each “short signature,” we retrieved correlated genes using Pearson and Spearman coefficients ≥0.9. Finally, we performed a functional analysis of the signatures using g:Profiler.³⁸ 

Statistics

Statistical analyses were performed using GraphPad Prism software (v10.0.3). Results are expressed as mean ± SD. P values <.05 were considered statistically significant.

Biorender (https://www.biorender.com/) was used to draw schematics and the graphical abstract.

Figure 1.

VIEW LARGEDOWNLOAD PPT

SARS-CoV-2 infection in the K18-hACE2 mouse model. (A) Schematic representation of the mouse experimental protocol. K18-hACE2 mice were inoculated intranasally with SARS-CoV-2 (250 TCID₅₀) or control media. Tissues (plasma [blood], bone marrow, and lung) were collected at 3 and 7 days after infection (DPI). (B) Survival curve of K18-hACE2 mice upon challenge with SARS-CoV-2 (250 TCID₅₀, red color) or with control media (Mock, black color). Mortality (end point requiring euthanasia) expressed as percentage of survival was monitored every day for 10 days (top; n = 17-19 per group). Survival curve statistical analysis was calculated with Log-rank (Mantel-Cox) test; ^∗∗∗∗P < .0001. Weight change (bottom; n = 4-40) was also reported every day for 10 days and expressed as mean (± SD) percentage of weight at day 0; the average weight before infection was 19.55 ± 2.95 g. Unpaired t–test with Welch correction was used to determine significance (^∗P < .05; ^∗∗P < .01; and ^∗∗∗P < .001). (C) Lung sections were stained with a Modified Carstairs method for fibrin and platelets detection: fibrin (bright red), platelets (gray-blue), collagen (bright blue), and red blood cells (orange-yellow). Platelet thrombus (arrow head) and cell infiltration (asterix) are indicated (scale bar, 100 μm). (D) Lung damage was evaluated using a histological score (left) and lung virus titer (right) was quantified at 3 and 7 DPI using lung homogenate on vero cells and expressed in TCID₅₀ per mg of tissue (n = 4-5 per group). Statistics: Mann-Whitney test; ^∗P < .05. (E) Heat map of cytokine and growth factor profiles in blood plasma (blood), bone marrow plasma (bone marrow), and lung during SARS-CoV2 infection (3 DPI, 7 DPI, or Mock). Cytokine levels (pg/mL) were subdivided into 3 groups: 0 to 80 pg/mL (low), 0 to 450 pg/mL (mid), and 0 to 75 000 pg/mL (high). The mean concentration for each cytokine is indicated (n = 4-5). Statistics: 1-way analysis of variance (ANOVA) with Dunnett multiple comparisons (compared with Mock); ^∗P <  .05; ^∗∗P <  .01; ^∗∗∗P < .001; and ^∗∗∗∗P < .0001; n = 4-5 per group.

Figure 2.

VIEW LARGEDOWNLOAD PPT

Megakaryocyte and SARS-CoV-2 detection in tissue during COVID-19. (A) Representative immunofluorescence images of BM MK during SARS-CoV-2 infection. Femurs of mice inoculated with SARS-CoV-2 (3 DPI or 7 DPI) or Mock-infected were stained with anti-CD41 (green) and the nuclei counterstained with Hoechst (blue; scale bar, 20 μm). (B) Megakaryocyte number was quantified in femurs and expressed as mean (± SD) × mm^–2 of tissue. Statistics: 1-way ANOVA with Tukey’s multiple comparisons, ^∗P < .05; n = 3-5 per group. (C) Frequency distribution of MK size (diameter in μm; left) and area (surface in μm² right panel) were analyzed, percentages of frequency are expressed as mean (± SD). Statistics, 2-way ANOVA, Mixed-effects analysis with Tukey’s multiple comparisons test; n = 4-5 per group. (D) Representative immunofluorescence images of BM or lung during SARS-CoV-2 infection. SARS-CoV-2 was detected using anti–SARS-CoV-2 nucleocapsid (green), and nuclei were counterstained with Hoechst (blue) (scale bar, 25 μm). (E) Tissues from mice uninfected (Mock) or infected with SARS-CoV-2 for 3 days or 7 days were analyzed for the presence of SARS-CoV-2. Results are expressed as mean (± SD) of SARS-CoV-2 coverage (% of tissue; top; dotted line indicates the mean value obtained for mock tissues). Viral RNA levels were determined in BM by RT-ddPCR and compared with those in the lungs (bottom). SARS-CoV-2 E gene copies were normalized with Gapdh mRNA copies (n = 4-5 mice per condition). Results are expressed as mean (± SD). Statistics: 2-way ANOVA, mixed-effects analysis with uncorrected Fisher’s LSD, ^∗P < .05; ^∗∗∗∗P < .0001; n = 3-5 per group. BM, bone marrow; DPI, day post-infection; LSD, least significant difference; MK, megakaryocytes; RT-ddPCR, reverse transcription digital droplet polymerase chain reaction.

Figure 3.

VIEW LARGEDOWNLOAD PPT

Spatial transcriptomic analysis of bone marrow megakaryocytes during SARS-CoV-2 infection. Femurs were analysed from K18-hACE2 mice infected with SARS-CoV2 (500 TCID₅₀) for 3 days (COVID_3d) or 7 days (COVID_7d) or with conditioned media (Mock) for 7 days. Ten-μm sections of paraffin-embedded femurs were used to perform spatial transcriptomic analysis. (A) Schematic illustration of NanoString’s GeoMx Digital Spatial Profiler (DSP) workflow. Different steps are indicated: (1) stain: femur sections were hybridized with UV photocleavable probes from the whole mouse genome and fluorescent morphology markers; (2) ROI were selected based on the fluorescence of interest (CD41) and whether MK were adjacent to sinusoid vessels (ASV-MK) or nonadjacent to sinusoid vessels (NASV-MK); (3) ROI were illuminated with UV light that released the barcodes; (4) Each ROI was collected independently using microcapillaries; (5) Sequencing libraries were generated followed by sequencing and counting. A total of 19 072 genes normalized by third quartile (Q3) were expressed above Limit of Quantitation (LOQ) in at least 10% of ROI. (B) Immunofluorescence staining of bone marrow (femur) using morphology markers. Anti-CD41 antibody (turquoise) was used to stain MK, anti-endomucin antibody (purple) to stain sinusoids, and SYTO dye (blue) to stain nuclei; scale bar, 125 μm. Representative ROI (6-8 megakaryocytes per ROI) are illustrated in the right panel. ROI with ASV-MK is yellow, and ROI with NASV-MK is white. (C, D) Volcano Plots comparing COVID_3d (C) or COVID_7d (D) vs Mock. The log2 fold change of each gene is plotted against its statistical significance (−log10 P value). Red dots represent genes significantly upregulated and blue dots, genes significantly downregulated in MK during COVID-19. Thresholds are indicated with dotted lines. (E) UMAP plot of the short (30) gene signature of MK-ROI colored by class (blue circles: Mock [n = 12]; red circles: COVID_7d [n = 12]). (F) Heat map showing differences of the 30-gene signature between Mock and COVID_7d. Each row represents a MK-ROI, and each column represents a gene obtained from the short signature. (G) Relevant enriched biological processes using the 30-gene signature and their P values are presented. The 30-gene signature was obtained using a machine learning approach. BioDiscML was used to classified “Mock” and “COVID_7d” groups, and the best model was variable feature importance (VFI) optimized with false discovery rate (FDR). MCC for the model was 0.997. (H) PF4 mRNA expression was evaluated in bone barrow by RT-ddPCR. PF4 mRNA copies were normalized with Gapdh mRNA copies (n = 4-5 mice per condition). Results are expressed as mean (± SD). Statistics: unpaired t–test with Welch correction; ^∗P < .05 (n = 4-5 per condition). MK, megakaryocytes; mPF4, mouse platelet factor 4; mRNA, messenger RNA; RT-ddPCR, reverse transcription digital droplet polymerase chain reaction.

All animal experiments were conducted under biosafety level 3 confinement and approved by the Université Laval animal committee (2021-90).

Results

SARS-CoV-2 infection in mice stimulates the production of cytokines in BM that differ from those found in the blood and lungs

SARS-CoV-2 requires interaction between viral Spike (S) protein and its counterreceptor, the human angiotensin I converting enzyme 2 (ACE2), to gain cell entry. Given that SARS-CoV-2 protein S does not bind murine ACE2, we used transgenic mice that expressed human ACE2 under the control of the human keratin 18 promoter (K18-hACE2 transgenic mice).³⁹ Mice housed in a NC3 respiratory facility were either intranasally injected with diluent or with SARS-CoV-2 Delta variant (250 TCID₅₀). In these experimental conditions, infected mice showed robust lung viral loads at day 3 and gradually showed signs of severe pathogenesis, as assessed by weight loss and a moribund appearance, until day 7 when the manifestations were the most severe (Figure 1A-D). By day 7, the mice had largely eliminated SARS-CoV-2 infectious particles from the lungs but were affected by systemic inflammation (Figure 1D-E), as expected in this established model.

SARS-CoV-2 infection increases the lung cytokine load, such as IL-6, TNF, and chemokines, in both humans and mice.^40,41 To identify cytokines in the BM and to determine whether they resembled those in the blood or the lungs, we collected BM from the murine femurs, from which we isolated the plasma after eliminating the marrow cells by centrifugation. We screened for the presence of 33 cytokines and chemokines in the BM plasma from uninfected mice and SARS-CoV-2–infected mice at days 3 and 7 and compared these with autologous blood plasma and lung homogenates. As expected, several cytokines were detected in the blood plasma and lungs, both at days 3 and 7 of infection (Figures 1E). Interestingly, while disease developed, cytokines were differentially regulated in the BM compared with blood plasma or lung (Figure 1E; supplemental Figures 1-3). For instance, PF4, which increases in the blood of patients with COVID-19,⁴¹ was also increased in the blood of mice after 3 days of infection, whereas its levels were decreased in the lungs (Figure 1E). Calprotectins (S100A8/A9) were increased in the circulating MKs and platelets in humans with COVID-19^32,33,42 but were reduced in the BM plasma (Figure 1E), potentially due to emergency hematopoiesis and the egress of granulocytes. IFN-γ, which would promote the expression of IFN-stimulated genes in the platelets of patients with COVID-19,³⁰ was increased in the blood and lungs but remained unchanged in the BM (Figure 1E). In contrast, IL-9 and IL-15 levels were increased in the BM, whereas it remained unchanged in the other locations studied (Figure 1E). The striking occurrence of a cytokine gradient in blood relative to the BM compartment prompted us to examine the transcriptome of BM MKs closely and to compare MKs distant from sinusoid vessels with those directly adjacent to blood vessels, because the latter may be more likely to be influenced by circulating inflammatory factors released in response to viral infection.

Use of spatial transcriptomics to examine BM MKs

We first characterized the BM in terms of MK density and occurrence during SARS-CoV-2 infection. We observed an average of 75 MKs per mm^–2 in the uninfected mice, a number which gradually decreased through the course of the disease, falling by 50% at 7 days after infection (Figure 2A,B). Despite a decrease in the number of MKs in the BM, the platelet count remained unchanged during infection at all the examined time points (supplemental Figure 4). Furthermore, the frequency distribution of MK size and area remained unaltered during infection (Figure 2C). Although SARS-CoV-2 was easily detected in the lungs of infected mice, the virus was not detected in the BM at any of the time points examined, using either antibodies directed against the nucleocapsid protein N or polymerase chain reaction to detect the viral genome (Figure 2D,E). These observations suggest that SARS-CoV-2 does not reach the BM at any of the time points studied.

We used spatial transcriptomics to assess the transcriptomes of cells in situ (Figure 3A). MKs in paraformaldehyde-fixed femurs were identified by their characteristic size (mean size in μm ± SD, Mock: 22.91 ± 5.243 [n = 90]; 3 days: 23.07 ± 5.302 [n = 90]; 7 days: 24.72 ± 5.513 [n = 71]), nuclear staining, and CD41 surface expression, whereas the blood vessels were identified by using the endomucin marker. Previous work that used a similar spatial transcriptomic approach to study BM MKs in healthy mice showed that as many as 1000 gene transcripts could be identified by the capture of a single MK.⁴³ To enhance the power of downstream analyses, we determined the transcriptome of 8 to 10 MKs that were directly ASV and compared this with equivalent number of MKs that were non-ASV. We then replicated these analyses in 6 distinct areas each (Figure 3A,B). Using this approach, a total of 19 072 genes were identified from all conditions using this approach.

The transcriptome of BM MKs is altered in COVID-19

We first evaluated the overall impact of infection on MK transcription, irrespective of MK localization within the BM. By day 3 of infection, a total of 32 and 1466 transcripts were significantly upregulated and downregulated, respectively (Figure 3C). The use of machine learning to identify a minimal signature did not generate a reliable signature. Principal component analysis and uniform manifold approximation and projection (UMAP) on the obtained signature failed to exhibit clean separation, and although the Matthew correlation coefficient (MCC = 0.915) was high, the elevated SD (0.095) prevented to confidently identify the presence of a gene signature that would efficiently distinguish the 2 MK populations (supplemental Figure 5). However, pathogenesis had more segregating effects on the MK transcriptome by day 7, because 28 genes were upregulated and 1291 genes were downregulated at this time point (Figure 3D). In addition to the UMAP separation on all genes of the 2 conditions (supplemental Figure 6A), a short predictive signature derived from the machine learning analyses contained 30 gene transcripts (MCC = 0.997; SD = 0.008), among which numerous genes had statistically significant differences (supplemental Figure 6B) that could efficiently distinguish the 2 populations of MKs (Figure 3E). Gene ontology and Kyoto encyclopedia of genes and genomes enrichment analyses determined that pathways relating to biological processes such as histone modifications, MK differentiation, NETosis, COVID-19, and SLE were affected in the MKs in the COVID-19–infected mice (Figure 3F,G). PF4, identified by the signature and generally considered to be specifically expressed in the BM by MKs, was confirmed to decrease after 7 days of infection in an independent cohort of mice treated with the same experimental conditions (Figure 3H), further validating our spatial transcriptomic analyses. Notably, IFN-stimulated genes were not modulated in MKs (supplemental Figure 7A), which may be explained by the absence of IFN increase (Figure 1F) and SARS-CoV-2 (Figure 2D,E) in the BM. Moreover, calprotectin levels, which increase in circulating MKs in humans with COVID-19,³² was not altered in MKs in the BM (supplemental Figure 7B). Thus, in the murine model of lung infection by SARS-CoV-2, there were significant and unique alterations to the transcriptome of BM MKs.

While selecting MKs to process for transcriptome analyses, we selected a balanced number of MKs that were either ASV or non-ASV to the marrow sinusoids. We thus aimed to verify whether MK localization within the BM was a factor that could influence the overall transcription pattern. We identified 528 differently expressed genes in healthy mice (Figure 4A). UMAP plot did not distinguish the 2 populations of MKs (Figure 4B), and machine learning analyses failed to generate a model capable of accurately predicting MK localization (MCC = 0.769; SD = 0.186) in the absence of disease (supplemental Figure 8), suggesting that the spatial localization of MKs relative to sinusoid vessels has no or only a limited impact on MK transcriptomes under healthy conditions.

Figure 4.

VIEW LARGEDOWNLOAD PPT

Spatial transcriptomic analysis of megakaryocytes within the bone marrow. MK-ROI were selected depending on their position: ASV- or NASV-MK. (A) Volcano plots (Mock ASV-MK [n = 6] vs Mock NASV-MK [n = 6]) shows the log2 fold-change of each gene plotted against its statistical significance (−log10 P value). Red dots represent genes significantly upregulated and blue dots, genes significantly downregulated in megakaryocytes depending on their position. Thresholds are indicated with dotted lines. (B) UMAP plot of the short gene signature of MK-ROI colored by class in Mock condition (blue circles: Mock ASV [n = 6]; red circles: Mock NASV [n = 6]). (C) Volcano plot of COVID_7d ASV-MK (n = 6) vs Mock ASV-MK (n = 6) (left) and COVID_7d NASV-MK (n = 6) vs Mock NASV-MK (n = 6) (right) showing the log2 fold-change of each gene plotted against its statistical significance (−log10 P value). Red dots represent genes significantly upregulated and blue dots, genes significantly downregulated in megakaryocytes during COVID-19. Thresholds are indicated with dotted lines. (D) UMAP plot of the short gene signature of MK-ROI colored by class in Mock ASV-MK (blue circles [n = 6]; and COVID_7d ASV-MK red circles [n = 6]). (E) Volcano plot (COVID_7d NASV-MK [n = 6] vs Mock NASV-MK [n = 6]) show the log2 fold-change of each gene plotted against its statistical significance (−log10 P value). Red dots represent genes significantly upregulated and blue dots, genes significantly downregulated in megakaryocytes depending on their position. Thresholds are indicated with dotted lines. (F) UMAP plot of the short gene signature of MK-ROI colored by class in Mock NASV-MK (blue circles [n = 6]) and COVID_7d NASV-MK (red circles [n = 6]). (G) Heat map showing differences of the 15-gene signature between ASV-MK in Mock and COVID_7d. Each row represents a MK-ROI and each column represents a gene obtained from the short signature. The 15-gene signature was obtained using a machine learning approach. BioDiscML was used to classify “Mock ASV-MK” and “COVID_7d ASV-MK.” The best model was the Functions Logistic optimized with BER (Balanced Error Rate), and MCC for the model was 0.922 ± 0.113. (H) Heat map showing differences of the 30-gene signature between NASV-MK in Mock and COVID_7d. Each row represents a MK-ROI, and each column represents a gene obtained from the short signature. The 15-gene signature was obtained using a machine learning approach. BioDiscML was used to classified “Mock NASV-MK” and “COVID_7d NASV-MK.” The best model was the complement Naïve Bayes with FDR; MCC for the model was 0.922 ± 0.114.

Given the effects of SARS-CoV-2 on MKs 7 days after infection, we examined at this time point whether the pathogenesis would preferentially alter MKs in a specific location relatively to BM sinusoids. Infection led to the respective upregulation and downregulation of 83 and 1084 genes in ASV MK (Figure 4C), whereas 40 and 747 were upregulated and downregulated in non-ASV MKs (Figure 4E). Non-ASV and ASV MKs from healthy and day 7 COVID-19 mice were distinguished by machine learning with a 15-gene signature that sufficed to specifically identify ASV MKs (Figure 4D and F-H), although a high SD was measured due to the presence of outliers (MCC = 0.922; SD = 0.113). A longer 30-gene signature was needed for the non-ASV MKs (MCC = 0.922; SD = 0.114), which is attributable to model overfitting and points to the absence of significant difference between cells in this location. Gene ontology and Kyoto encyclopedia of genes and genomes enrichment analyses using both signatures highlighted pathways such as histone modifications, MK differentiation, NETosis, and SLE, both in ASV and non-ASV MKs (supplemental Figure 9). Because these results are reminiscent of those found in comparisons that did not take into account position relative to sinusoid vessels, it suggests that most of the enriched pathways occur in MKs, irrespective of their location.

Discussion

Patients with COVID-19 have higher risk of thrombosis,^27,28 and postmortem investigations as well as intravital imaging in SARS-CoV-2–infected K18-hACE2 transgenic mice revealed platelet thrombi in the lung and brain vasculature.^40,44 Given the alterations to the platelet transcriptome in patients with COVID-19,^30,42 which may give rise to hyperactivated platelets,^29-31,45 we aimed to determine whether MKs in the BM were capable of sensing the disease and thereby modifying their transcriptome. Using a murine model of SARS-CoV-2 infection, we identified a set of cytokines and chemokines that are differentially expressed in the BM vs the blood and lungs at different stages of the disease. We show that the MK transcriptomic landscape is profoundly affected at the time of systemic inflammation, when the virus is essentially cleared from the lungs and is absent in the BM. Moreover, we suggest that the spatial localization of the cells within the marrow is not a variable that would affect MK transcriptomic alterations in this disease.

Respiratory tract infections are highly contagious and are responsible for millions of deaths worldwide each year, often following seasonal cycles. During influenza or Sendai virus infection, type I IFN produced in the lungs can affect cells in the BM,⁴⁶ thereby driving antiviral responses in emerging immune cells despite the absence of viral particles in the BM. IL-1 and IL-6, produced during influenza infection can also trigger emergency megakaryopoiesis, driving MK production directly from CD41⁺ (MK-committed) hematopoietic stem cells.⁴⁷ Megakaryopoiesis, however, was unlikely implicated in COVID-19 mice, because the number of MKs was reduced by the infection. SARS-CoV-2 primarily infects the host lung epithelium to cause respiratory tract inflammation, although its persistence and distribution in the human body, such as in kidneys, testis, and brain, varies greatly between individuals.^48-51 We show that MKs localized in the BM, more specifically in the femur, respond to inflammation through changes to their transcriptome during the course of lung infection. The transcriptional changes are unlikely to be due to SARS-CoV-2 itself, because the virus was mostly eliminated at this time and was undetectable in the BM. IFN was not induced in the BM, and the MKs in the BM examined at all stages of the disease never overexpressed IFN-stimulated genes (supplemental Figure 7A), perhaps due to the low IFN-inducing potential of SARS-CoV-2.⁵² It remains to be determined what molecule(s) stimulates MKs in the BM, but our survey of cytokines in both the BM and blood suggest that cytokines such as RANTES, known to induce proplatelet formation by MKs, or IL-9 and IL-15 might be implicated.⁵³ Damage-associated molecular patterns generated in the pathogenesis⁵⁴ or the reduced levels of S100A8/A9, a trigger of toll-like receptor 4, observed in the BM, may also play a role.⁵⁵ Moreover, there may be a role for extracellular vesicles, reportedly increased in COVID-19 and capable of reaching MKs in the BM,^56,57 or lipid mediators of inflammation, also increased in COVID-19.⁴¹ SARS-CoV-2 stimulates tissue factor expression from lung epithelial cells and thrombin activity,⁵⁸ which too may activate MKs⁵⁹ and might affect the MK transcriptome.

Given the occurrence of both intravascular and extravascular MKs in lungs,²⁵ we cannot rule out that lung MKs might preferentially interact with SARS-CoV-2. Such interactions may involve ACE2,⁶⁰ although its expression in human MKs is not definitive,⁶¹ or other molecules such as CD147.^62,63 The transcript encoding murine ACE2, which does not interact with SARS-CoV-2, was detected in the BM MKs in our analyses (supplemental Figure 10). In humans, MKs in the blood circulation of patients with severe COVID-19 overexpress IFN-stimulated genes and S100A8/9 and contain SARS-CoV-2 RNA and particles.^30,32,64 In our investigation, we did not observe these features in BM MKs. Although this might be due to intrinsic differences in the pathogenesis between humans and murine models, it may also suggest that a different pool of MKs from elsewhere than the BM, for example, the lungs, makes contact with SARS-CoV-2. These would be the MKs that generate the platelets containing SARS-CoV-2 (particles and RNA) reported in these patients.^29,30,65,66 MKs may also encounter viral particles after they egress the BM sinusoids during their transit toward the lung vasculature. Although this is unlikely, given the low number of viral particles in blood and the expected short duration of this event, studies suggest that MKs may potentially repeatedly circulate through the lung vasculature, which may contribute to longer exposure to viral and inflammatory mediators in lungs.⁶⁷ In sepsis, there is expansion of the MK population in the spleen.²⁶ Whether this occurs in the murine model is currently unknown.

Valuable information on the occurrence of potential subpopulations of MKs in BM was provided by recent single-cell RNA sequencing approaches.^21,68 A major challenge in the study of MKs is maintaining the integrity of these rare cells, and this is even more true for those imbricated in sinusoids because they are less likely to be captured by the harvest of the marrow or to survive the process of single-cell RNA sequencing analysis. This is an issue, because the single-cell RNA sequencing approach depends on cell integrity and excludes apoptotic cells. Moreover, most common, commercial single-cell sequencing approaches (eg, 10XGenomics, BD Rhapsody) use a 10 μm diameter nozzle, suggesting that the larger MKs may be lost through the process and may be underrepresented. To overcome this and for safety reasons given the use of virus-infected materials, we examined cells in femurs fixed in paraformaldehyde using Nanostring GeoMx. In addition to conserving the hard-to-get MKs imbricated in sinusoids, a notable advantage of this sort of analysis is the richness in terms of RNA content of polyploid cells, which would theoretically allow for the study of very few MKs while obtaining sufficient number of RNA molecules for downstream bioinformatic analyses. Little transcriptional changes were found in MKs at day 3 of infection, despite the large viral load at this time and the noticeable changes to cytokine composition in the BM. In contrast, by day 7, we observed numerous changes, attributed by pathway enrichment analyses to histone modifications, MK differentiation, NETosis, and autoimmunity. These characteristics may be congruent with changes to subpopulations of MKs, such as cycling and platelet-producing MKs. While separating MKs into the 4 subclasses on the basis of their transcripts signature identified by Sun et al²¹ (cell cycling, platelet generation, hematopoietic stem cell niche, and inflammatory), we find that, with the exception of pf4 and itgb2, none of the transcript changes comprised a particular subclass (supplemental Figure 11). The data suggest that modifications of histones may favor DNA replication, whereas MK differentiation may result from emergency megakaryopoiesis. As for NETosis, it is reportedly enhanced in COVID-19, whereas SLE shares common features with COVID-19 such as dysregulated immune cells and autoantibody production. Whether emperipolesis, the process by which neutrophils migrate into MKs,^69,70 is enhanced in COVID-19 and explains these changes to transcriptome is also a hypothesis to consider. More mechanistic investigations are needed to confirm the role of these MKs in inflammation and infection.

There are limitations to our study. Human ACE2 expression in K18-hACE2 transgenic mice may not fully replicate its natural distribution in humans. For instance, MKs may not interact directly with SARS-CoV-2 in these mice, although studies of human platelets suggest that SARS-CoV-2 is internalized despite ACE2 blockade.⁶⁵ Tissue immunofluorescence (2-dimensional imaging) has previously been used to define 4 main subtypes of MKs, and it highlights that both platelet-producing and hematopoietic stem cell–niche MKs were localized in the proximity to sinusoid vessels, while immune subpopulations were more distant.²¹ Similarly, we used 10-μm thick bone sections in our selection of MKs relative to sinusoid proximity. However, there might be vessels that were not efficiently visualized and might have been missed due to the 2-dimensional imaging approach, as demonstrated elsewhere.⁷¹ We did not identify any impact of sinusoid vessel proximity on MK transcriptome in healthy conditions. Even when the mice showed severe signs of disease, direct access of MKs to blood circulation was not a significant factor, even if a longer gene signature was necessary to segregate non-ASV MKs from infected and healthy mice. Although these findings need to be replicated, they are consistent with a very limited intervascular space in the BM and the suggestion that MKs may thus be vessel-biased in this tissue.⁷¹ 

The appreciation of MK function in immunity, their diversity in term of subpopulations, and their localization in various organs and tissues motivated broad research. We found that during SARS-CoV-2 infection, particularly during systemic inflammation, the transcriptomic profile of BM MKs may change. This suggests that the localization of MKs in the BM does not hinder their response to inflammatory stimuli in COVID-19, even in the absence of persistent infection. Such responsiveness may have implications for their platelet progeny and disease manifestations, as observed in conditions such as long COVID