The advent of CITE-Seq technology allows for capturing single cell RNA sequencing along with cell surface protein levels with antibodies conjugated a DNA-barcode, analogous to the fluorophore of flow cytometry antibodies. This allows for the quantitative and qualitative information on surface proteins with available antibodies on a single cell level, with matched RNA-seq, and B cell receptor (BCR) and T cell receptor (TCR) VDJ information (16). This allows us to perform a flow cytometric-style gating of the B cell populations (using CITE-seq) and compare this to the gene expression (GEX) and BCR sequencing (BCR-seq) information. We used data from the COMBAT study (17), which represents a comprehensive single cell multi-omic blood atlas encompassing acute patients with varying COVID-19 severity (18 critical, 20 severe, 12 mild and 12 convalescent), 10 influenza patients, 15 hospitalised sepsis and 10 healthy controls (sampled pre-pandemic). Integrative multi-omics analysis of scRNA-seq, CITE-seq and BCR/TCR-seq allowed for high confidence and quality annotations of B cell, T cell and myeloid populations, as outlined in (17) and characterised in Supplementary Figure S1 . Briefly, we first performed separate clustering of gene expression, clustering of surface protein expression, and analyses of T and B cell receptor V(D)J sequences [described fully in (17)]. Cell types and subsets were further refined using information from the BCR-seq, CITE-seq and GEX layers for each GEX cluster phenotype led by expert understanding of each immune cell subset, considering a combination of marker genes and transcription factors. Information from all three modalities was used to identify and exclude doublets from downstream analysis. In agreement with the literature, activation markers [CD69, CD80, CD86, CD70, and CD24 (18–20)], cytokines and cytokine receptors [IL-2R, IL-21R, and CXCR3 (21–25)] are elevated in memory and plasmablast populations compared to naïve (20), and IgD, CD21 and CD23 are downregulated (26–29) ( Supplementary Figures S1A–C ). Furthermore, plasmablast/plasma cell-specific transcription factors [IRF4, PRDM1 (BLIMP1), BCL2L1 and XBP1] are observed only in plasmablast populations, whereas early B cell stage TFs [BACH2, PAX5, MCL1 and BCL6 (30, 31)] are down-regulated in plasmablast populations. The GC-stage-specific TF, MYC, is seen highest in memory B cells as expected and decreases upon plasma cell differentiation (32). The transitional and naïve B cells contained no SHM and only unswitched BCRs (IgD/M), whereas the memory and plasmablast populations contained SHM and/or class-switched sequences ( Supplementary Figures S1D, E ). Finally, the level of expression of the heavy chain sequence (nUMIs) and expression of the J-chain was significantly elevated in the plasmablast population compared to the other B cell subsets ( Supplementary Figures S1B, D ) in agreement with elevated production of immunoglobulins (33).

Finally, a classical flow cytometric-style gating strategy was performed using multi-omics CITE-seq levels to define naïve, CD27+ IgM- (switched) memory, CD27+ IgM+ (unswitched) memory, IgD- CD27- B cell, CD27+ plasmablast, CD27+ IgM+ plasmablast, and IgD- CD27- plasmablast populations ( Figure 1A ; Supplementary Table S1 , see methods). These FACS-style gated B cell populations roughly overlaid the high-confidence multi-omics annotations ( Figures 1B, C ).

Using both the multi-omics and flow cytometric-style gating approaches for annotating the B cells, we were able to determine the concordance between labelling strategies, and thus characterise the level and nature of cellular contamination within standard flow cytometry gating ( Figure 2A , cells from all diseases stats and health). While the >99% of flow cytometric-style gated plasmablasts exhibited a plasmablast profile when using all the multi-omic information (multi-omic-annotation), only 69% of flow cytometric-style gated naïve B cells exhibited a multi-omics naïve cell profile. Similarly, only 5.13% of the unswitched memory B cells, as defined by the multi-omics annotation, were captured within the unswitched memory flow cytometric-style gate. Overall, we show the accuracy of flow cytometric-style gating ranged drastically between 77% to >99% depending on the populations of interest ( Figure 2B ). The same trend was observed when considering only cells from healthy individuals ( Supplementary Figure S2 ).

We next explored the cellular heterogeneity within the naïve and memory flow cytometric-style defined populations. 29.5% of flow cytometric-style gated naïve cells were defined as memory B cells via multi-omics information. These different phenotypes had distinct isotype distributions based on immunoglobulin RNA sequence expression ( Figure 3A ) that were significantly different between multi-omics labelled populations (p-values<0.05) and CD27 expression ( Figure 3B , p-value=2.2e-16)), albeit with low CD27 protein expression ( Figure 1A ). Likewise, 3.0% of the flow cytometric-style gated memory B cells were defined as plasmablasts via multi-omics information. These different phenotypes in the flow cytometric-style gated memory B cells also had distinct isotype distributions based on immunoglobulin RNA sequence expression ( Figure 3C ) that were significantly different between multi-omics labelled populations (p-values<0.05) and CD27 expression ( Figure 3D , p-value=1.7e-5). Indeed, we show that, while CD27 protein and gene expression is significantly correlated, the correlation is poor across B cell subsets ( Supplementary Figure S2C ). Finally, the flow cytometric-style gating of the switched memory B cells had 90% accuracy, however the gating of unswitched memory B cells was lower (77.3%), with the majority of impurities in this gate consisting of switched memory B cells.

Together this demonstrates that flow cytometric-style gating is successful at enriching particular cell groups such as plasmablasts, however classical naïve and memory B cell gates result in heterogenous B cell populations that can be clearly elucidated considering gene expression and VDJ information. These populations are functionally distinct, with different isotype usages, CD27 expression, and B cell repertoire features.

We next explored the phenotypic heterogeneity of B_ND cells (CD19⁺ IgD⁺ IgM⁻ CD27^- CD10^- CD24mid/low CD38mid/low) (4), CD21- Atypical B cells (CD19⁺ CD20⁺ CD10^- CD21^- CD27^-) (7), double negative B cells (DNB) (CD19⁺ CD27^- IgD^-) (1, 13–15), age-associated B cells (CD19⁺ CD21⁻ CD11c⁺) (34), anergic B cells (CD19⁺ CD21^−/low CD38^-) (5), autoreactive IgMlo naïve B cells (CD27^- IgD⁺ IgM^lo) (6) based on FACS gating strategies used in the literature ( Figure 4A ). Using the same gating strategies used in these studies on the CITE-seq values ( Supplementary Figure S3 ), we were able to capture each of these populations in the single-cell multi-omics data across health and disease states. Comparison of these populations with the multi-omics annotations revealed significant heterogeneity between these populations ( Figure 4B ). Indeed, age-associated, anergic B cells and CD21- atypical B cells were slightly enriched for unswitched and switched memory B cells, the double negative (DNB) B cells were enriched for switched memory and plasmablasts, whilst IgMlo naïve B cells were enriched for naïve B cells. The B_ND cells were not enriched for any specific multi-omics- phenotype, suggesting that CD19⁺ IgD⁺ IgM⁻ CD27^- is not a specific gating strategy. Overall, classically-gated anergic, age-associated, autoreactive IgMlo naïve, B_ND, CD21- atypical, and DNB B cells capture highly heterogenous populations representing diverse gene expression and protein expression patters.

To explore this further, we quantified the overlap between classically-gated anergic, age-associated, autoreactive IgMlo naïve, B_ND, CD21- atypical, and DNB B cells using the Jaccard Index, where a higher value indicates higher levels of overlap between two populations (Figure 4C). Whilst some populations were highly distinct with low Jaccard Indices, such as between autoreactive IgMlo naïve B cells and age-associated B cells, we show that there is high overlap between many of the other populations. This is most notable between CD21- atypical B cells (CD19+ CD21− CD11c+) and anergic B cells (CD19+ CD21−/low CD38-) in which the overlap was 0.98 and with the DNB cells (CD19+ CD27- IgD-) in which the overlap was 0.8. Likewise, B_ND B cells (CD19+ IgD+ IgM− CD27-) overlap highly with autoreactive IgMlo naïve B cells (overlap = 0.98). Together, this exemplified that these B cell populations are not mutually exclusive, and that functional studies on these populations are often measuring partially overlapping groups of B cells.

Next, we considered whether the cellular composition of flow cytometric-style gated B cell populations differed between disease states ( Supplementary Table S2 ). Using the flow cytometric-style gating of naïve, memory and atypical B cells, we observed significant differences in the proportions of some multi-omics defined populations between disease states ( Figure 5 ; Supplementary Figure S4 , p-values<0.05). Indeed, we show that the proportion of flow cytometric-style gated naïve B cells that have multi-omics profiles of memory unswitched B cells are significantly variable between patient groups, with COVID-19 mild patients with the lowest level of contamination, and influenza patients with the highest level of contamination from unswitched memory switched B cells ( Figures 5A, B , p-values<0.05). The proportion of flow cytometric-style gated IgD-CD27- B cells comprising multi-omics defined switched memory were significantly associated with disease status (p-values<0.05). Together this demonstrates that the compositions of classical flow cytometric-style defined B cell populations significantly differ between disease state.

To overcome the heterogeneity of flow cytometric-style gating of B cell populations, we performed a data-driven analysis of which markers would most appropriately enrich for purer populations of naïve, unswitched memory, switched memory, and IgD- CD27- B cells, that provide a lower level of transcriptional heterogeneity. Here, we considered only the addition of up to 2 additional markers to reflect the constraints of standard FACS sorting or flow-cytometric experiments. It is expected that increasing the number of cell markers will improve the separation of B cell subsets, as compared to gating with a lower number of markers. To achieve this, we performed an unbiased marker selection using a machine learning approach, named AlliGateR (All marker enrichment for additional flow-cytometric GATEs for purer populations) ( Figure 6A ). The Maximum Mean Discrepancy (MMD), which is a measure of dissimilarity between two probability distributions, was used to identify markers that are more adept at discriminating the true populations from impurities identified from the multi-omics data. Finally, the choice of markers needs to be biologically relevant and reflecting lineage definitions, rather than activation status.

For each flow-cytometric gated population, we identified true positive (TP) cells (those that were correctly identified as defined by multi-omics annotations), and false positive (FP) cells (those that were incorrectly identified by flow-cytometric gating but annotated differently by multi-omics analysis), and for each cell surface protein marker, we trained a sigmoidal support vector machine (SVM), and used its predicted annotation (TP or FP) to determine the sensitivity, specificity and accuracy for additional marker selection ( Supplementary Table S3 ). From this analysis, a combination of CD20, CD21, and CD24 were found to be the best markers for discriminating TP naïve, unswitched memory, and switched memory B cells from impurities (of which CD24 was already included in the classical gating strategy) ( Supplementary Figure S5A ). However, the multi-omics comparison demonstrated the need for non-linear gates through assessing the highest density of purer B cell populations. These markers were supported by the highly significant differences in protein levels between TP and FP populations ( Supplementary Table S4 ). Therefore, additional gating was performed using these markers ( Figures 6B, C ). These additional gates did increase the overall purity of naïve and unswitched and switched memory populations when compared to the original standard gating, but only by between 7.62%, 4.70% and 2.92% respectively. The additional gates for the naïve population (CD20lo/mid CD24lo/mid CD21mid) predominantly reduced the unswitched and switched memory B cell impurities and the plasmablast impurities (by 39.4%, 31.1% and 54.3% respectively, Supplementary Figure S3B ), which would likely significantly impact the functional readouts of any downstream experiments. The majority of the residual impurities were from transitional B cells.

The additional gates for the switched memory (CD24+ CD20lo/mid CD21hi) removed plasmablast impurities. The unswitched memory B cells were divided into three populations, of which two improved the impurity rate (CD24+ CD20hi and CD24+ CD20lo). The additional CD24- gate for the unswitched and switched memory effectively captured the CD27+ plasmablast impurities ( Supplementary Figure S3B ). We note that additional gating did, however, significantly reduce the number of cells captured within each gate, albeit with lower levels of impurities.

We also investigated whether the IgD- CD27- B cells (also termed double negative B cells, DNB, in the literature, Figure 3A ) could be subsetted into more homogeneous groups using these markers. Indeed, separating the IgD- CD27- B cells into (a) CD24- CD20hi CD21-, (b) CD24- CD20lo, and (c) CD24+ CD20lo/mid CD21hi, we were able to enrich for specific multi-omics phenotypes. Indeed, 92.19% of the IgD- CD27- [CD24- CD20lo] gated B cells consisted of switched memory B cells as annotated by multi-omics, whereas 41.12%, 56.57%, and 2.14% of the IgD- CD27- [CD24+ CD20lo/mid CD21hi] gated B cells consisted of unswitched memory, switched memory and naïve B cells, respectively, as annotated by multi-omics. Finally, the IgD- CD27- [CD24- CD20hi CD21-] gated B cells inhibited the highest proportion of plasmablasts (5.20%), as annotated by multi-omics, however, the majority of these cells were switched memory B cells (61.05%). Together, we provide a tool to identify additional protein markers that may be used to provide purer populations by flow cytometry, and may be used more generally for other cell types.

Finally, we assessed whether flow-cytometric gated B cell populations with increased purity would reduce the association with disease status. Interestingly, we show that, although the purities of the naïve and memory populations are increased with the additional gates, we showed that there were more associations between disease status and impurity levels of the flow-cytometric gating ( Supplementary Figure S5C ), particularly naïve and switched memory B cells. Overall, we demonstrate a data-driven multi-omics approach to improving experimental purity of B cell populations, and quantify the increased purity of these extra gating approaches. However, this also provides caution on the implications of enumeration and functional readouts of gating strategies when comparing between diseases.

Data was taken from the COMBAT study (17), which included a multi-omic blood atlas encompassing acute patients with varying COVID-19 severity. This data included single cell gene expression, CITE-seq, BCR and TCR VDJ information on matched cells, in which we performed high-confidence annotations of all immune cells with considering all modalities. The full list of CITE-seq markers used are included in Supplementary Table S1 . This was used as the foundation of this study.

Data pre-processing and annotation using the multi-omics information was performed as described in (17). Briefly, following inspection of the QC metrics, the dataset was filtered to retain cells with ngenes > 300 and pct_mitochondrial < 10%. For the annotation of the immune cell subsets, we used expert immunological knowledge to guide a curated integration of the data from the different modalities (GEX, ADT and VDJ) to identify and label the cell sub-populations present. We first performed separate clustering of gene expression, clustering of surface protein expression, and analyses of T and B cell receptor V(D)J sequences [described fully in (17)]. Cell types and subsets were further refined using information from the repertoire and GEX layers, or in the absence of definitive ADT information were identified by GEX cluster phenotype led by expert understanding of each immune cell subset. Finally, the identified cell types and subsets were further divided by inferred functional state based on targeted assessment of information from all three modalities. For example, cell cycle phase was determined by GEX phenotype, while assignment of B cell maturation status involved use of information from all three modalities (including BCR mutational status). Information from all three modalities was used to identify and exclude doublets from downstream analysis.

To gate cell population in a flow cytometry-like style, we first identified the most commonly used markers (cell surface markers) to identify each targeted population. CITEViz (version 0.1) in R was used to visualise, set thresholds, and gate cells from the original multi-omics dataset based on ADT information. The negative thresholds were based on ADT level densities within populations of cells that are known not to express each marker.

The same methodology was applied to identify and validate markers from the literature. Here we gated B_ND cells (CD19⁺ IgD⁺ IgM⁻ CD27^- CD10^- CD24mid/low CD38mid/low) (4), CD21- Atypical B cells (CD19⁺ CD20⁺ CD10^- CD21^- CD27^-) (7), double negative B cells (DNB) (CD19⁺ CD27^- IgD^-) (1, 13–15), age-associated B cells (CD19⁺ CD21⁻ CD11c⁺) (34), anergic B cells (CD19⁺ CD21^−/low CD38^-) (5), autoreactive IgMlo naïve B cells (CD27^- IgD⁺ IgM^lo) (6) based on FACS gating strategies used in the literature ( Figure 4A ; Supplementary Figure S3 ).

Three methods were used to prioritise additional antibody markers for separating B cell populations. Firstly, we used FindAllMarkers function from Seurat (version 5.0.1) in R to find differentially expressed markers for each FACS-like cluster. Subsequently, we filtered out any markers with an average log2FC < 2 and adjusted p-value > 0.05. Secondly, we used support vector machine model (svm model) to find the precision, accuracy, sensitivity, recall score, F1 score, false positive and false negative values for each marker in each cluster. Finally, the Maximum Mean Discrepancy (MMD), which is a measure of dissimilarity between two probability distributions, was used. The plasmablast populations (CD27+IgM+ PB, CD27+ PBs, IgD-CD27- PBs) were grouped together for the classification, as they showed a strong overlap with multi-omics plasmablast population. The selected markers were then examined for their effectiveness in distinctly distinguishing their respective cluster from other populations using density map. Population gating with new additional markers was performed as described in the above section, shown in Figure 6 .

All analysis were conducted using R version (4.2.3). Jensen-Shannon divergence was calculated to measure the similarity between literature-based gating. ANOVA or t-tests was used to find the significance between two groups or more. NS, not significant; *P < 0.05, **P < 0.01, ***P < 0.001, ****P < 0.0001.