A Single-Cell Atlas of Lymphocyte Adaptive Immune Repertoires and Transcriptomes Reveals Age-Related Differences in Convalescent COVID-19 Patients

COVID-19 disease outcome is highly dependent on adaptive immunity from T and B lymphocytes, which play a critical role in the control, clearance and long-term protection against SARS-CoV-2. To date, there is limited knowledge on the composition of the T and B cell immune receptor repertoires [T cell receptors (TCRs) and B cell receptors (BCRs)] and transcriptomes in convalescent COVID-19 patients of different age groups. Here, we utilize single-cell sequencing (scSeq) of lymphocyte immune repertoires and transcriptomes to quantitatively profile the adaptive immune response in COVID-19 patients of varying age. We discovered highly expanded T and B cells in multiple patients, with the most expanded clonotypes coming from the effector CD8⁺ T cell population. Highly expanded CD8⁺ and CD4⁺ T cell clones show elevated markers of cytotoxicity (CD8: PRF1, GZMH, GNLY; CD4: GZMA), whereas clonally expanded B cells show markers of transition into the plasma cell state and activation across patients. By comparing young and old convalescent COVID-19 patients (mean ages = 31 and 66.8 years, respectively), we found that clonally expanded B cells in young patients were predominantly of the IgA isotype and their BCRs had incurred higher levels of somatic hypermutation than elderly patients. In conclusion, our scSeq analysis defines the adaptive immune repertoire and transcriptome in convalescent COVID-19 patients and shows important age-related differences implicated in immunity against SARS-CoV-2.

Introduction

T and B lymphocytes are crucial for protection from SARS-CoV-2 infection, viral clearance and the formation of persisting antiviral immunity (1, 2). Yet, adaptive immune responses have also been implicated in contributing to immunopathology during COVID-19, with higher mortality rates in elderly individuals (3–6). However, the exact determinants of a successful adaptive immune response against SARS-CoV-2 and its variability between different age groups remain to be fully elucidated.

Lymphocytes express either T cell receptors (TCR) or B cell receptors (BCR), which possess a highly diverse pair of variable chains [variable alpha (Vα) and beta (Vβ) for TCR and variable light (V_L) and heavy (V_H) for BCR] that are able to directly engage with antigen (e.g., viral proteins or peptides). Diversity in these variable chains are generated by somatic recombination of V-, D- and J-gene germline segments and along with combinatorial receptor chain pairing and somatic hypermutation (BCR only) results in an estimated human TCR and BCR diversity of 10¹⁸ and 10¹³, respectively (7, 8). Upon encountering cognate antigens, lymphocytes are phenotypically activated and undergo massive proliferation, also referred to as clonal selection and expansion. Deep sequencing of TCRs and BCRs has become a powerful strategy to profile the diversity of immune repertoires and to reveal insights on clonal selection, expansion and evolution (somatic hypermutation in B cells) (9–12) and has been instrumental in studying long term effects following vaccination, infection and ageing (13–17). In the context of COVID-19, immune repertoire sequencing has shown diminished TCR repertoire diversity and BCR isotype switching and respective expansion during early disease onset (18).

In recent years, single-cell sequencing (scSeq) of transcriptomes has progressed substantially through the development and integration of technologies such as cell sorting, microwells and droplet microfluidics (19, 20); most notably commercial systems like those of 10X Genomics have been established and are providing standardized protocols for wider implementation of scSeq. To find interactions across multiple genes and cells, analysis and visualisation of this high dimensional single-cell data is facilitated by clustering and nonlinear dimensionality reduction algorithms [e.g., t-distributed stochastic neighbor embedding (t-SNE) or Uniform Manifold Approximation and Projection (UMAP)] (21, 22). scSeq of transcriptomes has been used extensively to profile the gene expression signatures of T and B cells to identify novel cellular subsets and phenotypes as well as their response to vaccination, infection and cancer (23–26). Furthermore, clustering with scSeq data enables the unbiased identification of cellular states and analyses of the broad continuum of T and B cell populations as well as their differentiation trajectories (27). In the context of patients with severe symptoms of COVID-19, scSeq has revealed a dysfunctional T cell response of interferon expression combined with elevated levels of exhaustion (28).

In addition to transcriptome sequencing, a major advantage of scSeq is that it also enables information on the native pairing of TCR Vα and Vβ chains and BCR V_L and V_H chains (29–32), which was not previously possible with the standard bulk sequencing of lymphocytes as these receptor chains are expressed as unique transcripts from separate chromosomes (33). Coupling TCR or BCR sequence to the transcriptome within an individual cell enables phenotypic analyses of a clonal population of lymphocytes and their dynamics (34–36). scSeq of transcriptomes and immune repertoires in COVID-19 patients with severe symptoms has shown a high level of clonal expansion in specific T cell subsets (Th1, Th2, and Th17) and preferential germline gene usage in clonally expanded B cells (28, 34, 37); while a more recent study found a positive correlation between clonal expansion of effector-like CD8+ T cells and disease severity (38).

An important question that remains to be answered is whether there are age-related differences in mounting a successful adaptive immune response against SARS-CoV-2. Here, we perform scSeq on the immune repertoires and transcriptomes of T and B cells derived from eight convalescent COVID-19 patients of two different age groups (mean ages = 31 and 66.8 years) at one month of convalescence following mild to moderate disease. We observed preferential clonal expansion of effector CD8+ T cells across all patients, although a significantly higher CD8-to-CD4 T cell ratio was detected in young patients of our cohort. Further, clonally expanded B cells in young patients displayed significantly higher levels of somatic hypermutation and an increased immunoglobulin (Ig) class-switching compared to clonally expanded B cells from older patients. Our analyses serve as a valuable resource for future scSeq characterization of SARS-CoV-2 adaptive immunity and highlight important age-related differences in the adaptive immune status of convalescent COVID-19 patients.

Results

Study Design and Single-Cell Profiling of Convalescent COVID-19 Patient Lymphocytes

We performed scSeq of immune receptor repertoires and transcriptomes of lymphocytes from convalescent COVID-19 patients to characterize the adaptive immune response against SARS-CoV-2. For this purpose, we selected eight patients enrolled in the SERO-BL-COVID-19 clinical study (39), all of which fully recovered from COVID-19 without requiring hospitalization or the administration of supplemental oxygen. Patients tested positive for the presence of SARS-CoV-2 after RT-PCR of naso/oropharyngeal swab samples (day 0), displayed COVID-19 symptoms for 4-14 days, and showed positive seroconversion at the time of blood collection (mean sample collection time = 32.5 ± 4.1 days post-symptom onset) (Figure 1A and Supplemental Table 1). Since COVID-19 often affects older patients more severely (40), subjects were divided into two groups according to their age, namely Group 1 (mean = 66.75 ± 6.9 years) and Group 2 (mean = 31 ± 5.9 years), with the aim of investigating potential differences in their responses against SARS-CoV-2. In addition to older age, significant differences in Group 1 versus Group 2 included elevated IgA/IgG SARS-CoV-2-specific antibody levels and an increased duration of COVID-19 symptoms (Figure 1B and Supplemental Figure 1). Patients from Group 1 also experienced an increased severity of COVID-19 symptoms relative to Group 2 (Supplemental Table 1). Despite increased symptom duration in the older cohort, correlation of this parameter with age was only modest (R² = 0.4647), likely reflecting the small sample size (Supplemental Figure 1C). In summary, patient selection criteria included IgG+ patients according to PoC test, age (youngest and oldest) and similar symptoms/severity. Compared to the other patients enrolled in the SERO-BL-COVID-19 clinical study and using a cutoff set at 40 years of age, the selected older patients have had a longer disease duration (12.25 vs. 9.1 days), whereas the selected younger patients had a shorter disease duration (6.75 vs. 8.5 days).

FIGURE 1

Figure 1 Overview of single-cell transcriptome and immune receptor profiling of convalescent COVID-19 patient lymphocytes. Convalescent COVID-19 patients enrolled in the SERO-BL-COVID-19 study were selected according to their age for single-cell sequencing analysis of their T cells and B cells. (A) Timeline illustrates symptom onset, symptom resolution and collection of blood samples from individual patients relative to the time of positive SARS-CoV-2 RT-PCR test (day 0). (B) Graph displays the ages and duration of COVID-19 symptoms in individual patients. Dotted lines show the mean duration of symptoms in the young (y = 6.75 days) and old (y = 12.25 days) groups. A significant difference in symptom duration between groups is indicated with an asterisk (p = 0.0127; unpaired t-test). (C) Single-cell sequencing protocol. Whole blood was collected following the resolution of COVID-19 symptoms and subjected to density gradient separation for isolation of PBMC. T cells and B cells from individual patients were purified from PBMC using negative immunomagnetic enrichment, pooled (intra-patient) and prepared for droplet generation using the 10x Genomics Chromium system. Single cells were emulsified with DNA-barcoded gel beads and mRNA transcripts were reverse-transcribed within droplets, resulting in the generation of first-strand cDNA molecules labelled with cell-specific barcodes at their 3’ ends (added by template switching). Emulsions were disrupted and cDNA was amplified by means of PCR for further processing of transcriptome libraries. Transcriptome libraries from individual patients were indexed and multiplexed for deep sequencing using the Illumina NovaSeq platform. Targeted enrichment of recombined V(D)J transcripts was performed by PCR and the resulting products were processed for the generation of BCR and TCR libraries, which were then indexed, multiplexed and deep-sequenced.

To profile patient lymphocytes, we isolated peripheral blood mononuclear lymphocytes (PBMC) from blood and purified T cells and B cells by negative immunomagnetic enrichment. Plasma cells (PCs) were depleted from PBMC samples prior to this step for scSeq in a companion study (Ehling et al., manuscript in preparation), and thus were excluded from our analyses. After purification, T cells and B cells underwent the 10X genomics protocol for scSeq 5’ library preparation, which included gel encapsulation single-cell barcoding of mRNA, followed by cDNA generation through polydT reverse transcription. Finally, after full-length V(D)J segment enrichment, construction of TCR and BCR V(D)J and transcriptome sequencing libraries was done according to the V(D)J enrichment and 5’ library construction kits, respectively. Deep sequencing of immune repertoires and transcriptomes was performed using the Illumina NovaSeq with paired-end 26 x 91 bp cycles per read. For TCR and BCR V(D)J and transcriptome libraries, we recovered on average 20.000 and 10.000 reads per cell, respectively (Figure 1C).

Single-Cell Transcriptome Analysis Defines Major T and B Cell Subsets

Bioinformatic filtering was performed to exclude the following: Cell doublets, cells with a very low or high number of genes, and T cells with no detectable expression of CD8 and CD4 (see Methods), which resulted in the identification of 30,096 cells in total from all eight patients. Cells were then split into CD8+ T cell (Figure 2A, 7,353 cells), CD4+ T cell (Figure 2B, 8,334 cells) and B cell (Figure 2C, 14,409 cells) datasets. In order to reduce the dimensionality of the data, while preserving the global structure, we used UMAP for better visualisation and interpretation purposes (41). UMAP and unsupervised clustering of these subgroups led to the identification of eleven dominant cell subsets (Figures 2A–C). CD8+ T cells clustered into naïve (SELL+, TCF7+), memory (IL7R+, CD40LG+) and effector cells (GZMB+, NKG7+) (Figures 2A, D), which also encompassed exhausted CD8+ T cells (Supplemental Figure 2). We identified four different CD4+ T cell subsets, namely naïve (SELL+, LEF1+), memory (S100A4+), effector (CCL5+, GZMK+) and regulatory cells (FOXP3+) (Figures 2B, E). The B cell compartment consisted of naïve (CD23A+), marginal zone (MZ) (FCRL3+, CD1C+), activated (CD83+) and memory cells (CD27+, TACI+) (Figures 2C, F). Notably, clustering of single cells based on transcriptome data revealed a trajectory that reflected a progression in lymphocyte differentiation from naïve to effector (or activated) subsets (Figures 2A–C). Pseudotime analysis of the dataset supports this differentiation trajectory in CD8+ and CD4+ T cells (Figures 2G, H). Interestingly, pseudotime analysis of B cell data not only showed naïve-MZ-memory and naïve-MZ-activated trajectories, but also a third MZ-memory-activated trajectory that suggests the presence of reactivated memory B cells, possibly through antigen encounter (Figure 2I).

FIGURE 2

Figure 2 Single-cell transcriptomic analysis delineates major T and B cell subsets. (A–C), Uniform manifold approximation and projection (UMAP) plots of major cellular subsets identified within the CD8+ T cell (A), CD4+ T cell (B) and B cell (C) populations. Cells from all patients are displayed in each plot.9 (D–F), UMAP plots showing the expression levels of selected genes used to delineate cellular subsets within the CD8+ T cell (D) CD4+ T cell (E) and B cell (F) populations. Cells from all patients are displayed in each plot. g-i, Graphs display pseudotime and trajectory inference analysis applied to CD8+T cell (G), CD4+ T cell (H) and B cell (I) clusters. (J–L), Bar graphs show the proportions of identified cellular subsets within the CD8+ T cell (J), CD4+ T cell (K) and B cell (L) populations in each patient. CCL5, C-C Motif Chemokine Ligand 5; CD27, TNFRSF7; CD40LG, CD40 ligand; FCER2, Fc Fragment of IgE Receptor II (also: CD23a); FCRL3, Fc Receptor Like 3; FOXP3, Forkhead Box Protein P3; GZMB, Granzyme B; GZMK, Granzyme K; IL7R, Interleukin-7 Receptor; LEF1, Lymphoid Enhancer Binding Factor 1; NKG7 , Natural Killer Cell Granule Protein 7; S100A4, S100 Calcium Binding Protein A4; SELL, Selectin L; TCF7, Transcription Factor 7; TNFRSF13B, TNF Receptor Superfamily Member 13B.

Having defined the major T cell and B cell subsets from pooled patient data, we next compared their proportions across patients (Figures 2J–L) and between different age groups (Supplemental Figure 3). We found that young patients had a significantly higher CD8-to-CD4 T cell ratio relative to older ones (Supplemental Figure 3A), which may reflect a previously reported age-dependent difference (42). Interestingly, despite this reduction, there was a trend that older patients had a higher proportion of effector CD8+ T cells relative to their younger counterparts (Supplemental Figure 3B). While this difference was not significant, it is consistent with the increased symptom severity experienced by older patients (Supplemental Table 1), a feature that has been associated with elevated proportions of effector/exhausted CD8+ T cells in the periphery (28). Of note, we found that older patients had a small but significant increase in CD4+ Tregs compared to young patients (Supplemental Figure 3H), and that increased proportions of MZ B cells occurred in two of the older patients (Supplemental Figure 3L). Taken together, our data highlights the diversity of elevated responses in specific patients across age groups, as exemplified by individuals with a high abundance of effector CD8+ T cells (e.g., Pt-2 and Pt-3) and/or activated B cells (e.g., Pt-2, Pt-5 and Pt-8).

Single-Cell Profiling of Immune Receptor Repertoires Identifies Highly Expanded TCR and BCR Clonotypes

We next determined the clonal expansion levels of T cells and B cells in convalescent COVID-19 patients by quantifying the number of cells expressing unique TCRs (Figure 3A) or BCRs (Figure 3C) (clonotype definition in the methods section). We found substantial heterogeneity in T cell clonal expansion levels across patients, with the highest number of expanded TCR clonotypes occurring in four patients, namely Pt-2 and Pt-3 (Group 1), Pt-7 and Pt-8 (Group 2). Within these patients, Pt-2 displayed the largest amount of expanded TCR clonotypes, which is consistent with the high abundance of effector CD8+ T cells in this subject (Figure 2G). Analysis of TCRα and TCRβ germline V-gene usage in the ten most expanded clonotypes per patient revealed a frequent occurrence of TRBV20-1 (7 out of 8 patients) and TRAV-29/DV5 genes (5 out of 8 patients), though pairing of these germline genes was not observed (Figure 3B and Supplemental Figure 4). In agreement with the overall higher expansion of CD8+ effector over CD4+ effector T cell subsets (Figures 2G, H), we found that the vast majority (85%) of the ten most expanded TCR clonotypes per patient originated from CD8+ T cells (Supplemental Figure 5). Based on this observation, we genotyped patient HLA class I alleles by means of amplicon deep sequencing (Supplemental Table 2). We found that the two patients with the highest levels of T cell clonal expansion (i.e., Pt-2 and Pt-8) shared the HLA-A*0201 allele, as well as a number of TRBV and TRAV genes in their most expanded clonotypes, which indicates a possible convergence towards germlines that may be related to SARS-CoV-2 specificity. Analysis of single-cell BCR repertoire sequencing data revealed that highly expanded BCR clonotypes occurred more frequently in younger patients, for example in Pt-5, Pt-6 and Pt-8 (Figure 3C). This is an unexpected finding, particularly as older patients in our cohort displayed significantly higher SARS-CoV-2-specific IgA/IgG levels in serum (Supplemental Figure 1E, F). Thus, this suggests that older patients may harbor a higher diversity of relatively unexpanded SARS-CoV-2-specific B cells. Supporting this observation, we found that older patients had a wider range of heavy chain complementarity determining region 3 (CDR3H) lengths relative to younger ones, indicating a possible larger degree of variability in their antibody paratopes (Supplemental Figure 6). Analysis of heavy chain and light chain germline V-gene pairing in the ten most expanded BCR clonotypes per patient revealed a frequent occurrence of the IGHV-3-23/IGKV-3-20 pairing (7 out 8 patients) (Figure 3D and Supplemental Figure 7). However, as this pairing is the most frequently found in healthy cohorts (43, 44), such antibodies may not necessarily be enriched for SARS-CoV-2 specificity.

FIGURE 3

Figure 3 Single-cell profiling of immune repertoires highlights differential levels of inter-patient T cell and B cell clonal expansion. (A, B), Analysis of T cell clonal expansion in convalescent COVID-19 patients. (A) Bar graphs show T cell clonal expansion, as determined by the number of cells identified per TCR clonotype. Each box represents the size of individual TCR clonotypes. TCR clonotypes present in more than one cell are shown. (B) Circos plots display V-gene usage in the top ten most expanded TCR clonotypes for each patient. The size and colour (dark to light) of outer bars reflect the relative abundance of T cells expressing specific V-genes on a per class basis (top: TCRα chain, bottom: TCRβ chain). (C–E), Analysis of B cell clonal expansion in convalescent COVID-19 patients. (C) Bar graphs show B cell clonal expansion, as determined by the number of cells identified per BCR clonotype. Each box represents the size of individual BCR clonotypes. BCR clonotypes present in more than one cell are shown. (D) Circos plots display V-gene usage in the top ten most expanded BCR clonotypes for each patient. The size and colour (dark to light) of outer bars reflect the relative abundance of B cells expressing specific V-genes on a per class basis (top: Ig light chain, bottom: Ig heavy chain). (E) Graph displays the levels of somatic hypermutation (SHM) in unexpanded (1 cell), expanded (2-4 cells) and highly expanded (5 cells) BCR clonotypes across patients. SHM levels are based on the percentage similarity between Ig heavy chain V-genes and their corresponding germlines. Data are displayed as median ± IQR. (F), Graph displays SHM levels in highly expanded BCR clones (≥5 cells) of old (n = 24 clones) and young (n = 29 clones) patients. Asterisks indicate a significant difference in SHM levels between groups (p = 0.0085; unpaired t-test). Data are displayed as median ± IQR. (G), Bar graphs show the Ig isotype distribution in unexpanded (1 cell), expanded (2-4 cells) and highly expanded (≥5 cells) BCR clonotypes across patients.

To further characterize BCR repertoires across different levels of clonal expansion we divided clonotypes into additional subsets: unexpanded (1 cell per clonotype), expanded (2-4 cells per clonotype) and highly expanded (≥5 cells per clonotype) and assessed their levels of somatic hypermutation (SHM). We found that the degree of SHM largely correlated with clonal expansion, with expanded and highly expanded clonotypes having higher SHM (i.e., more divergent from their germline V-genes) than unexpanded ones (Figure 3E). Strikingly, highly expanded BCR clonotypes from young patients had significantly higher SHM levels compared to older patients, potentially indicating more efficient affinity maturation had occurred in response to SARS-CoV-2 antigens (Figure 3F). In line with this, SHM levels of IgG⁺ expanded BCR clonotypes are significantly higher in young patients (Supplemental Figure 8). Finally, we examined the distribution of Ig isotypes across clonal expansion groups. As expected, IgM was the most frequent isotype in unexpanded BCR clonotypes, with the proportion of this isotype being reduced in expanded BCR clonotypes (2-4 cells) of all patients. Conversely, the proportions of IgG and IgA isotypes in expanded clonotypes increased for all patients, thus indicating class-switching in response to clonal expansion. Analysis of Ig isotype distribution in highly expanded BCR clonotypes (≥5 cells) revealed that a subset of patients harbored a vast majority of class-switched IgA (Pt-5, Pt-6 and Pt-7) or IgG (Pt-4) (Figure 3G). Notably, we found that Ig isotype class-switching in highly expanded clonotypes was correlated with SHM levels across patients (Figure 3E), highlighting the temporal connection between affinity maturation and class-switching processes in the germinal center (45).

Single-Cell Transcriptome and TCR Profiling Reveals Predominant Cytotoxic Programs in Highly Clonally Expanded CD4+ and CD8+ T Cells

We next investigated the patterns of clonal expansion in different T cell subsets by mapping single-cell TCR sequencing data onto individual CD8+ and CD4+ T cells visualized by UMAP (Figures 4A, B). For this analysis, we identified a total of 4,730 CD8+ T cells and 5,509 CD4+ T cells with available TCR clonotype and transcriptome information. Both CD8+ and CD4+ T cells showed increased levels of clonal expansion when progressing from naïve to effector phenotypes, with highly expanded TCR clonotypes (≥5 cells) almost exclusively expressed by effector T cells (Figures 4A, B, 2A, B). As previously observed (Supplemental Figures 3, 5), CD8+ T cells showed substantially higher levels of clonal expansion relative to CD4+ T cells, in which highly expanded clonotypes were rare. Notably, young patients (Group 2) had a markedly higher abundance of unexpanded CD8+ T cell clonotypes compared to older patients (Group 1), which could indicate an ongoing resolution of their CD8+ T cell response at the analyzed timepoint. Patients with high levels of CD8+ T cell clonal expansion, however, were identified across age groups (i.e., Pt-2, Pt-3, Pt-7 and Pt-8), with Pt-2 (Group 1) showing the highest abundance of highly expanded clonal T cells. We further explored the relationship between clonal expansion and T cell phenotype by performing differential gene expression analysis in unexpanded, expanded and highly expanded T cell clonotypes (Figures 4C, D). CD8+ T cells with high clonal expansion displayed elevated cytotoxicity (PRF1, GZMH, GNLY), activation (NKG7, CCL5), inflammation (NFGBIA, S100A4, S100A6) and type I interferon-induced (IFITM2) markers in all patients. Additionally, components of MHC class I (HLA-A, HLA-B, HLA-C and B2M) were also increased in this subgroup, indicating increased IFN-γ-induced activation (46). Conversely, unexpanded CD8+ T cell clonotypes across age groups displayed upregulated markers found in naïve and memory CD8+ T cell subsets (IL7R, LTB) (47), as well as markers likely associated with homeostatic proliferation (LDHB, NOSIP, EEF1B2, NPM1, TPT1, PABPC1).

FIGURE 4

Figure 4 Single-cell transcriptome and TCR sequencing reveals preferential clonal expansion in effector T cells. (A, B), UMAP plots display CD8+ (A) and CD4+ (B) T cells from specific patients according to their clonal expansion levels. T cells from other patients in each individual plot are shown in grey. (C, D), Heatmaps show differential gene expression (DGE) in unexpanded, expanded and highly expanded CD8+ (C) or CD4+ (D) T cells. Genes were filtered to include those with detectable expression in at least 50% of cells and that had a minimum 50% fold-change in expression level between groups.

Despite the low levels of clonal expansion observed in CD4+ T cells, we identified distinct gene expression signatures in the highly expanded clonotypes that were present in 5 of 8 patients. Similar to CD8+ T cells, highly expanded CD4+ T cell clonotypes showed upregulation of genes related to activation (CCL5), cytotoxicity (GZMA), inflammation (IL32, CD99, NFKBIA) and MHC class I molecules (HLA-A, HLA-B, HLA-C and B2M), while unexpanded CD4+ T cells displayed markers of naïve T cells (SELL) and proliferation markers (LDHB, NOSIP, PABPC1). Interestingly, both CD8+ and CD4+ memory T cells show higher expression of activation/proliferation markers (Granzyme A and K, NFKBIA, DUSP1) (48–50) in older compared to younger patients, potentially reflecting more recent viral clearance (Supplemental Figure 9). Taken together, integration of TCR sequencing and transcriptome data reveals clonal expansion as a hallmark of effector T cell subsets, and highlights a dominant role of CD8+ T cells in possible clonal responses against SARS-CoV-2 in convalescent COVID-19 patients.

Transcriptomic and BCR Profiling of Single B Cells Reveals Plasma Cell Transition, Class-Switching and SHM Patterns

We next integrated single-cell BCR sequencing data onto the B cell transcriptional landscape to relate clonal expansion, SHM and isotype distribution to different B cell phenotypes. For this analysis we identified 11,227 individual B cells with available BCR and transcriptomic information. We observed a generally low level of B cell clonal expansion, with preferential localization of expanded B cell clonotypes (2-4 cells) to the memory and MZ B cell regions and a rare occurrence of highly expanded (≥5 cells) clonotypes in most patients (Figure 5A). Analysis of differential gene expression showed that expanded clonotypes had increased expression of genes involved in cytoskeleton reorganization (VIM) and genes associated with the unfolded protein response (HSBP90, CALR, PPIB), indicating B cell activation and transition into plasma cells, respectively (51–53) (Figure 5B). Furthermore, expanded B cell clonotypes across patients showed downregulation of MHC class II genes (CD74, HLA-DR, -DQA1, -DRB1), further supporting their trajectory towards antibody-producing plasma cells (54).

FIGURE 5

Figure 5 Single-cell transcriptome and BCR profiling reveals elevated class-switching and somatic hypermutation levels in memory B cells. (A) UMAP plots display B cells from specific patients according to their clonal expansion levels. B cells from other patients in each individual plot are shown in grey. (B) Heatmap shows differential gene expression (DGE) in unexpanded, expanded and highly expanded B cells. Genes were filtered to include those with detectable expression in at least 50% of cells and that had a minimum 50% fold-change in expression level between groups. (C) Graph displays the levels of somatic hypermutation (SHM) in the BCRs of naïve, activated and memory B cells across patients. SHM levels are based on the percentage similarity between BCR heavy chain V-gene and its corresponding germline. Data are displayed as median ± IQR. (D) Graph shows the distribution of B cells expressing specific Ig isotypes relative to their location in transcriptome UMAP plots. B cells from all patients are shown. (E) Bar graphs show Ig isotype distribution of BCRs found in naïve, activated and memory B cells across patients. CALR, Calreticulin; CD74, HLA class II Histocompatibility Antigen Gamma Chain; DUSP1, Dual Specificity Phosphatase 1; HLA-D, Major Histocompatibility Complex, Class II; HSP90B1, Heat Shock Protein 90 Beta Family Member 1; PPIB, Peptidylprolyl Isomerase B; VIM, Vimentin.

We next analysed single-cell BCR sequencing data to assess SHM levels in different B cell subsets. BCRs from memory B cells showed the highest levels of SHM across patients, while BCRs extracted from naïve and activated B cells displayed similarly low median SHM values. Notably, however, activated B cells expressed a larger number of high-SHM outliers than naïve B cells, suggesting ongoing affinity maturation in this subset (Figure 5C). Mapping of Ig isotype information onto the B cell transcriptional UMAP space, revealed an even distribution of IgM expression across B cell subsets, rare occurrence of IgD-expressing B cells and, most notably, confinement of class-switched IgG- and IgA-expressing B cells to the memory and MZ B cell regions (Figure 5D). Finally, assessment of Ig isotype distribution across patients and B cell subsets revealed that as expected the vast majority of naïve B cells expressed the IgM isotype (Figure 5E), with minimal levels of class-switching observed in activated B cells but prominent class-switching to IgG and IgA isotypes in the memory B cell compartment across all patients (Figure 5E).

Together our results indicate ongoing transition of clonally-expanded B cells into antibody-producing plasma cells, as well as high levels of SHM and Ig isotype class-switching in the memory B cells of convalescent COVID-19 patients.

Computational Prediction of Shared Specificity Identifies Candidate SARS-CoV-2-Specific TCRs

Motivated by the high levels of CD8+ T cell clonal expansion and activation observed in convalescent COVID-19 patients, we further analyzed single-cell TCR repertoires for potential SARS-CoV-2 specificity. To this end, we applied GLIPH2, an algorithm developed by M. Davis and colleagues that clusters TCRs with a high probability of recognizing the same epitope into specificity groups (based on conserved motifs and similarity levels in CDR3β) (55). In addition, the provision of HLA typing data enables the prediction of HLA restriction in specific TCR clusters. Analysis of 23,010 paired TCRα and TCRβ sequences derived from the CD8+ T cells of eight patients led to the identification of a total of 552 specificity groups with attributed HLA restriction (seven alleles). We observed distinct proportions of shared specificity groups between pairs of patients, with Pt-1:Pt-8, Pt-3:Pt-7 and Pt-1:Pt-7 showing the highest overlap (Figure 6A). Furthermore, the vast majority of clusters were attributed with HLA-A*01:01, HLA-A*03:01, HLA-B*13:02 or HLA-C*03:04 restriction, and HLA-A*02:01, A*24:02 and C*04:01 were attributed to less than 15% of TCR clusters (Figure 6B). While some of these clusters may be defined by SARS-CoV-2 specificity, it is difficult to exclude reactivity to common human viruses (e.g., CMV, EBV). To further investigate potential for SARS-CoV-2 specificity, we analyzed the sequences of known HLA-A*02:01-restricted SARS-CoV-2-specific, CMV-specific and EBV-specific TCRs alongside those derived from patients expressing the HLA-A*02:01 allele (i.e., Pt-2, Pt-4 and Pt-8) (Supplemental Tables 3, 4). GLIPH2 identified 35 unique patient TCR sequences that clustered together with known SARS-CoV-2-specific TCRs, with the great majority originating from Pt-8 (Figure 6C and Table 1). Thus, such TCRs represent candidates for mediating CD8+ T cell immunity against SARS-CoV-2 infection.

TABLE 1

Table 1 Patient TCRs with predicted SARS-CoV-2-specificity, as determined by GLIPH2 clustering with known HLA-A*02:01-restricted SARS-CoV-2-specific TCRs.

FIGURE 6

Figure 6 GLIPH2 analysis of single-cell paired TCR repertoires reveals candidate SARS-CoV-2-specific TCRs. (A) Heatmap shows the proportion of TCR specificity groups containing sequences from specific pairs of patients, as determined by GLIPH2 analysis (total TCR clusters = 552). (B) Bar plot displays the proportions of predicted HLA class I alleles in HLA-attributed TCR specificity groups (total TCR clusters = 552). (C) Graph displays the proportions and numbers of candidate SARS-CoV-2-specific TCRs derived from HLA-A*0201-positive patients, as determined by GLIPH2 clustering with known SARS-CoV-2-specific TCR sequences.

Limitations of the Study

Our sample size of eight patients (four per age group) is small and reduces the number of conclusions that we can confidently make from the observed data. Nevertheless, single-cell data offers a deeper characterization of each patient than normal bulk transcriptome and repertoire studies. Furthermore, while the time between symptom onset and sample collection is highly uniform across patients, the time between symptom resolution and sample collection is significantly shorter in the old patient group due to prolonged symptom duration. This is an important variable that should be considered when interpreting the findings presented here.

Discussion

Here we apply scSeq for in-depth immune repertoire and transcriptomic analysis of T cells and B cells derived from non-severe COVID-19 patients at one month of convalescence. Our analyses of transcriptomic data defined eleven T cell and B cell subsets, of which the effector CD8+ T cell subset (GZMB, NKG7) showed the highest levels of expansion in specific patients, both in terms of proportion and clonality. These findings are in agreement with recent scSeq studies of convalescent COVID-19 patients (28, 56). For example, effector tissue-resident CD8+ T cells from bronchoalveolar lavage fluid were found to be highly clonally expanded in convalescent COVID-19 patients that experienced moderate but not severe infection (56). In addition, scSeq of PBMCs revealed increases in cytotoxic effector CD4+ and CD8+ T cell subsets in non-severe convalescent COVID-19 patients only (28). In our study, we observed high levels of clonal expansion of the effector CD8+ T cell subset but less evident expansion of CD4+ T cell subsets. Importantly, however, differential gene expression analysis revealed that both highly clonally expanded CD8+ and CD4+ T cells had elevated markers of cytotoxicity (CD8: PRF1, GZMH, GNLY; CD4: GZMA).

A recent study of functional T cell responses against SARS-CoV-2 reported significantly higher CD8+ T cell responses directed at spike, M/NP and ORF/Env epitopes in convalescent COVID-19 patients experiencing moderate symptoms compared to those recovering from severe infection (57). Further evidence supporting a potential role of CD8+ T cells in rapid viral clearance comes from the occurrence of SARS-CoV-2-specific T cells in asymptomatic seronegative family members of COVID-19 patients (58), as well as in samples from asymptomatic seronegative control subjects obtained during the COVID-19 pandemic but not prior to it (58, 59). In line with these findings, our analysis of clonal expansion in lymphocyte subsets suggest a key role of CD8+ effector T cells in the clearance and protection against SARS-CoV-2 in patients with moderate disease. Notably, the fact that younger patients in our cohort had significantly higher CD8-to-CD4 T cell ratios might have contributed to reduced symptom duration relative to older patients. In this context, we predict that methods for the identification of CD8-derived SARS-CoV-2-specific TCRs including functional assays (57, 58) and the application of motif clustering (59) or machine learning (60) to TCR repertoire data, as well as methods for the identification of their corresponding epitopes (61, 62) will become increasingly important tools for monitoring SARS-CoV-2 immunity.

In contrast to T cells, we found modest levels of B cell clonal expansion that was not exclusively restricted to a particular subset but spanned activated, memory and MZ B cells. These low levels of B cell clonal expansion are in contrast to the high IgG and IgA serum titers found in all patients of our cohort, particularly those in the older patient group. This suggests that the B cell compartment may experience clonal contraction at one month of convalescence, with the disconnect from serum titers of IgG likely explained by the long half-lives of secreted IgG (~3 weeks) (63). It should also be noted that antibody-producing plasma cells, which were not included in our analysis, may display clonal expansion levels that are more congruent with the observed levels of SARS-CoV-2-specific antibodies in serum. Although rare, highly expanded B cell clones showed elevated markers of plasma cell transition and activation, thus indicating possible ongoing differentiation into plasma cells, albeit at low levels, at one month after symptom onset.

Analysis of BCR repertoires and transcriptomes revealed that the highest levels of SHM and Ig class-switching occurred in the memory B cell subset. This finding is consistent with the generation of germinal center-derived memory B cells following antigen encounter (64, 65). Similarly, SHM and Ig class-switching levels were directly correlated with B cell clonal expansion. Remarkably, we found that highly expanded B cell clones from the young patient group had significantly higher SHM levels (median = 6.7% divergence from germline) than those derived from the old patient group (median = 2.5% divergence from germline). This occurred despite significantly higher levels of SARS-CoV-2-specific antibodies in the old patient group, and suggests that more effective affinity maturation may occur in younger patients. Initial studies characterizing SARS-CoV-2-specific antibodies reported low levels of SHM, with median divergence from germline ranging from 0.7-2% in convalescent patients of varied symptom severity analyzed at 20-40 days following symptom onset (66–68). Notably, a recent report has described ongoing affinity maturation of SARS-CoV-2-specific antibodies at six months following symptom onset in non-severe patients, with SHM levels rising to 3% divergence from germline (2). Interestingly, the cited study provides evidence of viable SARS-CoV-2 antigen in the gut of such patients, which has been proposed as a source of antigen for ongoing affinity maturation linked to elevated IgA serum titers. It is thus unclear whether increased SHM levels found in the young patients of our cohort are the result of a better capacity for affinity maturation upon initial antigen encounter or of ongoing affinity maturation resulting from longer exposure to antigen following symptom resolution. The clear dominance of IgA class-switched BCRs expressed by highly expanded B cells in 3 out of 4 patients in the young group appears to support the latter. It should be noted that we did not isolate SARS-CoV-2-reactive T or B cells and that our insights are based solely on the clonal and transcriptomic features of total T cells and B cells following scSeq.

In conclusion, our in-depth characterization integrating single-cell immune repertoire and transcriptome profiling of T and B cells represents a valuable resource to better understand the adaptive immune response and age-related differences in convalescent COVID-19 patients with moderate disease. Furthermore, it serves as an important point of reference for future single-cell characterization of lymphocytes at later time points of convalescence, or of lymphocytes isolated from patients experiencing long-term COVID-19 sequelae and can help in defining markers for clinical monitoring of disease progression.

Methods

Patient Samples

Patients were participants of the SERO-BL-COVID-19 study sponsored by the Department of Health, Canton Basel-Landschaft, Switzerland. All analyzed patients tested positive for SARS-CoV-2 after RT-PCR of naso- and oropharyngeal swab samples and experienced a resolution of COVID-19 symptoms without requiring hospitalization. Whole blood was collected 25 to 39 days following a positive RT-PCR test and subjected to density gradient centrifugation using the Ficoll Paque Plus reagent (GE Healthcare, #17-1440-02). After separation, the upper plasma layer was collected for semiquantitative ELISA detection of IgG and IgA SARS-CoV-2-specific antibodies (Euroimmun Medizinische Labordiagnostika, #EI2668-9601G, #EI2606-9601A). Peripheral blood mononuclear cells (PBMC) were collected from the interphase, resuspended in freezing medium (RPMI 1640, 10%(v/v) FBS, 10%(v/v) dimethyl sulfoxide) and cryopreserved in liquid nitrogen. Point-of-care lateral flow immunoassays assessing the presence of IgG and IgM SARS-CoV-2-specific antibodies (Qingdao Hightop Biotech, #H100) were performed at the time of blood collection. Both ELISA and point-of-care lateral flow tests employed Spike and Nuclear Capsid Protein for SARS-CoV-2-specific antibody detection.

Immunomagnetic Isolation of B Cells and T Cells

PBMC samples were thawed, washed in complete media (RPMI 1640, 10%(v/v) FBS) and pelleted by centrifugation. Cells were resuspended in 0.5 mL complete media, counted and treated with 10 U ml^-1 DNAse I (Stemcell Technologies, #) for 15 min at RT in order to prevent cell clumping. After DNase I digestion, cells were washed twice in complete media, pelleted by centrifugation and resuspended in 0.5 mL flow cytometry buffer (PBS, 2%(v/v) FBS, 2 mM EDTA). The cell suspension was filtered through a 40 μM cell strainer prior to immunomagnetic isolation. As a first step, plasma cells were isolated using the EasySep Human CD138 Positive Selection Kit II (Stemcell Technologies, #17877) for analysis in a companion study (manuscript in preparation). The negative fraction of the above selections was divided into two aliquots that were subjected to negative immunomagnetic isolation of either B cells (EasySep Human Pan-B cell Enrichment Kit, Stemcell Technologies, #19554) or T cells (EasySep Human T cell Isolation Kit, Stemcell Technologies, #17951). After isolation, B cells and T cells were pelleted by centrifugation, resuspended in PBS, 0.4%(v/v) BSA, filtered through a 40 μM cell strainer and counted. T cells and B cells originating from the same patient were pooled in equal numbers and the final suspension was counted and assessed for viability using a fluorescent cell counter (Cellometer Spectrum, Nexcelom). Whenever possible, cells were adjusted to a concentration of 1x10⁶ live cells/mL in PBS, 0.04%(v/v) BSA before proceeding with droplet generation.

Single Cell Droplet Generation and Preparation of Sequencing Libraries

Encapsulation of lymphocytes and DNA-barcoded gel beads was performed using the Chromium controller (10x Genomics, PN-110203). Briefly, 1.4x10⁴ to 1.7x10⁴ cells (in reverse transcription mix) were loaded per channel for a targeted recovery of 8x10³ to 1x10⁴ cells per sample. Reverse transcription and preparation of single-cell transcriptome, BCR and TCR libraries was performed according to the manufacturer’s instructions (CG000086 manual, RevM, 10x Genomics) and using the following kits: Chromium Single Cell 5’ Library & Gel Bead Kit (PN-1000006), Chromium Single Cell 5’ Library Construction Kit (PN-1000020), Chromium Single Cell V(D)J Enrichment Kit, Human T Cell (PN-1000005), Chromium Single Cell V(D)J Enrichment Kit, Human B Cell (PN-1000016), Chromium Single Cell A Chip Kit (PN-1000009), Chromium i7 Multiplex Kit (PN-120262).

Deep Sequencing

The quality and concentrations of transcriptome (i.e., cDNA), TCR and BCR libraries were determined using a fragment analyzer (Agilent Bioanalyzer) at specific steps of library preparation, as recommended in the 10x Genomics scSeq protocol (CG000086 manual, RevM). Following multiplexing (Chromium i7 Multiplex Kit, #PN-120262, 10x Genomics), transcriptome libraries were treated with free adapter blocking (FAB) reagent to prevent index switching (#20024144, Illlumina). Paired-end sequencing of multiplexed transcriptome libraries was performed using a NovaSeq 6000 sequencer (Illumina) and SP100-cycle kit (#20027464, Illumina). TCR and BCR libraries were multiplexed, FAB-treated and paired-end-sequenced using a second SP100-cycle kit in a separate run.

HLA Class I Typing

HLA class I transcripts were amplified and deep-sequenced from two overlapping RT-PCR reactions flanking exons 2 (PCR 1) or 3 (PCR 2) using barcoded primers designed to target conserved regions (Supplemental Table 5) (69). Total RNA was extracted from patient PBMCs by resuspension in TRIzol reagent (Invitrogen, # 15596018), and column-purified using the PureLink RNA Mini kit (Invitrogen, #12183025). For reverse transcription, 100 pmol of oligo dT, 10 nmol of each dNTP, 40 ng RNA and sufficient nuclease-free water for a final 14 μl volume were mixed, incubated at 65°C for 5 min and chilled on ice for 5 min. This was followed by addition of 4 μL of 5X RT buffer, 40 units of RiboLock RNAse inhibitor (Thermo Fisher, #EO0381) and 200 units of Maxima H-minus reverse transcriptase (Thermo Fisher, #EP0751) and mixing. Reverse transcription was performed at 50°C for 30 min, followed by inactivation at 85°C for 5 min. 5 μl of the resulting cDNA-containing reverse transcription mixes were then used as templates for 25 μL PCR reactions using the KAPA HiFi PCR kit with GC buffer (Roche Diagnostics, #07958846001) and the following thermal cycling conditions: 95°C for 3 min; 35 cycles of 98°C for 20 s, 61°C for 15 s, 72°C for 15 s; and final extension at 72°C for 30 s. HLA amplicons were purified by gel-extraction (QIAquick gel extraction kit, Qiagen #28704) and submitted for Illumina paired-end deep-sequencing (Amplicon-EZ, Genewiz). Unique sequences originating from specific patients were identified from their respective DNA barcodes and aligned using the ClustalOmega tool to cluster sequences arising from the same allele. Sequences with the highest amount of reads in each cluster were used as input for the basic local alignment search tool (BLAST; Nucleotide collection, Homo sapiens). Sequences returning matching or highly similar alleles across PCR 1 and PCR 2 in each patient were then assembled and queried against the IMGT/HLA database for final validation.

Transcriptome scSeq Alignment and Quality Control (QC)

Reads from transcriptome scSeq (FASTQ format) were aligned to the GRCh38 reference human genome and output as filtered gene expression matrices using the 10x Genomics Cell Ranger software (version 3.1.0). Subsequent data QC and analysis was performed using R (version 3.6.2) and the Seurat package (version 3.1.5). QC steps consisted of the exclusion of TCR and BCR genes (prevention of clonotype influence on subsequent clustering), the exclusion of cells with lower than 150 or greater than 3500 genes (low quality cells), and the exclusion of cells in which more than 20% of UMIs were associated with mitochondrial genes (reduction of freeze-thaw metabolic effects) (70).

Dataset Normalisation and Integration of Multiple Datasets

Patient datasets were merged into a Seurat object list using the merge and SplitObject function. Each patient dataset was then separately normalised using SCTransform. Variable integration features (3,000) were calculated using the SelectIntegrationFeatures function from the R package Seurat (71) and setting them as variable features after merging the normalised patient datasets. Principal component analysis for dimensionality reduction was performed using the RunPCA function with up to 50 principal components. Potential batch effects between patient samples were addressed with the Harmony R package (version 1.0) using the RunHarmony function (72). Finally, unsupervised clustering was performed using the FindNeighbours and FindClusters functions. Non-linear dimensionality reduction using the RunUMAP function was performed using the first 50 principal components to generate the final UMAP visualization of cell clusters.

Dataset Subsetting of CD8+ T Cells, CD4+ T Cells and B Cells

Initial T and B cell separation was performed by mapping of TCR and BCR (VDJ) cell-specific barcodes onto the scSeq transcriptome dataset. Double attribution of TCR and BCR to the same cell (i.e., barcode) was used to identify and exclude doublets. Separation of CD8+ and CD4+ T cells was performed using the WhichCells function, from the R package Seurat, based on the singular expression of CD8A and CD4, respectively. Additional filtering of B cells was done by discarding all B cells that showed expression of CD3E or SDC1 as well as excluding B cells whose cellular barcodes occured in the Plasma cell BCR (VDJ) cell barcodes (data not shown).

Cell State Annotation and Marker Identification

The expression of specific markers in identified clusters was determined using the FindAllMarkers function using the Wilcoxon Rank Sum test. Cluster-specific markers were thresholded by having a log2(fold-change) greater than 0.25 between cells in the respective cluster and remaining cells; with marker expression occurring in at least 25% of cells in the cluster. Clusters were then attributed with specific cell states based on the expression of canonical markers.

Differential Gene Expression Analysis

Differentially expressed genes between two groups of cells were identified using the FindMarkers function. Genes were thresholded by being expressed in more than 50% of the cells and by having a log2(fold-change) greater than 0.5 between cells of the different groups using the Wilcoxon Rank Sum test.

Paired TCR and BCR (VDJ) Single-Cell Sequencing Alignment and QC

TCR and BCR reads in FASTQ format were aligned with the VDJ-GRCh38-alts-ensembl reference using the 10x Genomics Cell Ranger VDJ software (version 3.1.0). This generated single-cell VDJ sequences and annotations such as gene usage, clonotype frequency and cell-specific barcode information. As a QC step, only cells with one productive alpha and one productive beta chain (T cells) or with one productive heavy and one productive light chain (B cells)were retained for downstream analysis.

Paired TCR and BCR (VDJ) Analysis

Clonotype definition was adjusted to count all sequences as clonal if they met the following criteria: (1) Same V and J gene usage in both chains, (2) Same CDR3 length in both chains and (3) 80% amino acid sequence similarity in the CDR3 region of the TCRβ (T cells) or BCR heavy chain (B cells). Shared cellular barcode information between TCR/BCR (VDJ) scSeq and transcriptome scSeq data was used to project TCR and BCR clonotypes onto the UMAP plots (colour-coded by clonal expansion level).

Somatic Hypermutation Analysis

SHM levels in individual BCR clonotypes were determined using the change-o toolkit from the Immcantation portal as a wrapper to run IgBlast on the Cell Ranger VDJ output. The IgBlast output enabled assessment of germline similarity of single-cell BCR (VDJ) sequences. Germline identity was used as a proxy for somatic hypermutation levels and was calculated from alignments of BCR clonotypes with their corresponding VH and VL germline sequences.

TCR Specificity Group Identification Using GLIPH2

GLIPH2 clusters TCRs into specificity groups predicted to share the same antigen specificity based on sequence similarity (55). We used this algorithm to cluster TCRs from HLA-A*0201 patients (Pt-2, Pt-4 and Pt-8) together with known SARS-CoV-2 binders as well as CMV and EBV binders, which were also from HLA-A*0201 background (obtained from VDJdb database). Specificity groups that were reported by GLIPH2, were filtered for groups that were significant according to the Fisher’s Exact test (significance level < 5%) and contained at least one patient TCR and one TCR of known specificity (i.e., SARS-CoV-2, EBV or CMV). Specificity groups were identified with either global (0-1 amino acid differences in CDR3β) or local similarities (CDR3β share a common motif that is rare in the reference dataset). CD4 expressing clonotypes were filtered out.

Pseudotime Analysis

Pseudotime and trajectory inference was applied to scSeq transcriptome data using the slingshot function with default parameters from the Slingshot package in R (73). The naive cluster from each CD8+ T cell, CD4+ T cell and B cell subgroup was set as the starting point for the minimum spanning tree. The previously generated UMAP clustering was set as the cellular embedding on which Slingshot performed trajectory inference computation.

Code Availability

The data analysis pipeline followed the standard procedures as outlined in Cell Ranger and Seurat documentations. Custom scripts and functions for easier downstream analysis and visualisation purposes are available upon request.

List of Utilized R Packages

Biobase (2.46.0), BiocGenerics (0.32.0), BiocParallel (1.20.1), Cell Ranger (3.1.0), Change-O (1.0.0), circlize (0.4.10), data.table (1.12.8), DelayedArray (0.12.3), dplyr (0.8.5), GenomeInfoDb (1.22.1), GenomicRanges (1.38.0), ggplot2 (3.3.2.9000), harmony (1.0), pheatmap (1.0.12), princurve (2.1.5), RColorBrewer (1.1-2), matrixStats (0.56.0), sctransform (0.2.1), Seurat (3.1.5), slingshot (1.4.0), stringdist (0.9.5.5), stringr (1.4.0), tibble (3.0.3), tidyr (1.1.0), tidyverse (1.3.0).

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: www.ebi.ac.uk/arrayexpress/, E-MTAB-10169.

Ethics Statement

SERO-BL-COVID-19 study sponsored by the Department of Health, Canton Basel-Landschaft, Switzerland. The patients/participants provided their written informed consent to participate in this study.

Author Contributions

FB, RV-L, AY, and SR designed the study. FB, RV-L, RE, DM, BW, EK, and RB performed experiments; FB, RV-L, AY, and CW analyzed data. FB, RV-L, AY, and SR wrote the manuscript with input from all authors. All authors contributed to the article and approved the submitted version.

Funding

This study is supported by funding from the Personalized Health and Related Technologies Postdoctoral Fellowship (to R.V.-L), the NCCR Molecular Systems Engineering (to S.T.R.), Helmut Horten Stiftung (to S.T.R.), Botnar Research Centre for Child Health (to S.T.R.).

Conflict of Interest

Author DM and CW were employed by company deepCDR Biologics AG.

The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

We thank the Genomics Facility Basel (D-BSSE, ETH Zurich), in particular Ms. Ina Nissen and Dr. Christian Beisel, for their excellent support with the scSeq protocol and for performing deep sequencing runs.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2021.701085/full#supplementary-material

References

1. Dan JM, Mateus J, Kato Y, Hastie KM, Yu ED, Faliti CE, et al. Immunological Memory to SARS-CoV-2 Assessed for Up to 8 Months After Infection. Science (2021) 371:(6529):eabf4063. doi: 10.1126/science.abf4063

CrossRef Full Text | Google Scholar

2. Gaebler C, Wang Z, Lorenzi JCC, Muecksch F, Finkin S, Tokuyama M, et al. Evolution of Antibody Immunity to SARS-Cov-2. Nature (2021) 591:639–44. doi: 10.1038/s41586-021-03207-w

PubMed Abstract | CrossRef Full Text | Google Scholar

3. Grifoni A, Weiskopf D, Ramirez SI, Mateus J, Dan JM, Moderbacher CR, et al. Targets of T Cell Responses to SARS-CoV-2 Coronavirus in Humans With COVID-19 Disease and Unexposed Individuals. Cell (2020) 181:1489–501.e15. doi: 10.1016/j.cell.2020.05.015

PubMed Abstract | CrossRef Full Text | Google Scholar

4. Liu L, Wang P, Nair MS, Yu J, Rapp M, Wang Q, et al. Potent Neutralizing Antibodies Against Multiple Epitopes on SARS-CoV-2 Spike. Nature (2020) 584:450–6. doi: 10.1038/s41586-020-2571-7

PubMed Abstract | CrossRef Full Text | Google Scholar

5. Juno JA, Tan H-X, Lee WS, Reynaldi A, Kelly HG, Wragg K, et al. Humoral and Circulating Follicular Helper T Cell Responses in Recovered Patients With COVID-19. Nat Med (2020) 26:1428–34. doi: 10.1038/s41591-020-0995-0

PubMed Abstract | CrossRef Full Text | Google Scholar

6. Altmann DM. Adaptive Immunity to SARS-Cov-2. Oxford Open Immunol (2020) 1:3–4. doi: 10.1093/oxfimm/iqaa003

CrossRef Full Text | Google Scholar

7. Rees AR. Understanding the Human Antibody Repertoire. MAbs (2020) 12:1729683. doi: 10.1080/19420862.2020.1729683

PubMed Abstract | CrossRef Full Text | Google Scholar

8. Rosati E, Dowds CM, Liaskou E, Henriksen EKK, Karlsen TH, Franke A. Overview of Methodologies for T-Cell Receptor Repertoire Analysis. BMC Biotechnol (2017) 17:61. doi: 10.1186/s12896-017-0379-9

PubMed Abstract | CrossRef Full Text | Google Scholar

9. Briney B, Inderbitzin A, Joyce C, Burton DR. Commonality Despite Exceptional Diversity in the Baseline Human Antibody Repertoire. Nature (2019) 566:393–7. doi: 10.1038/s41586-019-0879-y

PubMed Abstract | CrossRef Full Text | Google Scholar

10. Soto C, Bombardi RG, Branchizio A, Kose N, Matta P, Sevy AM, et al. High Frequency of Shared Clonotypes in Human B Cell Receptor Repertoires. Nature (2019) 566:398–402. doi: 10.1038/s41586-019-0934-8

PubMed Abstract | CrossRef Full Text | Google Scholar

11. DeWitt WS, Smith A, Schoch G, Hansen JA, Matsen FA, Bradley P. Human T Cell Receptor Occurrence Patterns Encode Immune History, Genetic Background, and Receptor Specificity. eLife (2018) 7:3–5. doi: 10.7554/elife.38358

CrossRef Full Text | Google Scholar

12. Greiff V, Menzel U, Miho E, Weber C, Riedel R, Cook S, et al. Systems Analysis Reveals High Genetic and Antigen-Driven Predetermination of Antibody Repertoires Throughout B Cell Development. Cell Rep (2017) 19:1467–78. doi: 10.1016/j.celrep.2017.04.054

PubMed Abstract | CrossRef Full Text | Google Scholar

13. Emerson RO, DeWitt WS, Vignali M, Gravley J, Hu JK, Osborne EJ, et al. Immunosequencing Identifies Signatures of Cytomegalovirus Exposure History and HLA-Mediated Effects on the T Cell Repertoire. Nat Genet (2017) 49:659–65. doi: 10.1038/ng.3822

PubMed Abstract | CrossRef Full Text | Google Scholar

14. Galson JD, Trück J, Fowler A, Clutterbuck EA, Münz M, Cerundolo V, et al. Analysis of B Cell Repertoire Dynamics Following Hepatitis B Vaccination in Humans, and Enrichment of Vaccine-Specific Antibody Sequences. EBioMedicine (2015) 2:2070–9. doi: 10.1016/j.ebiom.2015.11.034

PubMed Abstract | CrossRef Full Text | Google Scholar

15. Egorov ES, Kasatskaya SA, Zubov VN, Izraelson M, Nakonechnaya TO, Staroverov DB, et al. The Changing Landscape of Naive T Cell Receptor Repertoire With Human Aging. Front Immunol (2018) 9:1618. doi: 10.3389/fimmu.2018.01618

PubMed Abstract | CrossRef Full Text | Google Scholar

16. Reddy ST, Ge X, Miklos AE, Hughes RA, Kang SH, Hoi KH, et al. Monoclonal Antibodies Isolated Without Screening by Analyzing the Variable-Gene Repertoire of Plasma Cells. Nat Biotechnol (2010) 28:965–9. doi: 10.1038/nbt.1673

PubMed Abstract | CrossRef Full Text | Google Scholar

17. Vollmers C, Sit RV, Weinstein JA, Dekker CL, Quake SR. Genetic Measurement of Memory B-Cell Recall Using Antibody Repertoire Sequencing. Proc Natl Acad Sci USA (2013) 110:13463–8. doi: 10.1073/pnas.1312146110

PubMed Abstract | CrossRef Full Text | Google Scholar

18. Niu X, Li S, Li P, Pan W, Wang Q, Feng Y, et al. Longitudinal Analysis of T and B Cell Receptor Repertoire Transcripts Reveal Dynamic Immune Response in COVID-19 Patients. Front Immunol (2020) 11:582010. doi: 10.3389/fimmu.2020.582010

PubMed Abstract | CrossRef Full Text | Google Scholar

19. Macosko EZ, Basu A, Satija R, Nemesh J, Shekhar K, Goldman M, et al. Highly Parallel Genome-Wide Expression Profiling of Individual Cells Using Nanoliter Droplets. Cell (2015) 161:1202–14. doi: 10.1016/j.cell.2015.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

20. Han X, Wang R, Zhou Y, Fei L, Sun H, Lai S, et al. Mapping the Mouse Cell Atlas by Microwell-Seq. Cell (2018) 173:1307. doi: 10.1016/j.cell.2018.05.012

PubMed Abstract | CrossRef Full Text | Google Scholar

21. van der Maaten L, van der Maaten L, Hinton G. Visualizing non-Metric Similarities in Multiple Maps. Mach Learn (2012) 87:33–55. doi: 10.1007/s10994-011-5273-4

CrossRef Full Text | Google Scholar

22. McInnes L, Healy J, Saul N, Großberger L. Umap: Uniform Manifold Approximation and Projection. J Open Source Software (2018) 3:861. doi: 10.21105/joss.00861

CrossRef Full Text | Google Scholar

23. Singer M, Wang C, Cong L, Marjanovic ND, Kowalczyk MS, Zhang H, et al. A Distinct Gene Module for Dysfunction Uncoupled From Activation in Tumor-Infiltrating T Cells. Cell (2017) 171:1221–3. doi: 10.1016/j.cell.2017.11.006

PubMed Abstract | CrossRef Full Text | Google Scholar

24. Brummelman J, Mazza EMC, Alvisi G, Colombo FS, Grilli A, Mikulak J, et al. High-Dimensional Single Cell Analysis Identifies Stem-Like Cytotoxic CD8 T Cells Infiltrating Human Tumors. J Exp Med (2018) 215:2520–35. doi: 10.1084/jem.20180684

PubMed Abstract | CrossRef Full Text | Google Scholar

25. Cohn LB, da Silva IT, Valieris R, Huang AS, Lorenzi JCC, Cohen YZ, et al. Clonal CD4 T Cells in the HIV-1 Latent Reservoir Display a Distinct Gene Profile Upon Reactivation. Nat Med (2018) 24:604–9. doi: 10.1038/s41591-018-0017-7

PubMed Abstract | CrossRef Full Text | Google Scholar

26. Chen J, Tan Y, Sun F, Hou L, Zhang C, Ge T, et al. Single-Cell Transcriptome and Antigen-Immunoglobin Analysis Reveals the Diversity of B Cells in Non-Small Cell Lung Cancer. Genome Biol (2020) 21:152. doi: 10.1186/s13059-020-02064-6

PubMed Abstract | CrossRef Full Text | Google Scholar

27. Trapnell C, Cacchiarelli D, Grimsby J, Pokharel P, Li S, Morse M, et al. The Dynamics and Regulators of Cell Fate Decisions Are Revealed by Pseudotemporal Ordering of Single Cells. Nat Biotechnol (2014) 32:381–6. doi: 10.1038/nbt.2859

PubMed Abstract | CrossRef Full Text | Google Scholar

28. Zhang J-Y, Wang X-M, Xing X, Xu Z, Zhang C, Song J-W, et al. Single-Cell Landscape of Immunological Responses in Patients With COVID-19. Nat Immunol (2020) 21:1107–18. doi: 10.1038/s41590-020-0762-x

PubMed Abstract | CrossRef Full Text | Google Scholar

29. DeKosky BJ, Kojima T, Rodin A, Charab W, Ippolito GC, Ellington AD, et al. In-Depth Determination and Analysis of the Human Paired Heavy- and Light-Chain Antibody Repertoire. Nat Med (2015) 21:86–91. doi: 10.1038/nm.3743

PubMed Abstract | CrossRef Full Text | Google Scholar

30. Stoeckius M, Hafemeister C, Stephenson W, Houck-Loomis B, Chattopadhyay PK, Swerdlow H, et al. Simultaneous Epitope and Transcriptome Measurement in Single Cells. Nat Methods (2017) 14:865–8. doi: 10.1038/nmeth.4380

PubMed Abstract | CrossRef Full Text | Google Scholar

31. Howie B, Sherwood AM, Berkebile AD, Berka J, Emerson RO, Williamson DW, et al. High-Throughput Pairing of T Cell Receptor α and β Sequences. Sci Transl Med (2015) 7:301ra131. doi: 10.1126/scitranslmed.aac5624

PubMed Abstract | CrossRef Full Text | Google Scholar

32. Setliff I, Shiakolas AR, Pilewski KA, Murji AA, Mapengo RE, Janowska K, et al. High-Throughput Mapping of B Cell Receptor Sequences to Antigen Specificity. Cell (2019) 179:1636–46.e15. doi: 10.1016/j.cell.2019.11.003

PubMed Abstract | CrossRef Full Text | Google Scholar

33. Friedensohn S, Khan TA, Reddy ST. Advanced Methodologies in High-Throughput Sequencing of Immune Repertoires. Trends Biotechnol (2017) 35:203–14. doi: 10.1016/j.tibtech.2016.09.010

PubMed Abstract | CrossRef Full Text | Google Scholar

34. Han A, Glanville J, Hansmann L, Davis MM. Corrigendum: Linking T-Cell Receptor Sequence to Functional Phenotype at the Single-Cell Level. Nat Biotechnol (2015) 33:210. doi: 10.1038/nbt0215-210c

PubMed Abstract | CrossRef Full Text | Google Scholar

35. Stubbington MJT, Lönnberg T, Proserpio V, Clare S, Speak AO, Dougan G, et al. T Cell Fate and Clonality Inference From Single-Cell Transcriptomes. Nat Methods (2016) 13:329–32. doi: 10.1038/nmeth.3800

PubMed Abstract | CrossRef Full Text | Google Scholar

36. Lönnberg T, Svensson V, James KR, Fernandez-Ruiz D, Sebina I, Montandon R, et al. Single-Cell RNA-seq and Computational Analysis Using Temporal Mixture Modelling Resolves Th1/Tfh Fate Bifurcation in Malaria. Sci Immunol (2017) 2:3–4. doi: 10.1126/sciimmunol.aal2192

CrossRef Full Text | Google Scholar

37. Xu G, Qi F, Li H, Yang Q, Wang H, Wang X, et al. The Differential Immune Responses to COVID-19 in Peripheral and Lung Revealed by Single-Cell RNA Sequencing. Cell Discov (2020) 6:73. doi: 10.1038/s41421-020-00225-2

PubMed Abstract | CrossRef Full Text | Google Scholar

38. Su Y, Chen D, Yuan D, Lausted C, Choi J, Dai CL, et al. Multi-Omics Resolves a Sharp Disease-State Shift Between Mild and Moderate Covid-19. Cell (2020) 183:1479–95. doi: 10.1016/j.cell.2020.10.037

PubMed Abstract | CrossRef Full Text | Google Scholar

39. Kaltenbach H-M, Rudolf F, Linnik J, Deichmann J, Ruf T, Altamura R, et al. Initial Characterisation of ELISA Assays and the Immune Response of the Clinically Correlated SARS-CoV-2 Biobank SERO-BL-COVID-19 Collected During the Pandemic Onset in Switzerland. doi: 10.1101/2020.07.05.20145888

CrossRef Full Text | Google Scholar

40. Davies NG, , CMMID COVID-19 working group, Klepac P, Liu Y, Prem K, Jit M, et al. Age-Dependent Effects in the Transmission and Control of COVID-19 Epidemics. Nat Med (2020) 26:1205–11. doi: 10.1038/s41591-020-0962-9

PubMed Abstract | CrossRef Full Text | Google Scholar

41. Becht E, McInnes L, Healy J, Dutertre C-A, Kwok IWH, Ng LG, et al. Dimensionality Reduction for Visualizing Single-Cell Data Using UMAP. Nat Biotechnol (2018) 37:38–44. doi: 10.1038/nbt.4314

CrossRef Full Text | Google Scholar

42. Yan J, Greer JM, Hull R, O’Sullivan JD, Henderson RD, Read SJ, et al. The Effect of Ageing on Human Lymphocyte Subsets: Comparison of Males and Females. Immun Ageing (2010) 7:4. doi: 10.1186/1742-4933-7-4

PubMed Abstract | CrossRef Full Text | Google Scholar

43. Galson JD, Schaetzle S, Bashford-Rogers RJM, Raybould MIJ, Kovaltsuk A, Kilpatrick GJ, et al. Deep Sequencing of B Cell Receptor Repertoires From Covid-19 Patients Reveals Strong Convergent Immune Signatures. Front Immunol (2020) 11:605170. doi: 10.3389/fimmu.2020.605170

PubMed Abstract | CrossRef Full Text | Google Scholar

44. DeKosky BJ, Ippolito GC, Deschner RP, Lavinder JJ, Wine Y, Rawlings BM, et al. High-Throughput Sequencing of the Paired Human Immunoglobulin Heavy and Light Chain Repertoire. Nat Biotechnol (2013) 31:166–9. doi: 10.1038/nbt.2492

PubMed Abstract | CrossRef Full Text | Google Scholar

45. Wabl M, Steinberg C. Affinity Maturation and Class Switching. Curr Opin Immunol (1996) 8:89–92. doi: 10.1016/s0952-7915(96)80110-5

PubMed Abstract | CrossRef Full Text | Google Scholar

46. Rosa FM, Fellous M. Regulation of HLA-DR Gene by IFN-Gamma. Transcriptional and Post-Transcriptional Control. J Immunol (1988) 140:1660–4.

PubMed Abstract | Google Scholar

47. Upadhyay V, Fu Y-X. Lymphotoxin Signalling in Immune Homeostasis and the Control of Microorganisms. Nat Rev Immunol (2013) 13:270–9. doi: 10.1038/nri3406

PubMed Abstract | CrossRef Full Text | Google Scholar

48. Nowacki TM, Kuerten S, Zhang W, Shive CL, Kreher CR, Boehm BO, et al. Granzyme B Production Distinguishes Recently Activated CD8(+) Memory Cells From Resting Memory Cells. Cell Immunol (2007) 247:36–48. doi: 10.1016/j.cellimm.2007.07.004

PubMed Abstract | CrossRef Full Text | Google Scholar

49. Wang M, Windgassen D, Papoutsakis ET. A Global Transcriptional View of Apoptosis in Human T-Cell Activation. BMC Med Genomics (2008) 1:53. doi: 10.1186/1755-8794-1-53

PubMed Abstract | CrossRef Full Text | Google Scholar

50. Yamada T, Gierach K, Lee P-H, Wang X, Daniel Lacorazza H. Cutting Edge: Expression of the Transcription Factor E74-Like Factor 4 Is Regulated by the Mammalian Target of Rapamycin Pathway in CD8 T Cells. J Immunol (2010) 185:3824–8. doi: 10.4049/jimmunol.1000718

PubMed Abstract | CrossRef Full Text | Google Scholar

51. Tsui C, Maldonado P, Montaner B, Borroto A, Alarcon B, Bruckbauer A, et al. Dynamic Reorganisation of Intermediate Filaments Coordinates Early B-Cell Activation. Life Sci Alliance (2018) 1:e201800060. doi: 10.26508/lsa.201800060

PubMed Abstract | CrossRef Full Text | Google Scholar

52. Shaffer AL, Shapiro-Shelef M, Iwakoshi NN, Lee A-H, Qian S-B, Zhao H, et al. XBP1, Downstream of Blimp-1, Expands the Secretory Apparatus and Other Organelles, and Increases Protein Synthesis in Plasma Cell Differentiation. Immunity (2004) 21:81–93. doi: 10.1016/j.immuni.2004.06.010

PubMed Abstract | CrossRef Full Text | Google Scholar

53. Gass JN, Gunn KE, Sriburi R, Brewer JW. Stressed-Out B Cells? Plasma-Cell Differentiation and the Unfolded Protein Response. Trends Immunol (2004) 25:17–24. doi: 10.1016/j.it.2003.11.004

PubMed Abstract | CrossRef Full Text | Google Scholar

54. Shimoda M, Li T, Pihkala JPS, Koni PA. Role of MHC Class II on Memory B Cells in Post-Germinal Center B Cell Homeostasis and Memory Response. J Immunol (2006) 176:2122–33. doi: 10.4049/jimmunol.176.4.2122

PubMed Abstract | CrossRef Full Text | Google Scholar

55. Huang H, Wang C, Rubelt F, Scriba TJ, Davis MM. Analyzing the Mycobacterium Tuberculosis Immune Response by T-Cell Receptor Clustering With GLIPH2 and Genome-Wide Antigen Screening. Nat Biotechnol (2020) 38:1194–202. doi: 10.1038/s41587-020-0505-4

PubMed Abstract | CrossRef Full Text | Google Scholar

56. Liao M, Liu Y, Yuan J, Wen Y, Xu G, Zhao J, et al. Single-Cell Landscape of Bronchoalveolar Immune Cells in Patients With COVID-19. Nat Med (2020) 26:842–4. doi: 10.1038/s41591-020-0901-9

PubMed Abstract | CrossRef Full Text | Google Scholar

57. Peng Y, Mentzer AJ, Liu G, Yao X, Yin Z, Dong D, et al. Broad and Strong Memory CD4 and CD8 T Cells Induced by SARS-CoV-2 in UK Convalescent Individuals Following COVID-19. Nat Immunol (2020) 21:1336–45. doi: 10.1038/s41590-020-0782-6

PubMed Abstract | CrossRef Full Text | Google Scholar

58. Sekine T, Perez-Potti A, Rivera-Ballesteros O, Strålin K, Gorin J-B, Olsson A, et al. Robust T Cell Immunity in Convalescent Individuals With Asymptomatic or Mild Covid-19. Cell (2020) 183:158–68.e14. doi: 10.1016/j.cell.2020.08.017

PubMed Abstract | CrossRef Full Text | Google Scholar

59. Shomuradova AS, Vagida MS, Sheetikov SA, Zornikova KV, Kiryukhin D, Titov A, et al. Sars-CoV-2 Epitopes Are Recognized by a Public and Diverse Repertoire of Human T Cell Receptors. Immunity (2020) 53(6):1245–57.e5. doi: 10.1016/j.immuni.2020.11.004

PubMed Abstract | CrossRef Full Text | Google Scholar

60. Shoukat MS, Foers AD, Woodmansey S, Evans SC, Fowler A, Soilleux E. Use of Machine Learning to Identify a T Cell Response to SARS-Cov-2. Cell Rep Med (2021) 2(2):100192.. doi: 10.1016/j.xcrm.2021.100192

PubMed Abstract | CrossRef Full Text | Google Scholar

61. Ferretti AP, Kula T, Wang Y, Nguyen DMV, Weinheimer A, Dunlap GS, et al. Unbiased Screens Show Cd8 T Cells of COVID-19 Patients Recognize Shared Epitopes in SARS-CoV-2 That Largely Reside Outside the Spike Protein. Immunity (2020) 53:1095–107.e3. doi: 10.1016/j.immuni.2020.10.006

PubMed Abstract | CrossRef Full Text | Google Scholar

62. Tarke A, Sidney J, Kidd CK, Dan JM, Ramirez SI, Yu ED, et al. Comprehensive Analysis of T Cell Immunodominance and Immunoprevalence of SARS-CoV-2 Epitopes in COVID-19 Cases. Cell Rep Med (2021) 2:(2). doi: 10.1016/j.xcrm.2021.100204

PubMed Abstract | CrossRef Full Text | Google Scholar

63. Vidarsson G, Dekkers G, Rispens T. Igg Subclasses and Allotypes: From Structure to Effector Functions. Front Immunol (2014) 5:520. doi: 10.3389/fimmu.2014.00520

PubMed Abstract | CrossRef Full Text | Google Scholar

64. Suan D, Sundling C, Brink R. Plasma Cell and Memory B Cell Differentiation From the Germinal Center. Curr Opin Immunol (2017) 45:97–102. doi: 10.1016/j.coi.2017.03.006

PubMed Abstract | CrossRef Full Text | Google Scholar

65. Laidlaw BJ, Cyster JG. Transcriptional Regulation of Memory B Cell Differentiation. Nat Rev Immunol (2020) 21:209–20. doi: 10.1038/s41577-020-00446-2

PubMed Abstract | CrossRef Full Text | Google Scholar

66. Robbiani DF, Gaebler C, Muecksch F, Lorenzi JCC, Wang Z, Cho A, et al. Convergent Antibody Responses to SARS-CoV-2 in Convalescent Individuals. Nature (2020) 584:437–42. doi: 0.1038/s41586-020-2456-9

PubMed Abstract | Google Scholar

67. Rogers TF, Zhao F, Huang D, Beutler N, Burns A, He W-T, et al. Isolation of Potent SARS-CoV-2 Neutralizing Antibodies and Protection From Disease in a Small Animal Model. Science (2020) 369:956. doi: 10.1126/science.abc7520

PubMed Abstract | CrossRef Full Text | Google Scholar

68. Brouwer PJM, Caniels TG, van der Straten K, Snitselaar JL, Aldon Y, Bangaru S, et al. Potent Neutralizing Antibodies From COVID-19 Patients Define Multiple Targets of Vulnerability. Science (2020) 369:643–50. doi: 10.1101/2020.05.12.088716

PubMed Abstract | CrossRef Full Text | Google Scholar

69. Lank SM, Golbach BA, Creager HM, Wiseman RW, Keskin DB, Reinherz EL, et al. Ultra-High Resolution HLA Genotyping and Allele Discovery by Highly Multiplexed cDNA Amplicon Pyrosequencing. BMC Genomics (2012) 13:378. doi: 10.1186/1471-2164-13-378

PubMed Abstract | CrossRef Full Text | Google Scholar

70. Luecken MD, Theis FJ. Current Best Practices in Single-Cell RNA-seq Analysis: A Tutorial. Mol Syst Biol (2019) 15:e8746. doi: 10.15252/msb.20188746

PubMed Abstract | CrossRef Full Text | Google Scholar

71. Stuart T, Butler A, Hoffman P, Hafemeister C, Papalexi E, Mauck WM 3rd, et al. Comprehensive Integration of Single-Cell Data. Cell (2019) 177:1888–902.e21. doi: 10.1016/j.cell.2019.05.031.

PubMed Abstract | CrossRef Full Text | Google Scholar

72. Korsunsky I, Millard N, Fan J, Slowikowski K, Zhang F, Wei K, et al. Fast, Sensitive and Accurate Integration of Single-Cell Data With Harmony. Nat Methods (2019) 16:1289–96. doi: 10.1038/s41592-019-0619-0

PubMed Abstract | CrossRef Full Text | Google Scholar

73. Street K, Risso D, Fletcher RB, Das D, Ngai J, Yosef N, et al. Slingshot: Cell Lineage and Pseudotime Inference for Single-Cell Transcriptomics. BMC Genomics (2018) 19:477. doi: 10.1186/s12864-018-4772-0

PubMed Abstract | CrossRef Full Text | Google Scholar