# Dual role of plasmablasts as immune regulators or amplifiers in COVID-19.

Emmanuel Ledoult, Thomas Guerrier, Sylvain Dubucquoi, Martin Figeac, Celine Villenet, Blanche Daunou, Helene Behal, M. R. Pokeerbux, Thomas Machet, Vincent Koether, et al.

► To cite this version:

Emmanuel Ledoult, Thomas Guerrier, Sylvain Dubucquoi, Martin Figeac, Celine Villenet, et al. Dual role of plasmablasts as immune regulators or amplifiers in COVID-19. Clinical Immunology, 2025, Clinical Immunology, 281, pp.110609. 10.1016/j.clim.2025.110609. hal-05371308

HAL Id: hal-05371308

https://lilloa.hal.science/hal-05371308v1

Submitted on 18 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
remodeling of the B-cell compartment. Severe COVID-19 has been associated with a reduction in classical (CD19⁺CD27⁺) memory B cells, which are essential for long-term immunity, as well as an expansion of double-negative (CD19⁺CD27⁻IgD⁻, DN) B cells. These DN B cells, particularly the CD11c⁺CD21⁻ subset (DN2), have been implicated in inflammatory responses and autoimmune-like conditions [11–13]. Notably in the context of primary infection, DN2 B cells also represent an extrafollicular source PB from naive B cells, which may suggest that upstream dysregulation may shape downstream PB function.

Here, we aimed to dissect the functional diversity of PB in COVID-19 using a combination of longitudinal immunophenotyping and cytokine assays, mathematical modeling, and transcriptomic profiling in two complementary patient cohorts. We identify two distinct PB trajectories, a rapidly resolving versus a persistent and amplified expansion, associated with divergent clinical outcomes and transcriptional programs. These findings suggest that PB may play dual roles in COVID-19, acting either as regulators of antiviral immunity or amplifiers of inflammation in severe disease.

## 2. Results

### 2.1. Patient demographics and characteristics

Cohort 1 prospectively included 50 patients with COVID-19 who were enrolled during the second wave of COVID-19 in France. The median (IQR) age was 62.0 (53.2; 71.8) years. Patients were predominantly male (74 %). The median body mass index (BMI) was 30.2 (25.4; 33.1) Kg/m². The median time from symptom onset to inclusion was 7.0 (5.0; 9.0) days. The maximum severity was critical in 36 cases (76 %), and 10 patients (20 %) died. Dexamethasone was used in 41 cases (82 %). Visits and blood sampling (soluble serum factors assay and immunophenotyping) were scheduled at D1, D7, D14, and D30 post-diagnosis. One hundred and eleven extended B-cell phenotyping were performed, respectively in 39 patients at D1, in 31 patients at D7, 14 patients at D14 and 27 patients at D30. Samples from healthy donors were analyzed as controls (n = 18) (Table 1).

### 2.2. Circulating B-cell subsets in COVID-19 in cohort 1

We prospectively applied a multiparameter flow cytometry approach to identify and longitudinally quantify peripheral blood B-cell subsets after COVID-19 diagnosis in cohort 1 at D1 (n = 39), D7 (n = 31), D14 (n = 14) and D30 (n = 27).

At baseline (D1), absolute B-cell (CD19⁺) counts from patients and controls were not significantly different (Fig. 1a). Compared to controls (n = 18), we observed in COVID-19 patients that: (i) PB were elevated (both in absolute and relative values) (p < 0.001); (ii) memory B cells were decreased (both in absolute and relative values) (p < 0.001); (iii) switched memory B cells were decreased (p < 0.001); and (iv) marginal zone-like (MZ) B-cell proportions were decreased (p = 0.01). The full distribution of B-cell subsets at D1 is shown in Supplementary Results-Table 1 and in Fig. 1b-j. We investigated whether baseline (D1) circulating B-cell subsets correlated with maximal severity expressed on a 7-point ordinal clinical scale. Among the B-cell subsets assessed, we only identified significant positive correlations between the maximal disease severity and the baseline (D1) subset size of PB (r = 0.35, p = 0.02) and of DN2 B cells (r = 0.33, p = 0.03) (Fig. 2).

After observing global dynamic changes in B-cell subsets in COVID-19 (Fig. 1), we analyzed the B-cell subset distribution over time stratified by the maximal severity (non-severe, n = 4; severe, n = 8; and critical, n = 27). At D1, the median PB proportion was increased in all groups (14 % [10; 21], 6 % [4; 14], 13 % [8; 19] in the non-severe, severe, and critical groups, respectively) compared to controls (2 % [1; 3]). At D7, there was an overall reduction in PB proportion compared to D1. The PB expansion persisted longer and at higher levels in the critical group (8 % [4; 13]) than in the severe (4 % [3; 5]) and non-severe (1 % [1; 1])

<table><caption>**Table 1**

Patient demographics and characteristics.</caption><thead><tr><th></th><th>Cohort 1 (n = 50)</th><th>Cohort 2 (n = 282)</th></tr></thead><tbody><tr><td><strong>Demographics</strong></td><td></td><td></td></tr><tr><td>Male; n (%)</td><td>37 (74)</td><td>185 (66)</td></tr><tr><td>Age; median (IQR); y</td><td>62.0 (53.2; 71.8)</td><td>60.6 (50.8; 70.5)</td></tr><tr><td><strong>Medical background</strong></td><td></td><td></td></tr><tr><td>BMI; median (IQR); kg/m²</td><td>30.2 (25.4; 33.1)</td><td>28.4 (25.1; 32.7)</td></tr><tr><td>Arterial hypertension; n (%)</td><td>28 (56)</td><td>127 (45)</td></tr><tr><td>Diabetes; n (%)</td><td>13 (26)</td><td>66 (23)</td></tr><tr><td>Heart disease; n (%)</td><td>8 (16)</td><td>34 (12)</td></tr><tr><td>Chronic lung disease; n (%)</td><td>8 (16)</td><td>54 (19)</td></tr><tr><td>Cancer; n (%)</td><td>0 (0)</td><td>30 (10)</td></tr><tr><td>Severe liver disease; n (%)</td><td>0 (0)</td><td>13 (5)</td></tr><tr><td><strong>COVID-19 characteristics</strong></td><td></td><td></td></tr><tr><td>Days since disease onset; median (IQR)</td><td>7.0 (5.0; 9.0)</td><td>7.0 (4.0; 11.0)</td></tr><tr><td>Severity at first record<sup>a</sup>; n (%)</td><td></td><td></td></tr><tr><td>non-severe</td><td>12 (24)</td><td>98 (35)</td></tr><tr><td>severe</td><td>11 (22)</td><td>106 (38)</td></tr><tr><td>critical</td><td>27 (54)</td><td>78 (28)</td></tr><tr><td>0-7<sup>b</sup>; median (IQR)</td><td>5.0 (4.0; 5.0)</td><td>4.0 (3.0; 5.0)</td></tr><tr><td><strong>Biological results at inclusion</strong></td><td></td><td></td></tr><tr><td>Ferritin; ng/mL; median (IQR)</td><td>956 (641; 1758)</td><td>904 (358; 1640)</td></tr><tr><td>Fibrinogen; g/L; median (IQR)</td><td>6.6 (5.9; 7.4)</td><td>6.3 (4.9; 7.4)</td></tr><tr><td>Plasmablasts, %; median (min; max)</td><td>12.3 (5.4; 19.7)</td><td>9.5 (0; 58.3)</td></tr><tr><td>Plasmablasts, c/µL; median (min; max)</td><td>15.9 (5.8; 28.5)</td><td>14.9 (0; 496.3)</td></tr><tr><td><strong>Follow-up</strong></td><td></td><td></td></tr><tr><td>Maximum severity<sup>a</sup>; n (%)</td><td></td><td></td></tr><tr><td>non-severe</td><td>6 (12)</td><td>102 (36)</td></tr><tr><td>severe</td><td>8 (16)</td><td>59 (21)</td></tr><tr><td>critical</td><td>36 (72)</td><td>121 (43)</td></tr><tr><td>0-7<sup>b</sup>; median (IQR)</td><td>5.0 (4.0; 6.0)</td><td>4.0 (3.0; 6.0)</td></tr><tr><td>Death; n (%)</td><td>10 (20)</td><td>22 (8)</td></tr><tr><td>Time to maximum severity; days; median (IQR)</td><td>9.0 (6.0; 14.0)</td><td>10.0 (6.0; 14.0)</td></tr><tr><td><strong>Treatment</strong></td><td></td><td></td></tr><tr><td>Dexamethasone; n (%)</td><td>41 (82)</td><td>41 (15)</td></tr><tr><td>Tocilizumab; n (%)</td><td>2 (4)</td><td>3 (1)</td></tr><tr><td>Antiviral therapy; n (%)</td><td>15 (30)</td><td>41 (15)</td></tr></tbody></table>

groups. At D30, the PB proportion was similar between the 3 groups of severity (Fig. 1b, f-g and Supplementary Results-Tables 2 to 5).

These results demonstrate an expansion of PB from the first days of COVID-19, with a subsequent normalization over time (among the survivors), occurring rapidly in the non-severe and severe patients and more slowly in critical patients. When focusing on the two main sources of PB—the memory and DN2 B cells, no difference was observed according to the severity groups; however, there was a gradual normalization of this subset over time (Fig. 1c, h and Supplementary Results-Tables 2 to 5).

### 2.3. Soluble serum factors in COVID-19

To characterize the dynamic changes of soluble serum factors in response to COVID-19, we prospectively measured 31 soluble serum factors in cohort 1 at D1 (n = 39), D7 (n = 31), D14 (n = 14) and D30 (n = 27). Samples from healthy donors served as controls (n = 18).

#### 2.3.1. Longitudinal analysis of soluble serum factors in COVID-19 patients compared to controls

At baseline (D1), BAFF (p < 0.0001) and APRIL (p < 0.01) blood levels were significantly increased in patients with COVID-19 and remained higher than in controls at D30 (BAFF, p < 0.05; APRIL, p < 0.01). Interestingly, we observed in patients with COVID-19 a rapid transient fall in BAFF blood levels between D1 and D7 before rising again at D14, which may be related to the short course of steroids [14]

<footnote>a Severity according to WHO classification (e.g., non-severe/severe/critical).</footnote>

<footnote>b 0–7 scale used for severity scoring (details in Methods).</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Longitudinal analysis of B-cell subset in patients from cohort 1. Longitudinal quantification of peripheral blood B-cell subset in cohort 1 at D1 (n = 39), D7 (n = 31), D14 (n = 14), D30 (n = 27) and in healthy controls ((HC, n = 18). Data are expressed by median. Comparisons of classes were performed using Kruskal-Wallis test. In the case of significant results, pairwise comparisons were performed, and a Bonferroni correction was applied. \*: p < 0.05, \*: p < 0.01, \*: p < 0.001, \*\*\*p < 0.0001, \*\*\*\*: p < 0.00001. (a) B cells (absolute number (μL). (b/g) plasmablasts (absolute number (μL) and proportion). (c/h) memory B cells (absolute number (μL) and proportion). (d/i) switched memory B cells (absolute number (μL) and proportion). (e/j) Marginal Zone-like B cells (absolute number (μL) and proportion). (f) IgM+ among plasmablasts (proportion).</figcaption></figure>

recommended for COVID-19 (82 % of patients in cohort 1 received 10 days of dexamethasone) (Fig. 3a). sCD27 levels were significantly higher in patients with COVID-19 from D1 (p < 0.0001) to D30 (p < 0.001) than in controls. We did not observe any significant changes in IL-10, IL-21, IL-4, IL-5, and TSLP (Supplemental Results-Fig. 1a). CXCL13 blood levels remained elevated throughout the follow-up in patients with COVID-19 (p < 0.05). We noted an early, marked and significant elevation of other markers related to neutrophil activation and inflammation—such as IL-8, CXCL1, IL-6, MPO, PR3, lactoferrin, and MMP-8—in patients with COVID-19, followed by a slow decrease over time (Fig. 3b, d and Supplemental Results-Fig. 1e). Regarding the IFN-γ-induced chemokines—CCL2, CXCL10 and CXCL11—, the blood levels were higher from D1 to D30 in COVID-19 patients compared to controls (Figs. 3b). We found no significant differences about the serum levels of CCL11, CCL24, CCL26, IL-3 and IL-13 (Supplemental Results-Figs. 1be).

### 2.3.2. Correlations between serum levels of soluble factors and disease severity

We assessed correlations between the serum levels of soluble factors and the maximum severity. At D1, neither APRIL ($r = 0.07, p = 0.7$) nor BAFF ($r = -0.28, p = 0.12$) blood levels correlated, whereas at D7, APRIL blood levels negatively correlated ($r = -0.49, p < 0.001$) and BAFF blood levels correlated positively ($r = 0.74, p < 0.001$). The maximal severity correlated positively with baseline (D1) levels of inflammatory and neutrophil-activation factors such as MMP-8 ($r = 0.53, p = 0.009$), IL-8 ($r = 0.42, p = 0.015$), and CCL24 ($r = 0.49, p = 0.004$),

and negatively with CCL27 ($r = -0.39, p = 0.02$) (Fig. 3e, Supplemental Results-Fig. 2). At D7, we observed strong positive correlations between the maximal severity and CCL2 ($r = 0.73, p < 0.001$), CXCL2 ($r = 0.47, p = 0.04$), and CXCL10 ($r = 0.76, p = 0.004$) blood levels. At D14, higher levels of IL-6 ($r = 0.73, p = 0.05$) and CXCL13 ($r = 0.81, p = 0.02$) correlated positively with the maximum severity (Fig. 3e, Supplemental Results-Fig. 2b-d).

These results indicate a strong correlation between the maximum severity and the blood levels of soluble factors involved in the inflammatory storm, particularly between D7 and D14 rather than at baseline (D1). A strong correlation was also observed between the maximal severity and BAFF blood levels at D7.

### 2.3.3. Correlation between soluble serum factors and circulating B-cell subsets

We next analyzed the correlations between the serum levels of soluble factors and the proportion of PB, DN2, and classical memory B cells. At D1, we found strong positive correlations between the PB proportion and sCD27 ($r = 0.53, p = 0.03$), IL-13 ($r = 0.45, p = 0.02$), and CCL17 ($r = 0.46, p = 0.022$). At D7, we observed a positive correlation of BAFF ($r = 0.59, p = 0.002$) and IL-10 blood levels ($r = 0.43, p = 0.04$) with DN2 B-cell proportion, but not with PB proportion. PB and DN2 B-cell proportions showed positive correlations with IL-8 ($r = 0.45, p = 0.033; r = 0.42, p = 0.001$, respectively). At D30, PB proportion correlated positively with sCD27 ($r = 0.77, p = 0.007$) and CXCL10 ($r = 0.90, p = 0.001$). Similar results were observed for the DN2 B-cell subset (Supplemental Results-Fig. 3).
<figure><img src="image_3.png" /></figure>

<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Correlations between circulating B cell subsets at D1 with the maximal severity related to COVID-19 in cohort 1. (a) Spearman correlation. (b) p-value.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>**Fig. 3.** Longitudinal analysis of soluble serum factors in COVID-19 and correlations with the maximum severity in cohort 1. Longitudinal analysis of soluble serum factors assessed in at D1 (n = 33), D7 (n = 20), D14 (n = 9) and D30 (n = 14), in patients with COVID-19 and in healthy controls (n = 18). Data are expressed by median. Comparisons of classes were performed using Kruskal-Wallis test. In the case of significant results, pairwise comparisons were performed, and a Bonferroni correction was applied. \*: *p* < 0.05, \*: \*: *p* < 0.01, \*: \*: \*: *p* < 0.001, \*: \*: \*: *p* < 0.0001. Soluble factors were grouped according to their dominant immunological roles: (a) B-cell support / plasmablast differentiation; (b) Inflammatory chemokines and cytokines; (c) Mixed / regulatory cytokine; (d) Neutrophil degranulation / tissue damage. (e) Spearman correlation between the maximal severity related to COVID-19 expressed on a 7-point ordinal clinical scale (see Methods section) and the blood levels of soluble factors at D1 (n = 33), D7 (n = 20), D14 (n = 9) and D30 (n = 14). Soluble factors were grouped according to their dominant immunological roles: B-cell support / plasmablast differentiation; Inflammatory chemokines and cytokines; Mixed / regulatory cytokine; Myeloid and hematopoietic growth factors; Neutrophil degranulation / tissue damage. Several markers exhibit pleiotropic roles and are annotated accordingly (see *Supplemental Methods, Table 4*).</figcaption></figure>
Overall, these results suggest a positive correlation between PB
proportion and DN2 B cells and the levels of soluble serum factors
involved in the inflammatory storm in COVID-19, particularly from D7
onward.

## 2.4. Modeling of plasmablast trajectories

To further decipher the role of PB in the severity of COVID-19, we established a cohort before the widespread use of steroids in COVID-19. A total of 282 patients had documented COVID-19 infection and at least one lymphocyte phenotyping (cohort 2) (Table 1 and Supplemental Results-Table 6). At that time, few patients were treated by steroids (15 %) or tocilizumab (1 %) allowing us to explore the longitudinal behavior and functional profiles of PB in COVID-19. We identified a high mean (±SD) PB proportion (13.7 % ±12.4) with wide variability between patients (range min to max: 0–58.3 %) (Fig. 4a). Baseline PB proportions (both in % and cells/µL) correlated positively with the maximal disease severity (r = 0.53, p < 0.0001) and with numerous biological variables relative to the inflammation storm including IL-6 (r = 0.53, p < 0.0001), fibrinogen (r = 0.54, p < 0.0001), C-reactive protein (r = 0.61, p < 0.0001), and ferritin (r = 0.50, p < 0.0001) (Supplemental Results-Table 7).

### 2.4.1. Modeling of the plasmablast trajectories in COVID-19

To determine whether longitudinal changes in PB (%) were associ-
ated with disease severity in COVID-19, we selected patients who had at
least two records of PB proportion over time (Cohort 2b, n = 166)
(Fig. 4a, Supplemental Results-Table 6). The median time of the first
immunophenotyping was 12.0 days (IQR: 8.0; 18.0) since disease onset.
Models with one to four latent classes were sequentially performed
(Supplemental Results-Fig. 4). The 2-class model represented the best fit
to the data (Supplemental Results-Tables 8–9, Supplemental Results-
Figs. 4b-c). The results were similar when including only patients with at
least 3 records of PB overtime (Supplemental Results-Table 10, Sup-
plemental Results-Fig. 6) and when excluding patients who received
steroids (data not shown).

Class 1 showed a gradual decrease from a mean estimated baseline PB proportion of 16.3 % [95 % CI: 14.4; 18.1] to 5.0 % [95 % CI: 3.5; 6.5] at 30 days. One hundred and fifty patients were assigned to this class. Class 2 was characterized by a quadratic trajectory of PB proportions, with a median estimated baseline of 37.3 % [95 % CI: 14.4; 18.1], followed by an incremental phase to 40.9 % [95 % CI: 35.3; 46.2] at 8 days of follow-up, then a subsequent resolutive phase (mean estimated PB proportion at 30 days: 15.6 % [95 % CI: 9.7; 21.5]). Sixteen patients were assigned to this class. Compared to class 1, class 2 included a higher proportion of male patients (94 % vs. 69 %) and older patients (median age of 67.0 (59.5; 70.8) vs. 61.6 (51.1; 69.7) years old) (Table 2). The maximum disease severity was higher in class 2 (5.5 (4.0; 7.0)) than in class 1 (4.0 (4.0; 5.0)). Patients in class 2 also had a higher proportion of death at 30-day from ICU admission (31 % vs. 5 %, p < 0.001; Fig. 4d), even after adjusting for sex and age (p = 0.008). Taken together, the massive and prolonged expansion of PB (%) was associated with higher mortality related to COVID-19.

## 2.5. Transcriptomic analysis of plasmablasts in COVID-19

To explore the hypothesis that PB exhibit different functional profiles
in COVID-19 according to disease severity, we selected 41 samples from
23 patients in cohort 2b according to maximum severity (either non-
severe, n = 24; or severe [severe or critical], n = 17) and time since
onset of disease ("early" [n = 21] or "late" [n = 20] time points, less than
14 days or more than 20 days from symptom onset, respectively). In the
non-severe group, 9 males and 1 female were included, whereas in the
severe/critical group, 13 males and no females were included. The
median age was 70.0 (IQR: 56.5; 78.4) and 60.5 (53.2; 67.5) years-old,
respectively in the non-severe group and in the severe/critical group.

After isolating PB from other B cells using FACS, we investigated the transcript expressions by bulk RNA sequencing. Ten samples were eliminated from the analysis due to poor RNA quality. A total of 7096 genes were identified from the 72 samples (Fig. 5a). Principal component analyses distinguished PB from other B cells but did not clearly distinguish neither PB nor other B cells according to disease mortality or time of sampling (T<sub>early</sub> or T<sub>late</sub>) (Fig. 5b-e). We then focused on PB according to both maximal disease severity (severe or non-severe) and time of sampling (T<sub>early</sub> or T<sub>late</sub>). We identified 733 differentially expressed genes (DEG) at T<sub>early</sub>, between severe vs non-severe groups, respectively 154 DEGs were overexpressed in severe group (compared to non-severe group, q < 0.1 and absolute fold change (FC) > 1.5) and 579 DEGs were overexpressed in non-severe group (compared to severe group) (Figs. 5b and 6). At T<sub>late</sub>, 45 DEGs were identified, with 37 DEGs overexpressed in the severe group and 8 DEGs overexpressed in the non-severe group (Fig. 6, Supplementary Results-Fig. 5).

### 2.5.1. PB transcriptomic signatures at T<sub>early</sub> time point according to disease severity

Enrichment analysis of DEGs overexpressed in the severe group (compared to the non-severe group) revealed enrichment for pathways related to oxidative phosphorylation (e.g. hsa00190: oxidative phosphorylation), aerobic respiration (e.g. R-HSA-1428517: aerobic respiration and respiratory electron transport, WP4922: Mitochondrial complex IV assembly), metabolism of RNA (R-HSA-8953854: Metabolism of RNA, R-HSA-72163: mRNA Splicing), nucleoside triphosphate metabolic process (e.g. GO:0009141: nucleoside triphosphate metabolic process, GO:0046034: ATP metabolic process, GO:0009205: purine ribonucleoside triphosphate metabolic process), cellular response to starvation (R-HSA-9711097: cellular response to starvation) (Fig. 6b,d-e and **Supplementary material 2**). TRRUST-based transcription factor enrichment analysis did not reveal any significant regulatory link (Fig. 6h).

Enrichment analysis of DEGs overexpressed in the non-severe group (vs. severe) revealed enrichment for RNA metabolic processes (e.g. R-HSA-8953854: metabolism of RNA, GO:0008380: RNA splicing, GO:0006397: mRNA processing), ribonucleoprotein complex biogenesis (GO:0022613: ribonucleoprotein complex biogenesis), SARS-CoV infection (e.g. R-HSA-9679506: SARS-CoV Infections, R-HSA-9679506: SARS-CoV Infections, R-HSA-913531: Interferon Signaling, R-HSA-128021: Cytokine Signaling in Immune system, R-HSA-909733: Interferon alpha / beta signaling), signaling by Rho GTPases (R-HSA-194315: Signaling by Rho GTPases) (Fig. 6c,f,g and **Supplementary material 3**). TRRUST-based transcription factor enrichment analysis revealed a regulatory link involving CIITA (p < 0.0001) (Fig. 6i).

Transcriptomic profiling of PB at T<sub>early</sub> revealed distinct transcriptional programs between severe and non-severe patients. PB from severe cases upregulated genes involved in cell metabolism and immune modulation, whereas PB from non-severe cases were enriched in antigen presentation and interferon response pathways.

### 2.5.2. PB transcriptomic signatures at T<sub>late</sub> time point according to disease severity

At T<sub>late</sub>, enrichment analysis of the 37 DEGs overexpressed in the severe group (vs. non-severe) revealed enrichment in SUMO E3 ligases (R-HSA-3108232: SUMO E3 ligases SUMOylate target proteins), TNF-α signaling via NFKB (M5890: HALLMARK TNFA SIGNALING VIA NFKB), and negative regulation of cell differentiation (GO:0045596). TRRUST-based transcription factor enrichment analysis did not reveal any significant regulatory links involved (Supplementary Results-Fig. 4). At T<sub>l</sub>, 8 DEGs were overexpressed in the non-severe group (vs. severe): PRR3, S1PR1, CBX7, RUBCNL, YEATS4, RIPOR2, FBXL17, and ARL14EP, but no enriched term was found.
<figure><img src="image_6.png" /><figcaption>b</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>e</figcaption></figure>

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /><figcaption>f</figcaption></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /><figcaption>g</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>c</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>d</figcaption></figure>

<figure><img src="image_14.png" /></figure>

(caption on next page)
**Fig. 4.** Modeling of the expansion of plasmablasts in COVID-19 in cohort 2. (a) Individual plasmablast (PB) trajectories (n = 166). Time 0 was defined by the date of the first PB proportion record. (b-c) The 2-class latent class mixed model of PB trajectories. (d) Survival curves according to PB trajectory (p < 0.001). Time 0 was defined by the date of the hospital admission. (e-f) CD39 expression (mean fluorescence intensity, MFI) on PB and other B cells according to trajectory classes (n = 6 for class 1 and n = 5 for class 2), and disease stage (early: <14 days; or late: >20 days from symptom onset), in patients with paired early and late data. Comparisons were performed using the paired Mann-Whitney test in (e) and the unpaired Mann-Whitney test in (f). (g). CD39 expression (MFI) on PB and other B cells at the early time point, stratified by disease severity (n = 12 for non-severe and n = 7 for severe/critical). Comparisons were performed using the unpaired Mann-Whitney test. *: p < 0.05, **: p < 0.01.

<caption>**Table 2**</caption>

Patient demographics and characteristics assigned to class 1 and class 2 in Latent Class Mixed Model (Cohort 2b).

<table><thead><tr><th></th><th>Class 1 (n = 150)</th><th>Class 2 (n = 16)</th></tr></thead><tbody><tr><td><strong>Demographics</strong></td><td></td><td></td></tr><tr><td>Male; n (%)</td><td>104 (69)</td><td>15 (94)</td></tr><tr><td>Age; median (IQR); y</td><td>61.6 (51.1; 69.7)</td><td>67.0 (59.5; 70.8)</td></tr><tr><td><strong>Medical background</strong></td><td></td><td></td></tr><tr><td>BMI; median (IQR); kg/m²</td><td>28.3 (25.5; 34.6)</td><td>27.0 (24.7; 30.5)</td></tr><tr><td>Diabetes; n (%)</td><td>34 (23)</td><td>4 (25)</td></tr><tr><td>Heart disease; n (%)</td><td>16 (10)</td><td>4 (25)</td></tr><tr><td>Chronic lung disease; n (%)</td><td>29 (19)</td><td>5 (31)</td></tr><tr><td>Cancer; n (%)</td><td>17 (11)</td><td>3 (19)</td></tr><tr><td>Severe liver disease; n (%)</td><td>8 (5)</td><td>2 (13)</td></tr><tr><td><strong>COVID-19 characteristics</strong></td><td></td><td></td></tr><tr><td>Days since disease onset; median (IQR)</td><td>13.0 (8.0; 19.0)</td><td>7.0 (6.0; 10.0)</td></tr><tr><td>Severity at first record<sup>a</sup>; n (%)</td><td></td><td></td></tr><tr><td>non-severe</td><td>40 (27)</td><td>6 (38)</td></tr><tr><td>severe</td><td>65 (43)</td><td>5 (31)</td></tr><tr><td>critical</td><td>45 (30)</td><td>5 (31)</td></tr><tr><td>0-7<sup>b</sup>; median (IQR)</td><td>4.0 (4.0; 5.0)</td><td>4.0 (3.0; 5.0)</td></tr><tr><td><strong>Biological results at inclusion</strong></td><td></td><td></td></tr><tr><td>PB trajectories</td><td></td><td></td></tr><tr><td>Mean baseline PB (%) [95CI]</td><td>16.3 [14.4; 18.1]</td><td>37.3 [14.4; 18.1]</td></tr><tr><td>Mean peak PB (%) [95CI]</td><td></td><td>40.9 [35.3; 46.2]</td></tr><tr><td>Mean day 30 PB (%) [95CI]</td><td>5.0 % [3.5; 6.5]</td><td>15.6 [9.7; 21.5]</td></tr><tr><td>Ferritin; ng/mL; median (IQR)</td><td>770.0 (351.0; 1614.0)</td><td>819.0 (428.0; 2139.0)</td></tr><tr><td>Fibrinogen; g/L; median (IQR)</td><td>6.2 (5.0; 7.3)</td><td>6.9 (5.5; 7.9)</td></tr><tr><td>LDH; U/L; median (IQR)</td><td>309.0 (279.0; 397.0)</td><td>415.0 (346.0; 475.0)</td></tr><tr><td><strong>Follow-up</strong></td><td></td><td></td></tr><tr><td>Maximal severity<sup>a</sup>; n (%)</td><td></td><td></td></tr><tr><td>non-severe</td><td>46 (31)</td><td>5 (31)</td></tr><tr><td>severe</td><td>38 (25)</td><td>2 (13)</td></tr><tr><td>critical</td><td>66 (44)</td><td>9 (56)</td></tr><tr><td>0-7<sup>b</sup>; median (IQR)</td><td>4.0 (4.0; 6.0)</td><td>5.5 (4.0; 7.0)</td></tr><tr><td>ICU admission; n (%)</td><td>72 (48)</td><td>8 (50)</td></tr><tr><td>ICU; median of days (IQR)</td><td>11.5 (5.5; 28.0)</td><td>9.5 (2.0; 32.5)</td></tr><tr><td>Death at 30 days from ICU admission; n (%)</td><td>7 (5)</td><td>5 (31)</td></tr><tr><td><strong>Treatment</strong></td><td></td><td></td></tr><tr><td>Steroids; n (%)</td><td>26 (17)</td><td>3 (19)</td></tr><tr><td>Tocilizumab; n (%)</td><td>0 (0)</td><td>0 (0)</td></tr><tr><td>Antiviral therapy; n (%)</td><td>34 (23)</td><td>1 (6)</td></tr></tbody></table>

## 2.6. CD39 expression on PB and other B cells

After identifying that PB are enriched in genes involved in purine metabolism in severe patients (compared to non-severe), we used FACS to quantify CD39 expression on sorted PB in 11 patients with paired data at T_early and T_late, stratified by class assignment (n = 6 for class 1 and n = 5 for class 2). CD39 MFI on PB was higher at the T_early of COVID-19 than at T_late (p = 0.03 for class 1 and p = 0.01 for class 2) (Fig. 4e). CD39 expression on PB was higher in class 2 compared to class 1 (p = 0.03) only at T_early (Fig. 4f). Regarding the maximal disease severity, we quantify CD39 expression on sorted PB in 19 patients with unpaired data for early and late time points. CD39 MFI on PB was higher in severe COVID-19 (n = 7) than in non-severe COVID-19 (N = 12) (p = 0.04) at T_early. We did not find any significant difference at T_late (Fig. 4g).

These results suggest that PB differ in terms of longitudinal trajectory

and functions according to class assignment and disease severity, at least during the early stage of disease.

## 3. Discussion

This study provides a comprehensive analysis of the interplay between PB expansion, the cytokine storm, and disease severity revealing a striking functional duality of PB during COVID-19. We identify circulating PB as central players in COVID-19 immunopathology, exhibiting divergent functional programs and longitudinal dynamics. Using longitudinal immune profiling, transcriptomic analysis, and trajectory modeling across two independent cohorts, we show that PB responses in COVID-19 are not only quantitatively distinct according to disease severity, but also qualitatively different.

### 3.1. Alterations in B-cell subsets and expansion of plasmablasts

We confirm a robust and early expansion of the PB subset in COVID-19 patients, particularly prominent in severe cases before normalizing during the convalescence stage consistent with previous reports in COVID-19 [15,16] and in other severe infections [12,17,18]. However, we noticed considerable heterogeneity of the dynamic of PB expansion, the importance of PB expansion, and patient survival suggesting different functional profiles of PB.

Classical memory B cells (both switched and unswitched), critical for long-term immunity, were markedly depleted in severe cases, at least during the incremental phase of infection as previously reported [19]. This depletion reflects a shift from germinal center responses to extrafollicular B-cell responses, characterized by rapid but short-lived antibody production. As previously reported, we observed a positive correlation between the DN2 B-cell proportion and the maximal disease severity further underscores the role of B-cell dysregulation in driving severe COVID-19 [13,19]. DN2 B cells are associated with extrafollicular B-cell responses and have been shown to exhibit pro-inflammatory phenotypes [10] and their expansion may therefore represent a maladaptive immune response that exacerbates tissue damage and systemic inflammation in severe COVID-19 and potentially an extrafollicular pathway to generate PB [15,20]. Not only is there an increase of the extrafollicular response, but our results suggest that PB have a pro-inflammatory functional program in severe patients that could contribute to disorganizing the immune response by promoting the inflammatory storm but also by serving as a nutrient sink, where B-cell differentiation is biased towards metabolically hyperactive, short-lived PB responses that by their sheer numbers constrain the germinal center reaction limiting the generation of memory B-cell and long-lived plasma cell responses, especially in patients with large PB expansion [21].

### 3.2. Interplay between cytokines, purine metabolism, and interferon signaling

During the early stage of COVID-19, PB from severe patients showed upregulation of genes related to cell proliferation and purine metabolism as previously reported [15]. This metabolic reprogramming supports antibody production but also drives early inflammation through extracellular ATP production [22,23]. ATP is converted to ADP and AMP by enzymes belonging to three ectonucleotidase families, including CD39 [22,23] [24]. ATP/ADP/AMP releases act as a pro-

<footnote>a Severity according to WHO classification (e.g., non-severe/severe/critical).</footnote>

<footnote>b 0-7 scale used in severity scoring (details in Methods).</footnote>
<figure><img src="image_15.png" /><figcaption>**Fig. 5.** Transcriptomic analysis of B cells during COVID-19 in cohort 2. Transcript expression profiles of plasmablasts (PB) and other B cells (BC) according to disease severity (NS: non severe, *n* = 24; SEV: severe or critical, *n* = 17) and to time of sampling (*T*<sub>early</sub>: Time point “early”, *n* = 21; *T*<sub>late</sub>: Time point “late”, *n* = 20). (a) Hierarchical Clustering Analysis of all the expressed genes. Hierarchical clustering was constructed with Euclidean distance using R software. (b) Number of significantly DEG (under-expressed and over-expressed) for q-value <0.1 and absolute fold change (FC) >1.5. (c-e) Principal component analysis of DEG (q-value <0.1 and FC >1.5) by time of sampling (*T*<sub>early</sub> or *T*<sub>late</sub>) and cell type (PB or BC) (c), and by cell type (PB or BC) and survival (d-e).</figcaption></figure>

inflammatory mediator by binding purinergic receptors (e.g., P2X7R), triggering inflammatory cytokine release, perpetuating the inflammatory storm, and contributing to the thromboinflammation [23,25]. In line with these results, Diaz-Garcia et al. previously reported higher concentrations of plasma ATP and ADP in severe COVID-19 and Influenza A patients compared to healthy controls, reinforcing a purinergic signaling dysregulation, which might be underlying COVID-19 thromboinflammation disorder. The subsequent phase of inflammation is associated with reduced ATP release and the induction of ectonucleotidases that decrease the ATP/adenosine ratio. Adenosine thus released contributes secondarily to limiting the scale and duration of the inflammatory response [22,25]. So, the massive expansion of PB in severe patients could contribute to the inflammatory storm during the early stage of COVID-19, but the prolonged expansion of PB observed in severe patients possibly contributing to the long-term immunosuppression, as observed in sepsis survivors [4] through the generation of large amounts of extracellular adenosine [7,22,23,25].

Conversely in non-severe cases, PB showed upregulation in genes related to interferon signaling pathways, particularly type I and II

interferons, promoting antiviral responses by enhancing antigen presentation and adaptive immunity. By longitudinal multi-omics modeling, Bernardes et al. also identified a strong increase in IFN-related genes during the early stage of COVID-19, which was suppressed in severe COVID-19 in consistency with our findings [15]. Furthermore, we identified a relative upregulation of MHC class II transactivator (CIITA) as a key transcription factor in non-severe patients. CIITA has been reported to induce resistance to SARS-CoV-2 by up-regulation of the p41 isoform of CD74, which blocks cathepsin-mediated cleavage of viral GPs, thereby preventing viral fusion [26]. Low blood CIITA expression levels were positively correlated with greater severity in a prospective observational study [27]. PBs are not targeted by SARS-CoV-2 and do not contribute to limiting viral propagation in this manner. But relative upregulation of CIITA in non-severe patients PB could reflect the cytokine environment that contributes to optimize the antiviral response in immune cells. PB themselves could also contribute to viral control and limiting inflammation by facilitating T-cell-mediated immunity.

These contrasting profiles raise an important question: are PB in
<figure><img src="image_16.png" /><figcaption>**Fig. 6.** Enrichment analysis on DEGs in plasmablasts from severe patients compared to non-severe patients at T<sub>early</sub> in cohort 2. Transcript expression profiles of PB (n = 15 samples) according to disease severity (NS: non severe or SEV: severe or critical) at time point “early” (T<sub>early</sub>). (a) Heatmap summarizing differentially expressed genes (DEG) in PB (q-value < 0.1 and absolute fold change >1.5) rows (genes) are centered and scaled. Hierarchical clustering was constructed with Euclidean distance using R software. (b-c) Enrichment analysis summarizing the most representative up-regulated pathways related to DEG upregulated (n = 154 genes) in PB from severe patients (compared to non-severe patients) (b) and to DEG upregulated (n = 579 genes) in PB from non-severe patients (compared to severe patients) (c) using Metascape (http://metascape.org/gp/index.html#/main/step1). (d-e) All significantly enriched terms in PB from severe patients (compared to non-severe patients) were hierarchically clustered into a tree based on Kappa-statistical similarities among their gene memberships. A subset of representative terms from these clusters have been converted into a network layout. Each term is represented by a circle node, whose size is proportional to the number of input genes associated to that term, and whose color represents a cluster identity. Terms with a similarity score > 0.3 are linked by an edge. Each cluster was named according to the terms memberships. The same enrichment network with nodes colored following a *p*-value order. (f-g) The same analysis was performed on the significantly enriched terms in PB from non-severe patients (compared to non-severe patients). The enrichment analysis was performed using Metascape (http://metascape.org/gp/index.html - /main/step1). (h-i) Enrichment analysis of DEGs upregulated in PB from severe patients (h) and of DEG upregulated in non-severe patients using TRRUST (i).</figcaption></figure>

severe cases hyperactivated as a consequence of the inflammatory
milieu, or do they actively drive disease progression? While cytokine
feedback undoubtedly influences PB activity, the ability of PB to amplify
inflammation via purine metabolism suggests they may be a key
amplifier of disease severity.

3.3. *Correlations between plasmablasts, cytokines and disease severity*

The interaction between expanded PB and cytokines emerged as a
critical factor in shaping the immune response to SARS-CoV-2. PB pro-
portion correlated with key inflammatory mediators, especially from D7
onward, emphasizing that PB might play a pivotal role in perpetuating
the inflammatory cascade [19,28,29], which could in part be driven by
the increased purine metabolism observed in PB from severe patients at
the early stage of COVID-19 [9]. This highlights a central paradox: while
PB are critical for rapid antibody production and viral clearance, their
hyperactivation (in both intensity and duration) might contribute to
systemic inflammation. Furthermore, the association between PB
expansion and markers of endothelial damage, such as D-dimer, suggests

that PB-driven inflammation may extend beyond immune dysregulation
to include vascular pathology—a hallmark of severe COVID-19.

One of the most striking findings is the tight correlation between BAFF levels and maximal disease severity at D7 in the prospective cohort. Blood BAFF levels—a cytokine essential for PB survival [30]— were higher in COVID-19 patients compared to controls from baseline onward, as previously reported in COVID-19 patients with high numbers of PB (>20 %) [9]. Interestingly, the blood BAFF levels at D7 correspond to the transient fall before rising again at D14, which can be explained by the use of dexamethasone [14], although we cannot rule out that BAFF may have been consumed by PB expansion. BAFF-mediated B-cell activation appears to fuel PB expansion. Targeting BAFF may help mitigate PB-driven inflammation in early stage, PB-driven immunosuppression during convalescence stage of COVID-19 and thus reduce disease severity. Belimumab, a monoclonal antibody targeting BAFF, was developed to inhibit B cell survival. It was the first biologic therapy approved specifically for systemic lupus erythematosus. While no data are available on belimumab to treat COVID-19, data from the COVID-19 Global Rheumatology alliance bring safety issues. A total of 1606
patients with systemic lupus and COVID-19 were included. In the age-
adjusted and sex-adjusted model, belimumab was associated with a
good outcome (OR (95 %CI) = 0.51 (0.31 to 0.85), *p* < 0.001) [31].
Given BAFF's role in PB activation and survival, and its association with
disease severity, we propose that a dysregulated BAFF-PB axis may be a
key driver of pathogenic PB responses in COVID-19. Targeting this axis
could represent a novel therapeutic strategy to prevent the inflamma-
tory storm when PB expansion is observed. Further studies are needed to
determine whether targeting BAFF in patients with PB expansion pre-
vents disease severity.

Quantification of PB at baseline together with the assessment of CD39 expression on PB may provide valuable information for evaluating the risk of severity of COVID-19. Patients exhibiting a marked and sustained PB expansion and/or baseline high levels of CD39 expression appear to represent a subgroup at higher risk of severe outcomes in COVID-19 but also potentially in other severe infectious diseases. For example, patients with sepsis exhibited expanded CD39hi PB [4]. These observations highlight the potential diagnostic and prognostic value of integrating PB dynamics and CD39 expression into patient stratification strategies.

3.4. *Limitations of study*

The widespread use of steroids in the prospective cohort may have
impacted our longitudinal assessment of B-cell subsets and serum levels
of soluble factors. Since steroids are currently recommended as a first-
line therapy for severe patients, their effects could not be eliminated.
This cohort also included a higher proportion of critically ill patients
compared to milder forms, potentially limiting statistical power. In
contrast, the retrospective cohort offered the opportunity to assess PB
function with minimal steroid exposure, but its design constrained the
scope of analyses. Sampling at early and late time points did not always
involve the same individuals due to variable sample quality, preventing
individual-level modeling of PB dynamics. The limited material also
precluded B-cell receptor sequencing, which could have clarified PB
ontogeny, particularly regarding potential severity-dependent origins.
In addition, the metabolic reprogramming of PB observed here requires
confirmation by orthogonal approaches. A limitation of our study is the
lack of validation in an independent, treatment-free cohort to confirm
the role of BAFF in PB expansion and disease severity, which warrants
further investigation in COVID-19 and other severe infections. Finally,
whether targeting PB could modulate hyperinflammation in severe
COVID-19 warrants dedicated interventional studies.

To conclude, PB may play a central role in the immunopathology of
COVID-19, with divergent profiles reflecting and influencing disease
severity. In severe cases, PB act as amplifiers of inflammation through
purine metabolism and cytokine signaling, while in non-severe cases,
they support effective antiviral immunity via interferon pathways. The
interplay between PB, cytokines, and B-cell subsets underscores their
dual role as both mediators and markers of immune dysregulation.
Understanding these dynamics is key to developing targeted therapies to
modulate PB activity and improve outcomes in severe viral infections.
These results underscore the importance of PB not only as biomarkers of
disease severity but also as potential therapeutic targets. In particular,
the strong correlation between BAFF levels and disease severity supports
the rationale for targeting the BAFF-PB axis in selected patients when PB
expansion occurs. Modulating PB function or differentiation could
represent a novel strategy to limit hyperinflammation without
compromising antiviral immunity.

4. Methods

4.1. *Patients and controls*

4.1.1. *Sex as a biological variable*

Study participants were selected without consideration of their sex.

4.1.2. Cohort 1 (“Immuno-Cov2 cohort”)

The cohort 1 prospectively included 50 patients with COVID-19 admitted to the University Hospital in Lille, France between July 2020, and December 2020. Patients were included if they had to have symptomatic SARS-Cov2 infection confirmed by positive nasopharyngeal swab polymerase chain reaction. Data on demographics, clinical manifestations, routine laboratory features, and therapeutic management, and outcomes were collected prospectively. Blood samples had been collected at D1, D7, D14, D30 et D120. The first sample was obtained before administration of corticosteroids. The alpha (B.1.1.7) SARS-CoV-2 variant was dominant in France during this period. As a control cohort, 18 healthy donor samples were used as controls.

4.1.3. Cohort 2 (“Plasmacovid-19 cohort”)

All patients with a laboratory-confirmed diagnosis of COVID-19 who were hospitalized at the CHU Lille and had available lymphocyte immunophenotyping performed between April 8th and May 1st, 2020, were retrospectively enrolled. Patients had to have symptomatic SARS-Cov2 infection confirmed by positive nasopharyngeal swab polymerase chain reaction. Samples used for lymphocyte immunophenotyping were stored for B-cell isolation and DNA extraction. Data on demographics, clinical manifestations, routine laboratory features, and therapeutic management, and outcomes were collected. The alpha (B.1.1.7) SARS-CoV-2 variant was dominant in France during this period. We selected patients from cohort 2 who had at least two records of PB proportion over time (cohort 2b) for Latent Class Mixed Modeling. Among them, we retrospectively selected 41 samples from 23 steroid-naive patients with blood samples available according to maximum severity, time since onset of disease ("early" or "late" time points, respectively less than 14 days or more than 20 days from disease onset), availability of samples, and sufficient cells number per sample to perform RNA sequencing on PB cells and other B cells.

4.1.4. *Definition of variables*

Epidemiological, clinical, laboratory and radiological characteristics
were obtained from patients' electronic medical records. The COVID-19
severity was assessed at each visit according to (a) a 7-point ordinal
clinical scale of the WHO Master Protocol [32] defined as: 1. Not hos-
pitalized – no limitation on activities; 2. Not-hospitalized, limitation on
activities; 3. Hospitalized, not-requiring supplemental oxygen; 4. Hos-
pitalized, requiring supplemental oxygen; 5. Hospitalized, on non-
invasive ventilation or high flow oxygen devices; 6. Hospitalized, on
invasive mechanical ventilation or ECMO; and 7. Death, and (b) the
COVID-19 WHO disease severity [33]. Disease severity was defined as (i)
critical COVID-19 defined as fulfilling the criteria for acute respiratory
distress syndrome, sepsis, septic shock, or other conditions that would
normally require the provision of life-sustaining therapies such as me-
chanical ventilation (invasive or non-invasive) or vasopressor therapy;
as (ii) as severe COVID-19 defined by any of oxygen saturation < 90 %
on room air, signs of pneumonia, signs of severe respiratory distress
(accessory muscle use, inability to complete full sentences, respiratory
rate > 30 breaths per minute); (iii) as non-severe COVID-19 defined by
the absence of any criteria for severe or critical COVID-19.

4.2. Flow cytometry

4.2.1. B-cell immunophenotyping (cohort 1)

B-cell immunophenotyping was performed on fresh whole blood. Blood samples were washed twice with PBS-FBS 5 % (GIBCO) and incubated with antibodies listed in *Supplementary Method-Table 1* for 15 min at room temperature. After addition of one volume of FBS, red blood cells were lysed, and white blood cells were fixed using Immunoprep (BECKMAN COULTER). Cells were washed twice with PBS-FBS 5 % (GIBCO) and were analyzed using a CytoFLEX (BECKMAN COULTER). To determine absolute counts of cells, we used a routine leukocyte immunophenotyping: unwashed blood samples were incubated with the
following antibodies: anti-CD45, anti-CD3, anti-CD4, anti-CD8, anti-CD19, anti-CD56 and anti-CD16 (BECKMAN COULTER), then, red blood cells were lysed, and white blood cells were fixed using Immuno prep (BECKMAN COULTER). After addition of Flow-count fluorospheres (BECKMAN COULTER), cells were analyzed using a Navios flow cytometer (BECKMAN COULTER). Data were analyzed using Kaluza software (BECKMAN COULTER). The gating strategy for B-cell subsets was described in **Supplementary Method**-Fig. 1. The different B-cell subsets and corresponding phenotypes are presented in **Supplementary Method-Tables 1–2**.

### 4.2.2. Blood plasmablast quantification (cohort 2)

Routine leukocyte immunophenotyping (the same as described in chapter B-cell immunophenotyping (cohort 1)) was reanalyzed (by TG, EL and SD) to quantify PB based on Forward Scatter (FSC) and CD19 expression: FSChighCD19low. The gating strategy was described in Method-Fig. 2A. This gating strategy was validated on cohort 1 by comparison to the classical and more complete phenotypic definition of PB (CD19+/lowCD24-CD38high) (**Supplementary Method**-Fig. 2b).

### 4.2.3. Plasmablast sorting and CD39 expression analysis (cohort 2b)

Peripheral blood mononuclear cells (PBMCs) were isolated from venous blood samples via standard density gradient centrifugation and used after cryopreservation. Cells were thawed and washed twice using RPMI-1640 FBS 10% (GIBCO) DNase 10Kunitz/mL (SIGMA-ALDRICH). Then, cells were resuspended in cytometry buffer: PBS FBS 2% (GIBCO) DNase 10Kunitz/mL (SIGMA-ALDRICH) and incubated with LIVE/DEAD Fixable Yellow Stain Kit (THERMO FISHER SCIENTIFIC) and antibodies listed in **Supplementary Method-Table 3** during 30 min on ice. Cells were washed in a cytometry buffer and then analyzed and sorted using a FACS Aria Cell Sorter (BD BIOSCIENCES). PB and non-PB B cells were analyzed and sorted directly in the lysing buffer (MACHEREY-NAGEL) according to the gating strategy presented in **Supplementary Method**-Fig. 3.

### 4.3. Cytokines' analysis

Serum levels of APRIL, BAFF, Lactoferrin, Myeloperoxidase (MPO), Matrix metalloproteinase-8 (MMP8), Proteinase 3 (PR3), CCL2, CCL11, CCL17, CCL24, CCL26, CCL27, CXCL1, CXCL2, CXCL8/IL-8, CXCL10, CXCL11, CXCL13, CX3CL1, soluble CD27 (sCD27), soluble CD44 (sCD44), GM-CSF, IL-3, IL-4, IL-5, IL-10, IL-13, IL-21, Lipocalin-2/NGAL and Thymic stromal lymphopoietin (TSLP) were assessed by multiplexed bead assay (Luminex, R&D SYSTEMS) according to manufacturer's instructions. Serum levels of IL-6 were assessed by Ella Automated Immunoassay (R&D SYSTEMS) according to manufacturer's instructions.

## 4.4. Transcriptomics

### 4.4.1. RNA extraction

Total RNAs were extracted using the NucleoSpin® RNA kit (MACHEREY-NAGEL) as per the manufacturer's instructions.

### 4.4.2. RNA sequencing data analysis

The samples exhibited heterogeneity in total RNA quantity. To address this, we concentrated the samples with the lowest RNA concentrations using a SpeedVac and limited the input to a maximum of 50 ng for the most concentrated ones (n = 2) to avoid an excessive concentration range. After concentration, the estimated total RNA amounts ranged from 0.1 ng to 50 ng in a volume of 4 µL. From 4 µL of total RNA, we added 1 µL of ERCC spike-in control. Library preparation was then initiated using oligo(dT) priming, following the QuantSeq FFPE 3' mRNA FWD protocol with UDI (Unique Dual Indexes) from Lexogen, incorporating the low-input specific workflow. Unique Molecular Identifiers (UMIs) were introduced to facilitate the removal of PCR duplicates during analysis. After generating the double-stranded cDNA

library, we purified it using magnetic beads and amplified it. Library amplification was performed with 22 cycles. Libraries with the lowest yields (less than 600 pM in the 200 bp–1000 bp range) were re-amplified using the reamplification add-on kit for Illumina (LEXOGEN) with six to eight additional cycles. Final library concentrations ranged from 222 pM to 11,187 pM (200 bp–1000 bp interval). After equimolar pooling of the libraries, the pool was repurified to eliminate primer dimers. The libraries were sequenced on a NovaSeq 6000 system using a 100-cycle chemistry. All samples, except one, achieved at least 20 million reads. To eliminate low-quality regions and poly(A) tails from the reads, we used the fastp program with a quality score threshold of 20, removing reads shorter than 25 bp. Read alignments were performed using the STAR program with the human reference genome GRCh38 and reference gene annotations from Ensembl. UMIs were utilized to reduce errors and PCR bias through fastp and umi-tools. Based on deduplicated read alignments (UMI), gene-level molecule counts were determined using FeatureCounts. Additional programs such as qualimap, FastQC, and MultiQC were used for quality control of the reads.

### 4.4.3. Gene expression analysis

Differential gene expression analysis of RNA-seq data was performed using the R/Bioconductor packages DESeq2 and fdrtool [34,35]. Differentially expressed genes were identified using a *p*-value adjustment (FDR) threshold of <0.1. RNA-seq samples were obtained following cell sorting. In some cases, we had only limited cell numbers and consequently low RNA yields. A preliminary analysis was performed on 5 samples: the RIN values were more than 7. As a result, we chose to skip the qualification step in order to preserve as much material as possible. Samples exhibiting no detectable expression of housekeeping genes RPL13A and RPLP0 were deemed to contain insufficient RNA quantities for accurate transcriptomic analysis and were consequently excluded from the dataset. This strategy led to the exclusion of eight samples from the analysis. Moreover, due to the low number of molecules in some samples, normalization was performed using the following housekeeping genes: ACTB, HMBS, HPRT1, SDHA, TBP, RPL13A, GAPDH, B2M, RPLP0, IGKC, JCHAIN, RPS27A, and RPL27A. Additionally, bias correction was applied using the SVA package [36]. We did not perform any age- and sex-matching in the RNA-seq analysis.

### 4.4.4. Enrichment

To analyze Gene Oncology (GO), REACTOM, and KEGG pathways, we subjected differentially expressed genes (DEGs) (or proteins, DEPs) with cut-off *q*-value at <0.05 to the online tool Metascape (http://metascape.org/gp/index.html#/main/step1) [37]. All significant enriched terms (*p*-value <0.05) were hierarchically clustered into a tree based on Kappa-statistical similarities (threshold >0.3) among their gene membership. The clusters identified have been named according to the membership. Among the terms within each cluster, some of them (*q*-value <0.05) were selected for charting.

## 4.5. Statistics

The continuous variables were expressed by the mean ± SD or the median with interquartile range (IQR). The normality of distribution was checked graphically and using the Shapiro-Wilk test. The categorical variables were expressed as frequencies and percentages. Categorical variables were expressed as frequencies and percentages. Comparisons of groups were performed using the analysis of variance or the Kruskal-Wallis test for quantitative variables and the Fisher's exact test or the chi-square test for categorical variables. In the case of significant results, pairwise comparisons were performed, and a Bonferroni correction was applied. The correlations were assessed using the Spearman correlation coefficient. The overall survival was described using the Kaplan-Meier method and compared between the latent classes by using the Log-rank test. Then the comparison was adjusted for age and sex by using the Cox proportional hazards model.
To delineate groups of patients according to their PB trajectories (classes) as measured by PB proportion over time using latent class mixed models (LCMM) [38,39]. LCMM assumes that the population is divided into a finite number of groups called latent classes. Each latent class is characterized by a specific mean trajectory, which is described by a class-specific linear mixed model. In a given latent class, the individual trajectories of patients are close to each other, while the individual trajectories of patients of different classes tend to be dissimilar. In LCMM, latent classes correspond to an unknown categorical variable which is identified from data using a multinomial logistic model, and trajectories of PB are analyzed using the mixed model with random coefficients to take into account individual trajectories. The random effects (linear or quadratic) are determined from the analysis of residuals according to Verbeke and Molenberghs [40]. To identify the number of classes, several LCMMs are performed. Each model predicts the shape of the trajectory of each class, estimates the probability for each individual of class membership and assigns each of them to the class for which the likelihood is the highest. Time 0 was defined by the date of the first PB proportion measurement. Trajectories were censored 55 days after the first PB measurement, which corresponds to the 95th percentile of the time between the first and last PB measurements. The time between first symptoms and first measurement was considered as an adjustment factor. To determine the best number of latent classes, we considered both formal statistical criteria (such as Bayes information criterion (BIC)) and model adequacy. A low value of BIC and average posterior probabilities of class membership greater than 0.8 correspond to the better model [41,42]. Two sensitivity analyses were performed to identify latent classes: the first in patients with at least 3 measurements of PB; the second in patients with at least two measurements of PB and without steroids use.

All statistical tests were performed at the 2-tailed α level of 0.05, and data analysis was performed by using SAS software version 9.4 (SAS Institute Inc., Cary, NC, USA), GraphPad Prism version 10.4.2 (GraphPad Software, San Diego, CA, USA) and R software.

## 4.6. Study approval

This study was approved by the French institutional authority for personal data protection (Commission Nationale de l'Informatique et des Libertés (CNIL), registration number DEC20-086), and ethics committee (ID-CRB 2020-A00763-36). The protocol was registered as a clinical trial (registration numbers NCT-04327180 and NCT-04341792 for patients admitted to ICU and patients who consulted in the emergency department, respectively).

Cohort 1 patients were from the Lille COVID Research Network (LICORNE) cohort (ethics committee approval number ECH20/09). Informed consent with participation in the cohort 1 study was signed by the participants in accordance with the Declaration of Helsinki.

## CRediT authorship contribution statement

**Emmanuel Ledoult:** Writing – review & editing, Writing – original draft, Visualization, Validation, Project administration, Methodology, Investigation, Formal analysis, Data curation, Conceptualization.
**Thomas Guerrier:** Writing – review & editing, Writing – original draft, Visualization, Validation, Resources, Methodology, Investigation, Formal analysis, Data curation, Conceptualization.
**Sylvain Dubucquoi:** Writing – review & editing, Visualization, Validation, Resources, Methodology, Conceptualization.
**Martin Figeac:** Writing – review & editing, Software, Methodology, Formal analysis.
**Céline Villenet:** Software, Methodology, Formal analysis.
**Blanche Daunou:** Writing – review & editing, Formal analysis, Data curation.
**Helene Behal:** Writing – review & editing, Visualization, Software, Methodology, Formal analysis.
**Mohammad Ryadh Pokeerbux:** Writing – review & editing, Visualization, Software, Methodology.
**Thomas Machet:** Writing – review & editing, Validation, Investigation, Formal analysis.

**Vincent Koether:** Writing – review & editing, Visualization, Validation, Investigation.
**Aurore Collet:** Writing – review & editing, Visualization, Validation, Methodology, Conceptualization.
**David Launay:** Writing – review & editing, Visualization, Validation, Supervision, Software, Resources, Methodology, Investigation, Funding acquisition, Formal analysis, Conceptualization.

## Funding

European Union (REACT-EU, ERDF, File No. NP 0032488), and the French government (I-SITE ULNE/ANR-16-IDEX-0004 ULNE).

## Declaration of competing interest

The authors have declared that no conflict of interest exists.

## Acknowledgments

We would like to thank Dr. Frederic Lepetre for his help with preliminary data analysis of transcriptomic data, Jacques Trauet, Arnaud Dendooven and Dr. Julie Demaret for their help for immunophenotyping and the cytokine assays. This work was funded by the European Union's REACT-EU program, under the European Regional Development Fund (ERDF), as part of the European Union's response to the COVID-19 crisis (File No. NP 0032488). This work was supported by the French government through the Programme Investissement d'Avenir (I-SITE ULNE/ANR-16-IDEX-0004 ULNE) managed by the Agence Nationale de la Recherche.

## Appendix A. Supplementary data

Supplementary data to this article can be found online at https://doi.org/10.1016/j.clim.2025.110609.

## Data availability

Transcriptomics (RNA seq) datasets related to this article can be found at https://www.ebi.ac.uk, an open-source online data repository hosted at NCBI/GEO. Number: E-MTAB-15068. Additional information and detailed methods are available upon request to the corresponding author.

## References

[1] M. Akkaya, K. Kwak, S.K. Pierce, B cell memory: building two walls of protection against pathogens, Nat. Rev. Immunol. 20 (4) (2020) 229-238.

[2] K. Fink, Origin and function of circulating Plasmablasts during acute viral infections, Front. Immunol. 3 (2012) 78.

[3] C. Aggarwal, et al., Immunophenotyping and transcriptional profiling of human plasmablasts in dengue, J. Virol. 95 (23) (2021), https://doi.org/10.1128/jvi.00610-21.

[4] D.C. Nascimento, et al., Sepsis expands a CD39+ plasmablast population that promotes immunosuppression via adenosine-mediated inhibition of macrophage antimicrobial activity, Immunity 54 (9) (2021) 2024-2041.e8.

[5] A.K. McElroy, et al., Human Ebola virus infection results in substantial immune activation, Proc. Natl. Acad. Sci. 112 (15) (2015) 4719-4724.

[6] N.F. Ahsan, et al., The current understanding of the phenotypic and functional properties of human regulatory B cells (Bregs), Oxf. Open Immunol. 5 (1) (2024) iqae012.

[7] F. Figueiró, et al., Phenotypic and functional characteristics of CD39high human regulatory B cells (Breg), Oncolmmunology 5 (2) (2016) e1082703.

[8] M. Felisberto, et al., Lymphocyte B subtypes in peripheral blood: a prognostic biomarker for COVID-19 patients, J. Appl. Lab Med. 9 (3) (2024) 456-467.

[9] C. Schultheiß, et al., Maturation trajectories and transcriptional landscape of plasmablasts and autoreactive B cells in COVID-19, iScience 24 (11) (2021) 103325.

[10] M.C. Woodruff, et al., Extrafollicular B cell responses correlate with neutralizing antibodies and morbidity in COVID-19, Nat. Immunol. 21 (12) (2020) 1506-1516.

[11] C. Ruschil, et al., Specific induction of double negative B cells during protective and pathogenic immune responses, Front. Immunol. 11 (2020) 606338.

[12] Y. Sun, et al., Early B lymphocyte subsets in blood predict prognosis in sepsis, Front. Immunol. 15 (2024) 1437864.
[13] R. Cervantes-Díaz, et al., Severity of SARS-CoV-2 infection is linked to double-negative (CD27- IgD-) B cell subset numbers, Inflamm. Res. 71 (1) (2022) 131-140.

[14] V.A. Sosa-Hernández, et al., B cell subsets as severity-associated signatures in COVID-19 patients, Front. Immunol. 11 (2020) 611004.

[15] V. Hébert, et al., Modifications of the BAFF/BAFF-receptor Axis in patients with pemphigus treated with rituximab versus standard corticosteroid regimen, Front. Immunol. 12 (2021) 666022.

[16] J.P. Bernardes, et al., Longitudinal multi-omics analyses identify responses of megakaryocytes, erythroid cells, and Plasmablasts as hallmarks of severe COVID-19, Immunity 53 (6) (2020) 1296-1314.e9.

[17] L. Scharf, et al., Longitudinal single-cell analysis of SARS-CoV-2-reactive B cells uncovers persistence of early-formed, antigen specific clones, JCI Insight 8 (1) (2022) e165299.

[18] Butowt R, Pyrc K, Bartheld CS von. Battle at the entrance gate: CIITA as a weapon to prevent the internalization of SARS-CoV-2 and Ebola viruses. Signal. Transduct. Target. Ther. 2020;5(1):278.

[19] P. Nayak, K. Mukund, S. Subramaniam, The Janus face of proliferating plasmablasts in dengue and COVID-19 infections, Front. Immunol. 14 (2023) 1068424.

[20] C.E. Faliti, et al., Disease-associated B cells and immune endotypes shape adaptive immune responses to SARS-CoV-2 mRNA vaccination in human SLE, Nat. Immunol. 1-15 (2024).

[21] R. Vijay, et al., Infection-induced plasmablasts are a nutrient sink that impairs humoral immunity to malaria, Nat. Immunol. 21 (7) (2020) 790-801.

[22] C. Cekic, J. Linden, Purinergic regulation of the immune system, Nat. Rev. Immunol. 16 (3) (2016) 177-192.

[23] E. Díaz-García, et al., Role of CD39 in COVID-19 severity: dysregulation of purinergic signaling and Thromboinflammation, Front. Immunol. 13 (2022) 847894.

[24] A.L. Horenstein, et al., A CD38/CD203a/CD73 ectoenzymatic pathway independent of CD39 drives a novel adenosinergic loop in human T lymphocytes, OncoImmunology 2 (9) (2013) e26246.

[25] J. Linden, F. Koch-Nolte, G. Dahl, Purine release, metabolism, and signaling in the inflammatory response, Annu. Rev. Immunol. 37 (1) (2016) 1-23.

[26] A. Bruchez, et al., MHC class II transactivator CIITA induces cell resistance to Ebola virus and SARS-like coronaviruses, Science 370 (6513) (2020) 241-247.

[27] M. Girona-Alarcon, et al., Low levels of CIITA and high levels of SOCS1 predict COVID-19 disease severity in children and adults, iScience 25 (1) (2022) 103595.

[28] D.E. Leisman, et al., Cytokine elevation in severe and critical COVID-19: a rapid systematic review, meta-analysis, and comparison with other inflammatory syndromes, Lancet Respir. Med. 8 (12) (2020) 1233-1244.

[29] P. Liang, et al., Prognostic significance of serum interleukin-6 in severe/critical COVID-19 patients treated with tocilizumab: a detailed observational study analysis, Sci. Rep. 14 (1) (2024) 29634.

[30] G. Jego, R. Bataille, C. Pellat-Deceunynck, Interleukin-6 is a growth factor for nonmalignant human plasmablasts, Blood 97 (6) (2001) 1817-1822.

[31] D.T. Avery, et al., BAFF selectively enhances the survival of plasmablasts generated from human memory B cells, J. Clin. Invest. 112 (2) (2003) 286-297.

[32] M.F. Ugarte-Gil, et al., Characteristics associated with poor COVID-19 outcomes in individuals with systemic lupus erythematosus: data from the COVID-19 global rheumatology Alliance, Ann. Rheum. Dis. 81 (7) (2022) 970-978.

[33] F. Ader, et al., Protocol for the DisCoVeRy trial: multicentre, adaptive, randomised trial of the safety and efficacy of treatments for COVID-19 in hospitalised adults, BMJ Open 10 (9) (2020) e041437.

[34] WHO Clinical Management of Covid-19—Interim Guidance, 2025.

[35] B. Klaus, K. Strimmer, Signal identification for rare and weak features: higher criticism or false discovery rates? Biostatistics 14 (1) (2013) 129-143.

[36] K. Strimmer, fdrtool: a versatile R package for estimating local and tail area-based false discovery rates, Bioinformatics 24 (12) (2008) 1461-1462.

[37] J.T. Leek, J.D. Storey, Capturing heterogeneity in gene expression studies by surrogate variable analysis, PLoS Genet. 3 (9) (2007) e161.

[38] Zhou Y, et al. Metascape Provides a Biologist-Oriented Resource for the Analysis of Systems-Level Datasets. Nature.com. Doi: https://doi.org/10.1038/s41467-01-9-09234-6.

[39] C.E. McCulloch, et al., Discovering subpopulation structure with latent class mixed models, Stat. Med. 21 (3) (2002) 417-429.

[40] G. Verbeke, G. Molenberghs, Linear Mixed Models for Longitudinal Data, Springer Ser Stat, 2000, https://doi.org/10.1007/b98969.

[41] J.B. Schreiber, Latent class analysis: an example for reporting results, Res. Soc. Adm. Pharm. 13 (6) (2017) 1196-1201.

[42] R. van de Schoot, et al., The GRoLTS-checklist: guidelines for reporting on latent trajectory studies, Struct. Equ. Model. Multidiscip. J. 00 (00) (2017) 1-17.