Identifying “immunotypes” and relationships between circulating B and T cell responses with disease severity in COVID-19 patients
To further investigate the relationship between immune responses and COVID-19 disease trajectory, we stratified the COVID-19 patients (n = 125) into eight different categories according to the NIH Ordinal Severity Scale ranging from COVID 1 (death) and COVID 2 (requiring maximal clinical intervention) to COVID 8 (at home with no required care) (Fig. 6A). We then asked how changes in T and B cell populations defined above on D0 were related to disease severity. More severe disease was associated with lower frequencies of CD8 and CD4 T cells, with a greater effect on CD8 T cells in less severe disease (Fig. 6B). Taking all patients together, there were no statistically significant changes in the major T cell and B cell subsets related to disease severity though some trends were present (fig. S7, A to C). In contrast, HLA-DR+CD38+ CD8 T cells as well as both KI67+ and HLA-DR+CD38+ CD4 T cells were increased in patients with more severe disease (fig. S7, D and E).
Fig. 6  High dimensional analysis of immune phenotypes with clinical data reveals distinct COVID-19 patient immunotypes.
(A) NIH ordinal scale for COVID-19 clinical severity. (B) Frequencies of major immune subsets. Significance determined by unpaired Wilcoxon test with BH correction: *p < 0.05, **p < 0.01, ***p < 0.001, and ****p < 0.0001. (C) Heatmap of indicated immune parameters by row; donor type, disease severity, and mortality indicated across top. (D) UMAP projection of aggregated flow cytometry data. (E) Transformed UMAP projection; density contours drawn separately for healthy donor (HD), recovered donor (RD), and COVID-19 subjects (see Methods). (F) Bars represent mean of UMAP Component 1. Dots represent individual subjects; bars shaded by subject group and/or severity score. (G) Density contour plots indicating variation of specified immune features across UMAP Component coordinates. Relative expression (according to heat scale) shown for both individual patients (points) and overall density (contours). Spearman’s Rank Correlation coefficient (ρ) and p-value for each feature vs. Component 1 (C1) and Component 2 (C2) shown. (H) (Left) Spearman correlation between UMAP Components 1 and 2 and FlowSOM clusters. (Right) Select FlowSOM clusters and their protein expression. (I) Spearman correlation between UMAP Components 1 and 2 and clinical metadata. (J) Heatmap of immune parameters used to define Immunotype 3 indicated by row; disease severity and mortality indicated across top. (K) (Left) Transformed UMAP projection; patient status for Immunotype 3 indicated by color. (Right) Spearman correlation between Immunotype 3 and disease severity, mortality, and UMAP Components. There were two challenges with extracting meaning from these data. First, there was considerable inter-patient heterogeneity for each of these immune features related to disease severity score. Second, these binary comparisons (e.g., one immune subset versus one clinical feature) vastly underutilized the high dimensional information in this dataset. Thus, we next visualized major T and B cell subpopulation data as it related to clinical disease severity score (Fig. 6C). Data were clustered based on immune features and then overlaid with the disease severity score over time for each patient. This analysis revealed groups of patients with similar composite immune signatures of T and B cell populations (Fig. 6C). When individual CD8 T cell, CD4 T cell, or B cell populations were examined, a similar concept of patient subgroups emerged (fig. S7, F, G, and H). These data suggested the idea of “immunotypes” of COVID-19 patients based on integrated responses of T and B cells, though some individual cell types and/or phenotypes separated patients more clearly than others.
These approaches provided insight into potential immune phenotypes associated with patients with severe disease, but suffered from the use of a small number of manually selected T or B cell subsets or phenotypes. We therefore next employed Uniform Manifold Approximation and Projection (UMAP) to distill the ~200 flow cytometry features (see tables S5 and S6) representing the immune landscape of COVID-19 disease in two dimensional space, creating compact meta features (or Components) that could then be correlated with clinical outcomes. This analysis revealed a clear trajectory from HD to COVID-19 patients (Fig. 6D), which we centered and aligned with the horizontal axis (“Component 1”) to facilitate downstream analysis (Fig. 6E). An orthogonal vertical axis coordinate (“Component 2”) also existed that captured non-overlapping aspects of the immune landscape. We next calculated the mean of Component 1 for each patient group, with COVID-19 patients separated by severity (Fig. 6E). The contribution of Component 1 clearly increased in a stepwise manner with increasing disease severity (Fig. 6F). Interestingly, RD were subtly positioned between HD and COVID-19 patients. Component 1 remained an independent predictor of disease severity (P = 5.5 × 10−5) even after adjusting for the confounding demographic factors of age, sex, and race.
We next investigated how the UMAP Components were associated with individual immune features (tables S5 and S6). UMAP Component 1 captured immune features, including the relative loss of CD4 and CD8 T cells and increase in nonB:nonT cells (Fig. 6G). PB also associated with Component 1 (Fig. 6G). Other individual B cell features were differentially captured by UMAP Component 1 and 2. Component 1 contained a signal for T-bet+ PB populations (table S5) whereas Component 2 was enriched for T-bet+ memory B cells and CD138+ PB populations (table S6). Activated HLA-DR+CD38+ and KI67+ CD4 and CD8 T cells had contributions to both Component 1 and Component 2, with these features residing in the upper right corner of the UMAP plot (Fig. 6, G and H, and fig. S8, A to D). In contrast, T-bet+ non-naïve CD8 T cells were strongly associated with Component 2 whereas T-bet+ non-naïve CD4 T cells were also linked to Component 1 (Fig. 6G and tables S5 and S6). Eomes+ CD8 or CD4 T cells were both associated with Component 2 and negatively associated with Component 1 (Fig. 6G and tables S5 and S6).
We next took advantage of the FlowSOM clustering in Figs. 2 to 4 that identified individual immune cell types most perturbed in COVID-19 patients and linked these FlowSOM clusters to UMAP Components (Fig. 6H). For non-naïve CD8 T cells, FlowSOM Cluster 11 that contained T-bet+CX3CR1+ but non-proliferating effector-like cells was positively correlated with UMAP Component 2 and negatively correlated with Component 1 (Fig. 6H). In contrast, FlowSOM Cluster 14 contained activated, proliferating PD-1+CD39+ cells that might reflect either recently generated effector or possibly exhausted CD8 T cells (50) and was strongly associated with UMAP Component 1 (Fig. 6H). For CD4 T cells, FlowSOM Cluster 14, containing activated, proliferating CD4 T cells, was captured by both UMAP Components, whereas a second activated CD4 T cell population that also expressed CD95 (FlowSOM Cluster 13) was only captured by UMAP Component 1 (Fig. 6H). In addition, Component 1 was negatively correlated with CD4 T cell FlowSOM Clusters 2 and 3 that contained cTfh (Fig. 6H). Finally, for B cells, the FlowSOM Cluster of T-bet+CD138+ PB (Cluster 5) was positively correlated with Component 1 whereas the Tbet-CD138+ Cluster 3 was negatively correlated with Component 1 (Fig. 6H). Locations in the UMAP immune landscape were dynamic, changing from D0 to D7 for both Component 1 and Component 2 consistent with the data in Fig. 5 and fig. S9, A to F. The most dynamic changes in Component 1 were associated with the largest increases in IgM antibody levels (fig. S9G).
Given the association of the UMAP Component 1 with disease severity, we next examined the connections between UMAP Components with individual clinical features. UMAP Component 1 correlated with several clinical measurements of inflammation (e.g., ferritin, hsCRP, IL-6), co-infection, organ failure (APACHE III), and acute kidney disease and renal insufficiency (Fig. 6I). It was interesting, however, that, although D-dimer was elevated, this feature did not correlate with UMAP Component 1, but coagulation complication did (Fig. 6I). Several antibody features also correlated with Component 1 consistent with some of the immune features discussed above. In contrast, Component 2 lacked positive correlation to many of these clinical features of disease and rather was negatively correlated only to eosinophil count, NSAID use, and subsequent treatment with Remdesivir (Fig. 6I). UMAP Component 1, but not Component 2, also correlated with mortality, although there were clearly patients with high Component 2, but low Component 1 who succumbed to COVID-19 disease (Fig. 6E). These data indicate that the immune features captured by UMAP Component 1 have a strong relationship to many features of disease severity, whereas other features of immune dynamics during COVID-19 disease captured by UMAP Component 2 have a distinct relationship with clinical disease presentation.
More positive values in UMAP Components 1 or 2 captured mainly signals of change or differences in individual immune features in COVID-19 disease compared to HD and RD. UMAP Component 1 captured an immunotype (Immunotype 1) that was characterized by effector or highly activated CD4 T cells, low cTfh, some CD8 TEMRA-like activation, possibly hyperactivated CD8 T cells, and Tbet+ PB, whereas Component 2 or Immunotype 2 captured Tbetbright effector-like CD8 T cells, lacked some of the robust CD4 T cell activation but has some features of proliferating B cells (Fig. 6G and fig. S8). However, the data presented in Figs. 1 to 5 also suggested a subset of patients with minimal activation of T and B cell responses. To investigate this immune signature, we identified 20 patients who had responses more similar to HD and RD for five activated/responding B and T cell populations (Fig. 6J, middle, and fig. S10). If the UMAP Components 1 and 2 captured two distinct “immunotypes” of patient responses to SARS-CoV2 infection, this group of 20 patients represent a third immunotype. Immunotype 3 was negatively associated with UMAP Components 1 and 2 and negatively associated with disease severity, suggesting that a less robust immune response during COVID-19 was associated with less severe pathology (Fig. 6K and fig. S10), despite the fact that these patients were hospitalized with COVID-19 disease. These data further emphasize the different ways patients can present and possibly succumb to COVID-19 disease. These patterns may be related to pre-existing conditions in combination with immune response characteristics. It is likely that additional immune features, such as comprehensive serum cytokine measurements, will improve this model. Nevertheless, the current computational approach integrating deep immune profiling with disease severity trajectory and other clinical information revealed distinct patient immunotypes linked to distinct clinical outcomes (fig. S11).