Abstract
Broadening our understanding of the abundance and phenotype of B cell subsets that are induced or perturbed by exogenous Ags will improve the vaccine evaluation process. Mass cytometry (CyTOF) is being used to increase the number of markers that can be investigated in single cells, and therefore characterize cell phenotype at an unprecedented level. We designed a panel of CyTOF Abs to compare the B cell response in cynomolgus macaques at baseline, and 8 and 28 d after the second homologous immunization with modified vaccinia virus Ankara. The spanning-tree progression analysis of density-normalized events (SPADE) algorithm was used to identify clusters of CD20+ B cells. Our data revealed the phenotypic complexity and diversity of circulating B cells at steady-state and significant vaccine-induced changes in the proportions of some B cell clusters. All SPADE clusters, including those altered quantitatively by vaccination, were characterized phenotypically and compared using double hierarchical clustering. Vaccine-altered clusters composed of previously described subsets including CD27hiCD21lo activated memory and CD27+CD21+ resting memory B cells, and subphenotypes with novel patterns of marker coexpression. The expansion, followed by the contraction, of a single memory B cell SPADE cluster was positively correlated with serum anti-vaccine Ab titers. Similar results were generated by a different algorithm, automatic classification of cellular expression by nonlinear stochastic embedding. In conclusion, we present an in-depth characterization of B cell subphenotypes and proportions, before and after vaccination, using a two-step clustering analysis of CyTOF data, which is suitable for longitudinal studies and B cell subsets and biomarkers discovery.
Introduction
The humoral immune response is a major correlate of protection for poxvirus infection (1–8), as well as many infectious agents (reviewed in Ref. 9). Immunization of humans with live smallpox vaccines elicits very long-lived memory B cell and neutralizing Ab (nAb) responses (2, 10) that resulted in global eradication of the disease. In macaques, modified vaccinia virus Ankara (MVA), the current benchmark for vaccinia virus–based vaccines (reviewed in Ref. 11), also elicits robust nAb responses and confers protection against lethal monkeypox viral challenge via various routes (12, 13), yet it has a much improved safety profile as compared with vaccinia virus (14–17). Importantly, MVA is also being developed as a vaccine vector for many infectious and noncommunicable diseases (18–21).
The current model of humoral immunity includes memory B cells that are either short-lived or long-lived and can differentiate into plasmablasts (PBs) upon restimulation (22), and durable serum Ab responses (>1 y) that are maintained by specific PBs or plasma cells (PCs), including in the absence of memory B cells and Ag (23, 24). However, high-parameter single-cell phenotyping of the B cell response at steady-state or after immunological perturbation is limited to a few studies (25–28). Therefore, improving our understanding of the processes by which specific and durable Ab responses are elicited is highly relevant for the development and evaluation of new vaccines that also aim to induce effective humoral immunity.
As the discovery of molecules that are associated with B cell activation and memory differentiation continues to increase, so too has the requirement to assess greater numbers of markers at the single-cell level. Because of fluorescence overlap, conventional flow cytometry has reached its upper limit of ∼18 parameters. New technologies such as spectral flow cytometry (29) and mass cytometry (CyTOF) (30) are now being applied to increase the number of markers that can be investigated simultaneously. Mass cytometry has superseded fluorescent-based flow cytometry via the use of mass spectrometry to detect heavy metals conjugated to specific ligands of cell markers such as Abs, which can in theory allow the detection of up to a hundred parameters at the single-cell level. Today panels with >30 parameters have been used in immunological studies, including investigations of CD8+ T cells (31), CD4+ regulatory T cells (32), CD4+ follicular Th cells (33), B cells (27, 34–36), and vaccines (18, 37). The benefits of larger Ab panels include the discovery of new subsets (31, 35), the identification of coexpressed markers (27), and of signatures for differentiation, activation, and homing status of cells. Novel findings that are generated by high-dimensional data are also expected to facilitate the discovery of biomarkers or cell subsets with predictive capacity in terms of vaccine and immunotherapies efficacy, and unravel phenotypic progressions (27, 31, 38) and new mechanisms of induced immune responses that could guide future pharmaceutical developments (39–44).
In this study, we characterized the B cell immune response in the blood of nonhuman primates immunized twice with MVA, 2 mo apart. Because macaques have a high degree of similarity with humans in terms of markers and cell subsets (45), they represent a desirable model for translational studies, especially those requiring repeated tissue sampling. Importantly, the systemic B cell response is thought to reflect closely the overall humoral immune status of vaccinees (46). Therefore, blood sampling was chosen to preserve tissues associated with B cell priming, memory B cells, and PB/PCs reservoirs (47, 48), as well as increase the adaptability of these investigations to human trials. We used a 26-parameter CyTOF panel to identify MVA-elicited or altered cell phenotypes that could be differentially enriched at days 8 and 28 postvaccine boost (D8PB and D28PB, respectively) as compared with baseline, before any MVA immunization. The CyTOF panel was designed to include B cell–associated markers with the capacity to characterize vaccine-induced populations using high-dimensional data analysis. Whereas previous multidimensional studies have focused on comparing different tissue compartments (27, 33), healthy versus diseased samples (35), or Ag-specific responses (18, 31), we developed a data mining and analysis pipeline to identify and deeply characterize the phenotype of vaccine-induced B cells that changed in frequency over time and those that correlated with the magnitude of the anti-vaccine Ab response.
Materials and Methods
Experimental design
This work is an evaluation of the immune responses of five cynomolgus macaques after s.c. inoculation with MVA, an approved human vaccine. We hypothesized that changes in B cell subsets (subpopulation distribution and/or phenotypic changes) would be present in the periphery at D8PB and D28PB, and that single-cell CyTOF staining followed by high-dimensional data analysis would be able to identify these changes without the requirement for reagents to detect MVA-specific B cells via their BCR. Ab titer and neutralizing activity assays were performed. PBMCs were isolated from blood samples taken at baseline day zero (D0), D8PB, and D28PB to compare the cellular composition before Ag experience and around the time of a robust Ab response. Single-cell analysis was carried out using a 26-parameter mass cytometry Ab panel that had been designed to characterize B cells in detail. To mine the data, we developed an approach using cell clustering and statistical tools, with the aim of identifying proportional changes in B cell subset distribution that occurred after prime-boost vaccination.
MVA inoculation and macaque sampling
Five adult male cynomolgus macaques (Macaca fascicularis) were imported from Mauritius and housed in the animal facilities of the Infectious Disease Models and Innovative Therapies (IDMIT) Infrastructure, Fontenay-aux-Roses, France. Macaques were handled by veterinarians under anesthesia using an i.m. injection of 10 mg/kg ketamine (Rhone-Mérieux, Lyon, France) during the immunizations and blood sampling. On D0 (baseline) and D58, macaques were injected via the s.c. route with 2 × 1 ml inoculum in regions drained by the two inguinal lymph nodes, which contained a total of 4 × 108 PFU recombinant MVA that expresses HIV proteins Gag, Pol, and Nef (MVA HIV-B MVATG17401; Transgene, Illkirch-Graffenstaden, France) (49). Macaque blood was collected into either lithium-heparin tubes (Vacutainer BD) before immunization on D0 and also at D8PB and D28PB for PBMC isolation, into K3-EDTA tubes for complete blood count (Greiner Bio-One, Frickenhausen, Germany), or into serum clot activator tubes (Greiner Bio-One) at the indicated additional time points. Animals were monitored daily by staff for signs of disease, appetite, and lethargy. At each bleed and inoculation, clinical examination was performed. Although it would have been highly relevant to include lymphoid organs and tissues, such as lymph nodes, tonsils, spleen, and bone marrow, in addition to blood in our analysis of B cells, biopsies were not collected to avoid interference with the development and maintenance of immune memory. Nonhuman primates were handled in accordance to national regulations (Commisariat à l'Energie Atomique et aux Energies Alternatives Permit Number A 92-32-02), in compliance with Standards for Human Care and Use of Laboratory of the Office for Laboratory Animal Welfare under Office for Laboratory Animal Welfare Assurance number A5826-01, and the European Directive (2010/63, recommendation N°9). The study was approved by the Ministère de l’Éducation Nationale, de l’Enseignement Supérieur et de la Recherche (France) and the ethical committee Comité d’Éthique en Expérimentation Animale n°44 under the reference 2015031314518254 (APAFIS#319).02.
Sample processing and storage
The overall number of major leukocyte subsets per milliliter of EDTA tube macaque blood was determined using an automated blood counter (HmX; Beckman Coulter). PBMC isolation was performed by overlaying a 1:1 (v/v) mixture of Li-heparin–collected whole blood and Dulbecco’s PBS (Invitrogen) on a mixture of 95% Lymphocyte Separation Media 1077 (GE Healthcare, PAA, Austria) and 5% PBS, and centrifuging at 400 × g for 45 min at room temperature (RT). PBMCs at the Separation Media–plasma interface were collected and washed twice in complete media (complete RPMI), which comprised RPMI 1640 media (Invitrogen) supplemented with 10% heat-inactivated FBS (Lonza) and 100 U/ml penicillin, streptomycin, and neomycin (Invitrogen). Serum was isolated from blood by centrifuging the tubes at 3000 rpm for 10 min and stored at −80°C.
Measurement of total anti-MVA IgG and IgM in serum
In brief, wild-type MVA (B. Verrier, Biologie tissulaire et ingénierie thérapeutique, Institute of Biology and Chemistry of Proteins, Lyon, France) was used to coat 96-well MaxiSorp microplates (Nunc) at 105 PFU/well in coating buffer (200 mM NaHCO3, 80 mM Na2CO3, pH 9.5) overnight at 4°C. Wells were washed five times with wash buffer (PBS, 0.1% Tween 20, 10 mM EDTA) and blocked for 1 h at RT with 3% w/v BSA (Sigma). Plates were washed five times and incubated with 2-fold serial dilutions of macaque serum diluted in 1% w/v BSA in PBS for 2 h at RT, starting at 1:50 for IgG or 1:25 for IgM. Plates were then washed five times and 1:2500 or 1:1000 peroxidase-conjugated goat anti-monkey H+L chain IgG or IgM, respectively (AbD Serotec), in 1% BSA (w/v) PBS was added and incubated for 1 h at RT. Plates were washed five times and 100 μl 3,3′,5,5′-tetramethylbenzidine (Thermo Scientific) was added and incubated for 30 min at RT in the dark. The reaction was stopped by adding 2 N H2SO4. Absorbance was measured at 450 nm using a spectrophotometer (Multiskan FC, Thermo Scientific), and data were analyzed using SoftMax Pro software (version 4.6; Molecular Devices). Ab titers were calculated by extrapolation from the OD as a function of a serum dilution curve and defined as the dilution of the test serum reaching 2 OD of the corresponding preimmune serum tested at 1:50.
Ab neutralization assay
nAb titer was determined using a modified version of a standard plaque inhibition assay. In brief, wild-type MVA (1 PFU/cell) was mixed with an equal volume of 2-fold serial dilutions of serum in assay medium (DMEM, 2% FCS), starting at 1:20. After 60 min of incubation at 37°C, 0.1 ml of the serum–virus mixture was transferred, in duplicate, to a 96-well plate containing subconfluent UMNSAH/DF-1 chicken fibroblasts (ATCC). After 48 h of incubation at 37°C, cell viability was quantified using an MTS/PMS assay (CellTiter 96 AQueous Non-Radioactive Cell Proliferation Assay; Promega). Absorbance was measured at 492 nm using a spectrophotometer (Multiskan FC; Thermo Scientific), and data were analyzed using SoftMax Pro software (version 4.6; Molecular Devices). The sample dilution versus the percentage viability was plotted (four-parameter logistic curve) to calculate a neutralizing concentration, corresponding to the sample dilution resulting in 50% neutralization of virus-mediated cell mortality. Cell viability in uninfected control cells and in infected cells incubated with undiluted vaccinia immune globulin i.v. (human polyclonal anti–vaccinia virus IgG; BEI Resources) was equivalent as expected.
Ab–element conjugation
Pure (carrier-protein–free) mAbs or polyclonal Abs from various manufacturers were coupled to elements using MAXPAR lanthanide labeling kits (Fluidigm, San Francisco, CA), as indicated in the manufacturer’s preload method for 400 μg Ab, with the element–Ab combinations shown in Supplemental Table I. Element–Ab conjugates were adjusted to 1 μg/μl in Ab stabilizer buffer (Candor Bioscience, Wangen, Germany), supplemented with sodium azide (Santa Cruz Biotechnology) to a final concentration of 0.01%, and stored at 4°C in sterile conditions throughout the study. Element–Abs were titrated on cells from healthy macaques to obtain optimal staining concentrations.
CyTOF staining protocol and file preprocessing
Fresh PBMCs (5 × 106/well) were stained with element-tagged Abs. In brief, PBMCs were washed twice in staining buffer (SB; BD Biosciences, Franklin Lakes, NJ) and then incubated for 45 min on ice with the element-tagged surface-stain Abs detailed in Supplemental Table I in a total volume of 50 μl SB. Cells were washed twice with SB and then fixed in 2% paraformaldehyde (Electron Microscopy Sciences) in PBS for 45 min on ice. Cells were washed twice in SB and then three times in permeabilization buffer (eBioscience), and then stained with intracellular Abs (Supplemental Table I) in 50 μl permeabilization buffer on ice for 30 min. Cells were washed twice in SB and resuspended in 0.25 μM iridium DNA-intercalator (Fluidigm) in 200 μl PBS + 2% paraformaldehyde overnight at 4°C, before washing once in SB, once in PBS, and three times in ddH2O, and filtered through a 5 ml Polystyrene round-bottom tube with cell strainer cap (BD Biosciences). Fifty microliters of 4-Element EQ Beads (Fluidigm) was added to each sample, which were acquired as two replicates using an autosampler (Fluidigm) and mass cytometer (CyTOF; Fluidigm) following the standard procedure recommended by the manufacturer (50).
Flow cytometry standard files produced by the CyTOF were normalized with the MatLab Compiler software normalizer (51) using the signal from the 4-Element EQ beads (Fluidigm) as recommended by the software developers. The replicates of each sample were concatenated and modified to remove unwanted or unusable channels using R software (Version 6.2; R Foundation for Statistical Computing, Austria).
Spanning-tree progression analysis of density-normalized event
The spanning-tree progression analysis of density-normalized events (SPADE) algorithm was used to identify populations of cells having similar marker intensities in cytometry datasets. SPADE was performed using the publicly available R/Bioconductor package (52). To reduce bias during the SPADE down-sampling step caused by large variances between the numbers of CD20+ events in each flow cytometry standard file, which ranged from 3,750 to 37,983 CyTOF acquired B cells for macaque BC641 at D0 baseline and macaque BD620 at D28PB, respectively (Supplemental Fig. 1F), we first randomly selected 3000 cells uniformly for each sample and then performed the SPADE density-based down-sampling using a target of 750 cells from each pool of 3000 randomly preselected cells (i.e., 25% down-sampling). Individual SPADE trees were then generated using all CD20+ cells from each sample (up-sampling). The following SPADE settings were as follows: 10 clustering markers (CD20, CD21, CD22, CD23, CD27, CD40, CD80, HLA-DR, sIgG, and sIgM), 70 clusters, and an outlier density parameter of 0.01.
Identification of abundant clusters and SPADE differentially enriched clusters
Abundant clusters and SPADE differentially enriched clusters (DECs) were identified based on the percentage of cells in individual SPADE clusters, relative to total CD20+ B cells. For each SPADE cluster at each time point, the mean percentage of CD20+ cells and associated SD value for five macaques were computed. A one-sample Student t test (p < 0.01) was used to determine whether clusters were abundant, and whether the proportion of each cluster was similar between the five macaques. In addition, to identify SPADE DECs, we compared SPADE cluster sizes between each pair of the three time points using a paired Student t test corrected for multiple comparisons using local false discovery rate (q < 0.05).
Multidimensional scaling
Multidimensional scaling (MDS) methods aim to represent the similarities and differences among high-dimensional objects in a two- or three-dimensional space for visualization purposes (53). In MDS representations, pairwise distances between the data points are proportional to the Euclidean distances between the samples. In this study, we used MDS representations of: 1) the similarities between macaques in terms of their relative cell concentrations of different clusters, or 2) the similarities between SPADE clusters in terms of the means of their markers intensities. The Kruskal stress (KS) criterion (53) was calculated to quantify the fraction of the information lost during the dimensionality reduction procedure, as a measure of MDS quality.
Automatic classification of cellular expression by nonlinear stochastic embedding and identification of automatic classification of cellular expression by nonlinear stochastic embedding DECs
Automatic classification of cellular expression by nonlinear stochastic embedding (ACCENSE) was used to reduce the dimensionality of the mass cytometry data and to automatically identify cell clusters according to their distribution of markers (54). The 3000 CD20+ cells that were randomly selected from each sample during the initial step of our SPADE analysis were clustered using ACCENSE Version 0.4.3. The use of Barnes-Hut stochastic neighbor embedding (SNE), the k-means clustering method, the same 10 clustering markers used for the SPADE analysis, a perplexity value of 30, and a p value of 0.0001, resulted in the designation of 80 ACCENSE clusters. There was no up-sampling with all CD20+ cells from each sample, in contrast with SPADE. ACCENSE DECs were defined in a similar way to SPADE DECs by using a paired Student t test corrected for multiple comparisons using local false discovery rate and a q value <0.05 to compare the frequency of cells in ACCENSE clusters between different pairs of time points.
Phenotypic categorization of SPADE and ACCENSE clusters and manually gated cells
The phenotypic categorization of SPADE and ACCENSE clusters was calculated using the 5th and 95th percentile range for each marker intensity and then dividing this range into five uniform categories. The categorization was computed based on the means of the individual SPADE expression medians for each marker. The categories represented negative, low, medium, high, and very high marker expression using a color scale from white to dark red, which was used to produce a heat map. Double hierarchical clustering of the phenotypical heat map was performed using the Euclidean metric and the complete linkage method. All markers except for CD3, CD14, and Iridium DNA intercalators were included.
Phenotypic categorization using manual gating of CD20+ cells involved five successive bivariate gates that were defined using 10 markers: sIgM, sIgG, CD21, CD27, CD22, CD40, CD62L, CD45RA, Bcl-2, and HLA-DR, most of which show a large intensity range, and thus facilitating the gating of high, intermediate, low, or negative cells for a given marker. The positions of the gates were chosen based on the marker expression of SPADE DECs 9 and 16 from group L (Figs. 4, 5, 6A). For example, the manual gate was set on cells that were slightly positive for IgM and not IgM− because the majority of cells from DECs 9 and 16 showed this pattern of IgM expression (see Fig. 4A). The same gate drawing approach was applied to the other nine markers selected for the manual gating analysis.
Ab correlation
Associations among the frequency of SPADE clusters, ACCENSE clusters, or manually gated cells, expressed as a percentage of CD20+ B cells, and other biological variables, such as anti-MVA Ab and neutralization Ab titer, were quantified using a Pearson coefficient of correlation. A coefficient of correlation (R value) >0.80 or <−0.80 was considered significant.
Associations between SPADE and ACCENSE clusters and with manually gated cells
The similarity between SPADE and ACCENSE clusters and manually gated cells phenotypes were computed using a Spearman coefficient of correlation. Correlations were computed based on the 10 clustering markers mean signal intensity (MSI) expressed as categories as presented in the heat maps. Associations with an R coefficient >0.90 and p < 0.01 were considered significant.
Results
Multiparametric SPADE clustering to group B cells based on phenotype
Five cynomolgus macaques were primed and then boosted 2 mo later with MVA via the s.c. route (Fig. 1A) with the aim of inducing a robust B cell response that could be characterized in-depth using multiparametric clustering algorithms and used to generate hypotheses about changes in the peripheral B cell compartment after vaccination. The overall experimental approach (Fig. 1B) involved collection of blood samples before immunization (D0) and at two subsequent time points, D8PB and D28PB, followed by staining of PBMCs with a B cell–focused Ab panel (Fig. 1C, Supplemental Table I) and CyTOF acquisition. Two injections of MVA were used to increase the magnitude of the B cell response, and the two time points postboost were included in the study to analyze response kinetics and contrast an early and a later secondary response after vaccination. The CyTOF Ab panel was aimed at identifying 25 markers, including those associated with BCR signaling, activation, survival, exhaustion, homing, and interaction with T cells. In addition, two non–B cell exclusion markers, CD3 and CD14, were included to facilitate the gating of CD20+ cells (see gating strategy Supplemental Fig. 1A). Partial downregulation of CD20 occurs during differentiation into Ab-secreting cells, and human PB/PCs are well characterized as CD3−CD19+CD20−/low cells. CD19 expression could not be used to track the appearance of PB/PCs in our study because none of the anti-human CD19 clones we tested resulted in a positive signal, including clone J4.119 from Beckman Coulter (data not shown), which is the only reported clone with the ability to cross-react with cynomolgus macaque cells in flow cytometry, albeit with a low efficiency (http://www.nhpreagents.org/). In addition, in contrast with previous studies (26, 55), it was recently shown that rhesus macaque PBs are primarily CD19−. PB/PCs were defined as CD3−CD16−CD20−/intHLA-DR+CD14−CD11c−CD123−CD80+ cells (56). The cynomolgus macaque (Macaca fascicularis) model presented in this study could include phenotypes of PB/PCs that differ from other primates, and we expect that our Ab panel may identify a portion of PBs when gating HLA-DR+ CD20low cells. The specificity of the Ab staining for each marker was confirmed by comparing the signal on different PBMC subsets (Supplemental Fig. 1B), after in vitro stimulation of PBMCs with a polyclonal stimulus (Supplemental Fig. 1C), by analyzing B cells from lymphoid organs (Supplemental Fig. 1D) or by documenting interindividual differences in marker expression profiles (Supplemental Fig. 1E).
Analysis of the data consisted of: 1) multiparametric clustering of gated CD20+ cells using SPADE (57); 2) statistical identification of SPADE clusters that changed in size after vaccination; 3) characterization of the phenotype of these clusters using unsupervised double hierarchical clustering based on SPADE clusters and markers (58), which grouped the SPADE clusters to elucidate molecular signatures, coexpression trends, and predict trajectories of phenotypic development; 4) correlation between SPADE clusters and Ab responses; and 5) hypotheses generation (Fig. 1A). SPADE is a clustering method developed to analyze and visualize high-dimensional cytometry data (27, 57, 59–64). It includes an initial density-based down-sampling step performed to preserve cells of rare phenotypes, and it organizes cells into hierarchies of related phenotypes. The user-defined steps of the analysis pipeline involved gating CD20+ cells (Supplemental Fig. 1A), selection of the SPADE clustering markers (CD20, CD21, CD22, CD23, CD27, CD40, CD80, sIgG, sIgM, and HLA-DR), choice of the desired number of clusters (n = 70 clusters), and the level of down-sampling (see 2Materials and Methods and Supplemental Fig. 1F). SPADE analysis was performed using CD20+ events, rather than all cells in a PBMC sample, because the Ab panel was focused on B cell markers, and the exclusion of non-CD20+ events was expected to improve phenotypic separation of B cells. The choice of the 10 clustering markers was based on their wide use to delineate B cell subsets in humans and macaques in previous publications (26, 46, 55, 65, 66).
After the SPADE analysis, the sizes of the 70 clusters were calculated as a percentage of the parent population (CD20+ B cells), and the average cluster size for the combined five macaques at each time point is represented in Fig. 2A and listed in Supplemental Table II. The 70 clusters ranged in size from ∼0.1 to 5% of total CD20+ B cells. We found that 26, 68, and 24 clusters of 70 were statistically present in the group of 5 macaques on D0, D8PB, and D28PB, respectively (one-sample t test, p < 0.01; Fig. 2A, Supplemental Table II). This analysis highlighted an intermacaque variation in cluster size or the presence of empty clusters at baseline and D28PB that were transiently reduced early after immunization, at D8PB. In addition, it confirmed that 70 clusters were an appropriate choice for our dataset because for at least one time point, at D8PB, most of the clusters were homogeneously filled with cells in all five animals. Because absolute B cell counts per microliter whole blood did not undergo obvious changes for individual macaques between the three different time points (Supplemental Fig. 1G), these values were not factored into the subsequent analysis.
Identification of SPADE DECs
Cluster sizes within the combined SPADE trees at the different time points strongly suggested that vaccination caused changes in the proportion of cells in various clusters over time (Fig. 2A). To investigate this phenomenon, we performed paired t tests corrected for multiple comparisons to identify DECs, which we defined as clusters that had undergone either a significant expansion or a significant contraction relative to the total pool of CD20+ B cells, for any comparison of the three time points (q < 0.05). The majority of clusters were not impacted by immunization in terms of their size. However, for 14 of the 70 SPADE clusters (20%), a significant increase or reduction in size was observed after MVA prime boost (Fig. 2B, Supplemental Table II). It included 8 DECs increasing in size (with a fold change between 2.0 and 7.8 for DEC27 and DEC9, respectively) and 5 DECs decreasing in size (with a fold change between −5.0 and −1.4 for DEC9 and DEC31, respectively). A single cluster, DEC9, underwent both a significant expansion followed by a contraction in size over the course of the study.
When DECs were grouped according to their kinetic profile, based on the timing and directionality of population enrichment, a trend was observed. Nearly every case of significant expansion occurred between the baseline and either D8PB (DEC9, 16, 25, 27, 48, 62, and 67) or D28PB (DEC5), which was in contrast with reductions in cluster size that mainly occurred between D8PB and D28PB (DEC8, 9, 23, 50, and 66). Interestingly, vaccination exerted long-term effects on the proportion of selected B cell clusters, as evidenced by the increase in size of DEC5 and decrease in size of DEC8 and DEC10 that persisted for 1 mo after the second immunization.
MDS is used to represent the similarities between multiparametric data in two dimensions, allowing a global overview of complex datasets (53). An MDS based on the percentage of CD20+ parent values was computed using either all 70 SPADE clusters or only the 14 SPADE DECs (Fig. 3A). The MDS illustrated that when the proportion of cells within all 70 SPADE clusters are represented as a single point for each macaque, the three time points overlap in multidimensional space, whereas the dataset can be clearly separated by the three time points when only the 14 DECs are used. Strikingly, the intermacaque distance at D8PB was lower compared with other time points, again suggesting a more homogeneous B cell response (quantitatively) between animals early after immunization. Each MDS representation scored similarly low for KS, indicating that a high percentage of information was retained during the compression of the data and validating the comparison of the two MDS computations. The MDS results confirmed that an analytical approach based on DEC calculations was an effective means of highlighting major differences in the data in terms of cluster size over time.
The expansion and contraction of SPADE DECs are displayed in further detail in Fig. 3B, which shows DEC size for each macaque at each time point, as well as the average for five macaques, revealing a range of different expansion and contraction kinetic profiles over time. Of particular interest were the DECs that originally represented an insignificant proportion of total B cells (one-sample t test, p > 0.01; Supplemental Table II) at the baseline, such as DEC9 and DEC16, and significantly expanded in proportion after vaccination when comparing D0 with D8PB (paired t test, q < 0.05).
In-depth phenotypic characterization of SPADE DECs
After establishing the presence of clusters that expanded or contracted between pairs of time points, we performed the analysis to characterize DEC phenotype. The SPADE clustering process placed phenotypically similar cells into proximal or identical clusters, as evidenced, respectively, by SPADE trees that were colored based on marker intensity (Supplemental Fig. 2) or by bivariate dot plots where cells from individual DECs grouped relatively closely on the dot plot in comparison with total CD20+ cells, which displayed a larger range of marker intensity (Fig. 4). By using parallel coordinates to plot the mean intensity for each marker, for each macaque and time point, a phenotype signature relative to parent CD20+ cells was visualized for each DEC (Fig. 5). Taken together, these initial phenotypic analyses (Figs. 4, 5, Supplemental Fig. 2) demonstrated that the SPADE clustering was effective and that the 14 vaccine-altered DECs showed a wide variety of phenotypes and could be distinguished from each other based on subtle or obvious patterns of marker expression, including via the intensity of markers that are commonly used to delineate major B cell subsets such as sIgG, sIgM, CD21, and CD27 (25, 26, 46, 67).
To further stratify the SPADE clusters into phenotypic groups and compare DECs not only with each other but also with non-DECs, we performed double hierarchical clustering at the marker and SPADE cluster levels. The resulting heat map and dendrograms were used to define 15 phenotypic groups of clusters (from groups A to O). Each phenotypic group contained between 1 and 11 clusters, with 9 groups containing 1–3 DECs (groups A, B, C, E, H, J, K, L, and N) (Fig. 6A). This analysis highlighted the phenotypic complexity and diversity of circulating B cells at steady-state, as well as the multifaceted impact of vaccination on B cells when DEC phenotype and kinetic profile, including directionality, were taken into account (Fig. 6A). It appeared that DECs that expanded between D0 and D8PB (DEC 27, DEC9 and DEC16, DEC25, DEC48, DEC62 and DEC67) or contracted between D8PB and D28PB (DEC23and DEC8, DEC66, DEC50, and DEC9) belonged to distinct phenotypic groups (groups L, N, J, E, and C, and groups A, B, K, and L, respectively).
Some DECs that were part of to the same phenotypic group also showed a similar kinetic profile. Three DECs that increased in size between D0 and D8PB (DEC27, DEC9, and DEC16) were found in group L, whereas the remaining DECs with a similar kinetic profile (DEC25, DEC48, DEC62, DEC67) belonged to groups N, J, E, and C, respectively. Clusters in group L were sIgG+, sIgMlow, CD21hi, CD27lo, and CD195med, CD62Lmed, CD40hi, CD45RAhi, CD197med, and Bcl-2hi. Previously reported memory B cell subsets definitions based on sIg, CD21, and CD27 (25, 26, 67) suggest that group L comprised resting memory B cells (CD21+, CD27+). The three resting memory–like DECs could be differentiated from each other based on subtle differences in the intensity of clustering markers sIgG and also CD40 (Figs. 4A, 4C, 6A). They were generally similar to other clusters in group L, with a few notable features such as increased CD21 expression. Similarly, two DECs (DEC23 and DEC8), which both significantly decreased in size between D8PB and D28PB, were part of the same hierarchical clustering group (group A), whereas other DECs with the same kinetic profile (DEC66 and DEC50) were members of groups B and K, respectively. Group A was characterized as sIgGhi, sIgM−, CD21−, CD27+, and CD23−, CD195−, CD62L−, CD38lo, CD95med, and CD20high, suggesting they were activated memory–like B cells based on their CD21 and CD27 profile.
However, not all DECs within the same group displayed similar kinetics profiles. For example, in group J, DEC48 increased in size by D8PB, whereas DEC52 and DEC10 underwent expansion or contraction at other time points. DEC52 and DEC10 were unique among the vaccine-induced clusters in terms of sIg expression, because they were characterized by a clear sIgM-only phenotype and were otherwise phenotypically identical to each other except for slight differences in surface marker expression such as CD20, CD21, and CD22. Despite the absence of Ab directed against CD10 in our panel (26), the phenotypes of DEC52 and DEC10 cells were aligned with some subsets of peripheral immature B cells, which have been described as CD38hi, CD21lo, CD23lo, CD27−, CD62L−, and sIgMhi (46, 68). DEC52 also displayed a singular kinetic profile with enrichment occurring between D8 and D28PB.
In addition, some DECs had similar kinetic profiles and phenotypes, although they belonged to distinct groups, for example, DEC66 from group B and DECs 23 and 8 from group A. Both groups are likely to comprise activated memory–like cells (CD21−, CD27+) that are CD20hi, Bcl-2hi, CD62L−, CD195−, and CD23−. However, group B had a markedly higher expression of sIgM compared with group A, and among the highest levels of CD27, CD95, and CD22 of all B cell clusters.
In contrast, some DECs from distinct groups had similar kinetic profiles but distinct phenotypes. For instance, as DEC23, DEC8, and DEC66 (activated memory–like cells), DEC50 reduced in size between D0 and D28PB. However, it belonged to group K, which was mixed for sIgM expression, CD23−, CD195−, and shared many elements of the expression profile described for tissue-like memory B cells, that is, CD21−/lo, CD27−/lo, CD62Llo, and CD20hi (67). Notably, both DEC50 and DEC66 were low or negative for both sIgG and sIgM. Likewise, although DEC25 followed the same kinetics as DEC27, DEC9, and DEC16, which were all part of cluster group L and were reminiscent of resting memory B cells, it was found in cluster group N and contained cells that were intermediate memory–like (CD21hi, CD27−, IgG+) (69).
Finally, DEC62 and DEC67, the remaining DECs that expanded between D0 and D8PB, were found in groups E and C, which were composed of only one or two clusters, respectively. These DECs shared a similar expression profile for few markers such as sIgGhi and sIgM−, but differed greatly for other markers such as CD27, HLA-DR, CD23, CD40, and CD279. Notably, DEC62 (group E) cells showed no to very low expression of CD20, CD22, CD23, CD40, HLA-DR, and sIgM, whereas clusters in group C, including DEC67, contained cells with very high levels of CD23.
Hierarchical clustering performed on the relative MSI of each marker revealed several CD20+ B cell coexpression patterns. For example, the relative expression of CD38 was often mirrored with a similar relative expression level of CD45RA, and hence the hierarchical clustering placed these markers adjacent to one another (Fig. 6A). Similarities were observed for the relative MSI of some other adjacently placed pairs of markers such as CD195 and CD62L, CD20 and CD22, and CD80 and CD95. In contrast, markers placed distally on the heat map were generally characterized by disparate relative MSI values within an individual cluster.
Next, MDS was used in this instance to represent the phenotypic similarities between each cluster, based on the MSI of only the 10 clustering markers (Fig. 6B). Much like the SPADE tree, clusters with similar phenotypic profiles were placed in proximity. However, MDS allowed additional information to be extracted. The unique phenotype of DEC62 became more evident using the MDS representation, and the distance between DEC66 and DEC50 was smaller using MDS compared with using a spanning-tree organization of the clusters (Fig. 2). These observations demonstrate the potential of MDS to represent more accurately clustered data in terms of phenotype similarity, because minimum spanning tree diagrams can occasionally segregate phenotypically related clusters to opposite sides of the tree.
Serum Ab response and correlation with SPADE DECs
Serum Ab assays confirmed the presence of total vaccine-specific IgG at D8PB and D28PB (Fig. 7A) with some detectable neutralizing activity (Fig. 7B). Although low anti-MVA IgG titers were induced after a single immunization, the boost greatly enhanced the humoral response, which was the main rationale for selecting postboost time points for CyTOF evaluation. The peak of the serum Ab response occurred before D28PB; therefore, it is reasonable to assume that the contraction of the peripheral memory B cell response had commenced before this time point, as described in other reports of prime-boost vaccination of cynomolgus macaques (70). The expansion and contraction of cells in the peripheral CD20+ B compartment from a single SPADE DEC, DEC9, correlated with the total anti-MVA serum IgG response (R = 0.85, p = 0.00007) (Fig. 7C, 7D) and, to a lesser extent, with nAb titer (R = 0.79, p = 0.0004) (Fig. 7C). The size changes of DEC16 and DEC27 that belonged to clusters group L together with DEC9 also tended to correlate with the total anti-MVA IgG response (R = 0.67, p = 0.006 and R = 0.75, p = 0.001, respectively), although not as strongly as DEC9. There was also a trend for an inverse relationship between the size of either DEC10 or DEC52, both from clusters group J, with anti-MVA Ab (R = −0.64, p = 0.01 and R = −0.60, p = 0.01, respectively) (Fig. 7C).
Identification of B cell subphenotypes impacted significantly in frequency after vaccination using a different algorithm and manual gating
ACCENSE was used to analyze B cell responses in addition to SPADE. ACCENSE identifies density-based cell clusters within multidimensional data after a nonlinear dimensionality reduction (Barnes-Hut SNE) that retains single-cell resolution of the data (54, 71). Among the 80 ACCENSE clusters (Fig. 8A), 11 underwent a significant increase or decrease over time and were classified as ACCENSE DECs (q < 0.05) (Fig. 8A, 8B). t-SNE maps from each individual at the three time points were useful to visualize and realize the heterogeneity of cluster frequency between animals and time points and to intuitively identify likely DECs (Fig. 8A), before the actual statistical analysis (Fig. 8B). Similarly to the SPADE analysis (Fig. 3B), all ACCENSE DECs between D0 and D8PB expanded (ACCENSE DEC18, 20, 28, 51, 73), whereas most ACCENSE DECs between D8PB and D28PB contracted (ACCENSE DEC4, 20, 56). ACCENSE DEC20 showed the same expansion and contraction profile as SPADE DEC9. A heat map representing the double hierarchical clustering of ACCENSE clusters and markers was generated and allowed the visualization of complex and heterogeneous B cells subphenotypes, as with SPADE (Fig. 8C). Remarkably, the heat map–associated dendrogram of Ab markers after SPADE and ACCENSE clustering showed a similar structure (Figs. 6A, 8C), with sets of markers being grouped together (e.g., CD21, CD38, CD40, CD45RA) by hierarchical clustering regardless of the initial type of clustering used.
To assess whether ACCENSE and SPADE clustered cells into some common subphenotypes, we performed a correlation analysis that was illustrated using a circular graph (Fig. 9A). Of the 70 SPADE clusters, 45 correlated significantly for their phenotype with 53 of the 80 ACCENSE clusters (R > 0.9), with a total of 132 correlations between individual SPADE clusters and one or more ACCENSE clusters (Fig. 9A). The similarity between the two clustering analyses was even more obvious at the level of cluster groups defined by dendrograms on the respective heat maps (Figs. 6A, 8C). Importantly, the phenotype of 10 of the 14 SPADE DECs (SPADE DECs 8 and 23; 66; 5; 48 and 52; 50; 27, 9 and 16) correlated with 5 of the 11 ACCENSE DECs (ACCENSE DECs 4 and 56; 63 and 78; 20). Notably, ACCENSE DEC20, which showed the same kinetic profile as SPADE DEC9 (R = 0.91, p = 5.78E-10) (Figs. 3B, 8B), also tended to correlate with anti-MVA IgG titers (R = 0.75, p = 0.0014) (Fig. 9B) and was phenotypically highly similar to SPADE DEC9 (R = 0.93, p = 9.93E-05), as demonstrated in the selected dot plots (Figs. 4, 9C). Thus, similar results were generated using two distinct algorithms, SPADE and ACCENSE.
Finally, the computational analysis was confirmed by manual gating. Using the individual marker expression profile of cells from SPADE-clustered cells in group L, which was composed of SPADE DEC27, DEC9, cluster 35, DEC16 and cluster 17, we designed a gating strategy of five successive bivariate gates based on sIgM, sIgG, CD21, CD27, CD22, CD40, CD62L, CD45RA, Bcl-2, and HLA-DR expression (Fig. 10A). These manually gated SPADE clusters group L–like B cells showed the same kinetic profile as combined SPADE DEC9, DEC16, and DEC27, with a significant expansion between D0 and D8PB (Fig. 10B) and a correlation with anti-MVA IgG (R = 0.86, p < 0.0001) (Fig. 10C). The phenotypic features of the manually gated cells, represented on a radar chart with SPADE DEC9 and ACCENSE DEC20 cells (Fig. 10D), correlated with SPADE cluster 35 (R = 0.94, p < 0.0001), DEC16 (R = 0.97, p < 0.0001) and cluster 17 (R = 0.92, p < 0.0001), which belonged to SPADE cluster group L.
Discussion
In this study, mass cytometry and a two-step clustering (SPADE, followed by unsupervised hierarchical clustering of SPADE clusters with markers) were used for the first time, to our knowledge, to provide a comprehensive phenotypic characterization of vaccine-altered circulating B cells at multiple time points. Cynomolgus macaque blood B cells were separated into 70 SPADE clusters based on their expression profile of 10 markers. A phenotypically complex and diverse circulating B cell compartment was evidenced at steady-state. The changing proportion of cells in each cluster at various time points was the primary measure of vaccine immunogenicity. This approach allowed the identification of 14 clusters that significantly changed in size (percentage of parent CD20+ B cells) between two or more time points, which we defined as DECs. The identification of DECs was followed by in-depth characterization based on the expression of all 23 B cell–related markers using double hierarchical clustering of SPADE clusters and markers. The 14 SPADE DECs displayed diverse marker expression patterns, including multiple DECs with known memory B cell–like phenotypes such as activated (CD21−CD27+) or resting memory (CD21+CD27+) class-switched (IgG+) B cells, as well as novel B cell subphenotypes. Stratifying the clusters based on their kinetic profile or by using hierarchical clustering of the SPADE clusters and markers revealed marker coexpression patterns, and differences between individual clusters or between groups of clusters.
We compared two widely used and computationally different algorithms that can be applied to mass cytometry data (71). Concordant results were obtained using SPADE, a hierarchical agglomerative clustering analysis followed by a minimum spanning-tree projection, and using ACCENSE, a nonlinear dimensionality reduction analysis followed by clustering based on local maxima within the t-SNE map, both followed by double-hierarchical clustering, in terms of B cell subphenotypes that changed significantly in frequency after vaccination. SPADE analysis was designed to favor the identification of rare cell populations via density-based downsampling step before clustering, including from large datasets. In contrast, the current version of ACCENSE is designed for smaller datasets, and we therefore limited our ACCENSE analysis to the 3000 CD20+ cells that were randomly selected from each sample during the initial step of our SPADE analysis. In addition, the principal findings were confirmed using manual gating analysis, which was only possible after algorithmic identification of cluster phenotypes.
Mass cytometry is still a relatively new technology, and it remains a challenge to ensure high data quality and powerful unsupervised computational analysis, as recently highlighted (71, 72). Progress has been made and normalization beads are now used to overcome variation in mass cytometry instrument sensitivity, for instance. One of the main pitfalls resides in controlling for possible Ab batch effects. Batch effects are mainly related to Ab stability over time, variations in sample processing, staining efficacy, and Ab preparations, which are all critical issues for longitudinal studies that rely on the staining of fresh cells, such as this study. We believe that the possible batch effects were minimized, and thus that the significant impact of vaccination on the frequency of blood B cell clusters over time we report in this study is true and not artifactual. To our knowledge, there is no report in the literature about metal isotope-labeled Abs degradation over time, and the MSI for Ab batches that were used to stain cells at each of the three time points did not change over time, despite the 12-wk interval between baseline and 28 d postboost (data not shown). Blood sampling and processing to isolate PBMC and staining protocols were strictly adhered to throughout the study and measures were put into place to control factors that are well-known to affect staining efficacy, such as temperature and freshness of the fixation and permeabilization reagents. Finally, the same Ab batches were used for all three time points for most, but not all, Abs, including 6 of the 10 Abs that were used to cluster the data with SPADE or ACCENSE (Supplemental Table I). To ensure that different batches of the same Ab produced identical staining, we performed titrations of metal-labeled Ab batches and used these data to finely tune the amount of Ab that we used per well for the subsequent time points. Ideally, the same batch of Abs should be used for the entire longitudinal follow-up and a standardized control sample stained together with the tested samples should be included at each time point. A protocol where storage of stained cells is possible would also be preferable because this would allow CyTOF acquisition of all samples on the same date.
The methodology of vaccine evaluation presented in this analysis did not include the use of reagents to detect MVA-specific B cells, and instead focused on identifying phenotypes of cells that underwent significant proportional changes (DECs) in the B cell compartment. Therefore, we were unable to conclude whether the DECs contained MVA-specific cells, cells of other specificities, or a mixture of both. We hypothesize that at least some of the DEC cells were vaccine-specific, primarily because a robust anti-MVA serum IgG response was detected, and also because cells in some of vaccine-expanded clusters expressed CD80, a marker expressed after interactions between B and T cells and BCR ligation (55). It has been reported that a tetanus toxoid (TT) boost in humans increased TT-specific but also TT-nonspecific PCs of distinct phenotypes in the blood. TT-specific PCs represented precursors of bone marrow PCs, whereas TT-nonspecific PCs were composed of older PCs displaced from bone marrow (73, 74). Whether the same competition for survival niches occurs for memory B cells needs to be addressed. A future challenge will be to evaluate whether it is possible to use cluster characteristics to predict vaccine specificity or activation history of cells, because Ag-independent vaccine evaluation represents a method to investigate therapies for which the protective Ags have not yet been identified. The successful correlation of clusters of a particular phenotype with desirable clinical outcomes, such as long-lived protective immune responses after vaccination, would reduce drug development time and identify nonresponders or therapeutic formulations that require improvement. Research approaches that identify DECs, or other multidimensional analyses that lack specific-Ag detection reagents, are also suited to the evaluation of adjuvants (37) or vaccine-primed T cell responses.
Memory B cells are geared toward the recognition of previously encountered Ag via their sIg receptors and recall responses, and are not specialized for Ab secretion–like PB/PCs, which have downregulated or lost the expression of many surfaces molecules such as Ig, CD20, and HLA-DR (46, 75). Therefore, correlations between the proportion of memory B cells and serum Ab titer observed in this study and others (1, 10, 47) are indirect and suggest that circulating Abs arose from the differentiation of Ag-specific memory B cells in short-lived PCs after the boost. To further elucidate germinal center cell commitment to memory B cells or PB/PCs, future studies should incorporate markers that can characterize additional B cell subsets through their intracellular Ig levels, capacity to secrete Ab, or proliferate (56), as well as analyze additional tissues where most Ab-secreting cells preferentially reside (27, 47). In addition, reagents to detect specific BCRs as described for rotavirus, influenza, papillomavirus, and other Ags (27, 76–79) could be used to reveal SPADE clusters that are enriched with vaccine-specific cells.
This study has highlighted the phenotypic heterogeneity of circulating B cells by subdividing CD20+ cells into clusters, some of which were statistically enriched or impoverished after vaccination. The comparison of activated memory phenotype SPADE groups A and B with resting memory phenotype SPADE group L showed lower levels of expression of CD38, CD40, CD45RA, CD195, CD197, and higher levels of expression of CD95, CD20, and CD22 by cells from the activated memory phenotype groups. Biologically, this implies differences in the ability of these two memory subsets to home to secondary lymphoid organs via CD197 (reviewed in Ref. 80) and to survive to CD95 and CD40 triggering (81). SPADE DEC25 was CD21hi, CD27−, a profile normally defined as naive, but also IgGmed and IgMmed, indicating that class switching had occurred. Thus, DEC25 is more likely to be a CD27− IgG+ intermediate memory population described by Kardava et al. (69) or Fecteau et al. (82). The phenotypic signature of DEC25 was very similar to that of resting memory cells; however, the MSI of CD27 was below the value classified as CD27+, which demonstrates the utility of multiparametric clustering to categorize cells with intermediate phenotypes.
Modification of leukocyte cell subset distribution has been reported as a measure of disease progression. Skewing of some B cell subsets has been associated with repeated exposures to pathogens and chronic infections, including malaria and HIV. A CD21− CD27− atypical malaria-associated memory B cell subset is found in greater frequencies in patients who live in endemic areas (83, 84). In addition, many changes in the B cell compartment have been reported during HIV infection (reviewed in Refs. 85, 86), namely, a reduction in resting memory and an increase in activated and tissue-like exhausted memory B cells subsets. A recent study has also demonstrated that vaccination against influenza A induces IgG-secreting cells in parallel with increases in CD21− CD27+ activated memory cells, as well as a reduction in tissue-like CD21− CD27− memory cells that were exhausted (87). In contrast, our data demonstrate that vaccine-induced CD21− CD27+ IgG+ clusters (SPADE DEC23, DEC8, DEC66) were not enriched between the baseline and D8PB. However, these activated memory–like B cell clusters significantly contracted between D8PB and D28, suggesting either migration out of the blood or, alternatively, displacement via an increase in other clusters in the weeks after an MVA boost. These findings also indicate that DECs with the strongest correlation with Ab titer were all CD27+, a marker commonly used to define B cell subsets and functionality because of its upregulation at the end of the GC reaction (88, 89), and its association with an increased capacity to differentiate into Ab-secreting cells compared with CD27− cells (82). In reality, the use of one marker is unlikely to capture the spectrum of functionalities seen in biological systems. Instead, combinatorial molecular signatures are being used increasingly to assess the quality of the response in terms of clinically relevant lymphocyte subpopulations after vaccination or infection (41, 42). These findings reiterate the potential use of marker signatures obtained using single-cell cytometry, particularly high-parameter mass cytometry (37, 90) or genomic analysis (91), to assign multiple functional capabilities to lymphocytes.
SPADE and hierarchical clustering of SPADE clusters and markers are proving to be useful analytical approaches for identifying correlations of marker expression (27, 92), and which was confirmed with ACCENSE. As previously reported (25), CD21, a complement receptor involved in B cell Ag processing (93), showed patterns of coexpression with CD62L. This coexpression was observed in this study, for example, higher expression of CD21 and CD62L in the three resting memory–like clusters (SPADE DEC9, DEC16, and DEC27; Fig. 6A), in contrast with low or absent expression of both markers in the group A activated memory–like DECs. CD62L has previously been shown to be expressed at higher levels on CD21+CD27+ B cells in peripheral blood compared with CD27 single-positive cells, and drops in intensity on subsets in the mucosa (25), suggesting that cells from resting memory–like clusters defined in this study (CD62Lmed) are unlikely to migrate to or seed the periphery, contrary to activated memory–like DECs. Other correlations were observed between CD80 and CD95, and somewhat CD27, because these markers were grouped together via hierarchical clustering. This is in agreement with previous reports that >70% of rhesus macaque blood CD27+ B cells were also positive for CD80 and CD95 (55).
Multiparametric data analysis has also been proposed as a method to predict the trajectories of phenotypic development (27, 31, 33). Subtle phenotypic differences existed between the three vaccine-induced resting memory–like SPADE DECs (9, 16, 27), which all showed the highest degree of correlation with the Ab response in our dataset. If cells in these clusters shared a common progenitor, we speculate that the unique features of DEC9 compared with DEC16 and DEC27, such as higher levels of CD20 and HLA-DR and lower levels of CD40, CD45RA, and CD62L, are associated with the striking contraction of DEC9 cells between D8PB and D28PB, which was not evident for DEC16 and DEC27. Activated memory–like DEC23, DEC8, and DEC66, all members of phenotypic group A, showed overall phenotypic similarity but differed in their levels of certain markers such as CD22. CD22, alongside other markers such as CD69, have been used previously (26, 94) to identify recently activated cells, and could therefore be applied to derive information about the activation history of different clusters, especially in comparison with other clusters in the same phenotypic group. The high CD22 expression of DEC66 was similar to the majority of clusters in phenotype group A, and suggests these clusters were activated more recently compared with DEC23 and DEC8. In this respect, it would be informative for future studies to follow the size of memory-like clusters more frequently and after successive vaccinations, in an attempt to identify other markers that correlate with activation history or differentiation status.
To conclude, this proof-of-concept study used multiparametric clustering of single-cell mass cytometry data to identify groups of phenotypically similar B cells that expanded or contracted over time after MVA vaccination and those that correlate with the Ab response. This study can be used to guide the identification of B cell subsets at a more intricate level, and has applications for research into pathological or protective immune responses.
Acknowledgements
We thank the staff of IDMIT infrastructure and especially the FlowCyTech core laboratory at IDMIT (Fontenay-aux-Roses, France), François Huetz (Institut Pasteur, Paris) for recommendations on the Ab panel design, and Bernard Verrier (Biologie tissulaire et ingénierie thérapeutique, Institute of Biology and Chemistry of Proteins) for providing wild-type MVA and ELISA and plaque-reduction neutralization tests protocols. MVA HIV-B (MVATG17401) was obtained through a material transfer agreement between the Vaccine Research Institute, the Agence nationale de recherches sur le SIDA et les hépatites virales, INSERM, and the Commisariat à l'Energie Atomique et aux Energies Alternatives.
Footnotes
This work was supported by the Laboratories of Excellence–2011 Vaccine Research Institute, Créteil, France, funded by the Investment in the Future Program under Grant ANR-10-LABX-77; Facilities of Excellence–2010 FlowCyTech funded by the Investment in the Future Program under Grant ANR-10-EQPX-02-01; National Infrastructures in Biology and Health–2011 Infectious Disease Models and Innovative Therapies funded by Investment in the Future Program under Grant ANR-11-INBS-0008; and European Commission Advanced Immunization Technologies Grant FP7-HEALTH-2011-280873.
The online version of this article contains supplemental material.
Abbreviations used in this article:
- ACCENSE
automatic classification of cellular expression by nonlinear stochastic embedding
- D0
day zero
- DEC
differentially enriched cluster
- D8PB
day 8 postvaccine boost
- D28PB
day 28 postvaccine boost
- IDMIT
Infectious Disease Models and Innovative Therapies
- KS
Kruskal stress
- MSI
mean signal intensity
- MVA
modified vaccinia virus Ankara
- nAb
neutralizing Ab
- PB
plasmablast
- PC
plasma cell
- RT
room temperature
- SB
staining buffer
- SNE
stochastic neighbor embedding
- SPADE
spanning-tree progression analysis of density-normalized event
- TT
tetanus toxoid.
References
Disclosures
The authors have no financial conflicts of interest.