Transcriptional Basis of Mouse and Human Dendritic Cell Heterogeneity

Abstract

Dendritic cells (DCs) play a critical role in orchestrating adaptive immune responses due to their unique ability to initiate T cell responses and direct their differentiation into effector lineages. Classical DCs have been divided into two subsets, cDC1 and cDC2, based on phenotypic markers and their distinct abilities to prime CD8 and CD4 T cells. While the transcriptional regulation of the cDC1 subset has been well characterized, cDC2 development and function remain poorly understood. By combining transcriptional and chromatin analyses with genetic reporter expression, we identified two principal cDC2 lineages defined by distinct developmental pathways and transcriptional regulators, including T-bet and RORγt, two key transcription factors known to define innate and adaptive lymphocyte subsets. These novel cDC2 lineages were characterized by distinct metabolic and functional programs. Extending our findings to humans revealed conserved DC heterogeneity and the presence of the newly defined cDC2 subsets in human cancer.

Keywords: ATAC-sequencing; T-bet; dendritic cells; myeloid cells; single-cell RNA-sequencing; transcriptional regulation.

Graphical abstract

Figure 1

Single-Cell Survey Reveals Heterogeneity of cDC2s with Two Subsets Delineated by Expression of T-Bet (A) Representative contour plot showing gating strategy for splenic DCs in Tbx21^RFP-Cre mice. DCs defined as Lin(CD3,CD19,CD49b,Siglec-F)^–Ly6C^–CD64^–CD11c⁺MHCII⁺. (B) Frequency of T-bet⁺ cDC2s across tissues. Each circle represents one mouse. In the peripheral and mesenteric LN (PLN and MLN), migratory DCs were defined as MHCII^hiCD11c^int and resident DCs as MHCII^intCD11c^hi. Error bars represent mean ± SEM. (C) Analysis of RFP⁺ and YFP⁺ splenic cDC2s from Tbx21^RFP-CreERT2Rosa26^YFP mice, 3 days post tamoxifen gavage. (D) Percent RFP⁺ and YFP⁺ of cDC2 cells. Percent RFP⁺ of YFP⁺ cDC2s at indicated time points post tamoxifen gavage (right). Error bars represent mean ± SEM; n = 3–4 mice per time point. (E) t-SNE embedding of 4,464 DCs. Colors indicate unsupervised clustering by Phenograph (left panel) or classification based on expression of canonical markers (right panel). (F) Expression of canonical DC markers across the transcriptionally defined DC clusters from (E). (G) Proportion of T-bet (RFP⁺) cells in each cell cluster identified in (D). (H) Violin plot showing expression of the cell-cycle signature across the DC clusters from (E). (I) Similarity of bulk T-bet^– cDC2s, T-bet⁺ cDC2, and cDC1 transcriptomes to the reference single-cell DC clusters (E). Colors represent the correlation coefficient between the cell population identified in the row label and the DC cluster identified by the column label. See also Figures S1 and S7.

Figure S1

Single-Cell Survey Reveals Heterogeneity of cDC2s, Related to Figure 1 A. Representative histogram showing expression of T-bet (RFP) in splenic cells from Tbx21^RFP-cre mice. (B). Expression of T-bet in CD11b⁺XCR1⁺ DCs from the intestinal lamina propria. Data representative of > 5 independent experiments, with at least 3 mice per experiment. (C). Expression of T-bet in splenic myeloid cells. Cells were defined as: (i) Ly-6C^hi monocytes (Lin ^–Ly6C⁺Ly6G^–CD11b⁺CX3CR1⁺); neutrophils (Lin^–Ly6C⁺Ly6G⁺); macrophages (Lin^–CD64⁺Ly6C^–). Lineages (Lin) were defined as: CD3e, CD90.2, CD19, CD49b and Siglec F. Each circle represents an individual mouse, error bars represent mean ± SEM. (D). Left: Gating strategy for single-cell sorting. DCs were defined as Lin(CD3, CD19, CD90)^–Ly6C^–CD64^–CD11c⁺MHCII⁺. Two populations were sampled: RFP⁺ DCs and RFP^– DCs (encompassing XCR1⁺ cDC1s, CD11b⁺RFP^– and CD11b^–XCR1^– DCs). Right: Post-sort purity of RFP⁺ and RFP^– cells. Contaminating population of Ly6C⁺ cells identifiable on post-sort purity (lower panel). (E). Similarity of splenic CD11c⁺MHCII⁺ cells to reference myeloid cells (ImmGen Consortium) Colors represent the Pearson correlation between the mean gene expression from the dendritic cell cluster in the rows and the bulk reference transcriptome in the columns. (F). Top 20 positive and negative gene loadings of PC1 for T-bet⁺ cDC2 clusters after cell-cycle correction (left panel). Scatterplot of PC1 and PC2 for T-bet⁺ cDC2 clusters after cell-cycle correction (right panel).

Figure 2

T-bet and RORγt Expression Delineates Distinct cDC2 Subsets (A) Heatmap reporting scaled, imputed expression for the top 10 differentially expressed genes for each cDC2 cluster identified in Figure 1E. Genes of interest are shown on the right. (B) Expression of Rorc transcript in bulk sorted T-bet^– cDC2s, T-bet⁺ cDC2s, and cDC1s. (C) Expression of CLEC12A, Esam, and CLEC10A in cDC2s from Tbx21^RFP-Cre mice. Far right: CLEC12A and CLEC10A expression within T-bet^– cDC2s. (D) Summary graph for (C). Each symbol represents one mouse. (E) Enrichment of GO pathways in in T-bet⁺ versus T-bet^– cDC2 clusters. (F) Correlation between splenic and MLN cDC2 transcriptomes. Error bars represent mean ± SEM. See also Figures S2 and S7 and Tables S1, S2, and S3.

Figure S2

Distinct cDC2 Subsets Express Divergent Transcription Factors and Cell Surface Markers. Related to Figure 2 (A). Enrichment of GO pathways (biological processes) in cluster 14. (B). Graph showing AUC score (x axis) for genes differentially expressed between T-bet⁺ and T-bet^– cDC2 clusters. Earth movers distance (EMD) on the y axis. Dashed lines represent μ_EMD ± 3σ_EMD. (C). Graph showing AUC score (x axis) for genes differentially expressed between T-bet^– cluster 10 and all other cDC2 clusters. EMD on the y axis. Dashed lines represent μ_EMD ± 3σ_EMD. (D). t-SNE map of splenic DCs colored by imputed expression of Clec10a and Mgl2 demonstrating co-expression of these genes by cells in cluster 10. (E). Representative flow cytometry plots showing expression of CLEC12A, CLEC10A and T-bet (RFP) in cDC2s isolated from spleen, peripheral lymph nodes (PLN), mesenteric lymph nodes (MLN), large intestine lamina propria (LI), liver and lung (left). Frequency of CLEC12A⁺ and CLEC10A⁺ cells within all cDC2s (right). Each symbol represents a single mouse. Data representative of 2 independent experiments. (E). Cytospin analysis of Wright-Giemsa stained sorted cDC2 subsets.

Figure 3

T-bet^– cDC2s Are Distinct from Monocytes (A) Expression of Zbtb46 in bulk sorted DC subsets. (B) Representative plot showing expression of CSF1R in cDC2s and summary graph (right). Each symbol represents one mouse. (C) Representative plot showing expression of CX3CR1 in cDC2s and summary graph (right). Each symbol represents one mouse. (D) Expression of CCR2 (GFP) in splenic DCs. T-bet^– cDC2s identified as CLEC12A⁺Esam^– and T-bet⁺ cDC2s as CLEC12A⁺Esam^– cDC2s. Monocytes gated as Lin^–Ly6C⁺CD11b⁺CX3CR1⁺. Each symbol represents one mouse. (E) Frequency of monocytes and indicated cDC2 subset in CCR2^–/– mice and CCR2^+/– littermates. Error bars represent mean ± SEM; p values calculated using Student’s t test.

Figure 4

Transcriptional and Epigenetic Landscape of cDCs (A) Heatmap showing differentially accessible ATAC-seq peaks in cDC1s, T-bet⁺ cDC2s, and T-bet^– cDC2s. Color bar, accessibility Z score. (B) ATAC-seq analysis of chromatin accessible sites in cDC2s. Peaks shown in purple showed increased read count in T-bet⁺ cDC2s, peaks shown in red were enriched in T-bet^– cDC2s. (C) Predictive value of TF motifs in peaks more accessible in T-bet^– cDC2s (red) or T-bet⁺ cDC2s (green). (D) Heatmap reporting scaled expression of predicted transcriptional regulators in T-bet^– cDCs versus T-bet⁺ cDC2s. Color bar, accessibility Z score. (E) Correlation between average differential accessibility of peaks associated with a gene and gene expression in T-bet⁺ versus T-bet^– cDC2s. (F) Diamond plot showing gains and losses of regulatory elements for top 20 most differentially expressed genes in T-bet⁺ versus T-bet^– cDC2s. Each diamond represents a chromatin accessibility peak associated with the indicated gene. Green denotes ATAC-seq peaks that gained accessibility in T-bet⁺ cDC2s, red diamonds denote peaks that lost accessibility. The bottom-most peak on the y axis corresponds to the log₂FC in differential expression of the gene. See also Figure S3 and Table S4.

Figure S3

Transcriptional and Epigenetic Landscape of cDCs, Related to Figure 4 (A). Principal component analysis (PCA) for ATAC-seq of cDC1s, T-bet⁺ cDC2s and T-bet^– cDC2s. Each symbol represents a biological replicate for each cell-type. (B). Predictive value of TF motifs in peaks more accessible in cDC2s versus cDC1s. (C). Correlation between differential ATAC-seq peak accessibility and gene expression in cDC1s versus cDC2s.

Figure S4

Environmental Cues Drive Distinct DC2 Differentiation Pathways within the Spleen, Related to Figure 5 (A). Gating strategy for the identification of DC progenitors in the bone marrow (BM) (B). Palantir pseudo-time analysis of differentiation potential and branch probabilities from the Siglec-H⁺ pre-DC state to T-bet⁺ cDC2 and T-bet^– cDC2 terminal states. (C). Plots showing Palantir differentiation potential (y axis) along Palantir pseudo-time (x axis) for Siglec-H⁺ DC and T-bet⁺ cDC2s (top) or Siglec-H⁺ DC and T-bet^– cDC2 clusters (bottom) (D). Plots showing the top two diffusion component embeddings for Siglec-H⁺ DC and T-bet⁺ cDC2 clusters (top) or Siglec-H⁺ DC and Tbet^– cDC2 clusters (bottom). Black arrow indicates Siglec-H⁺ DC cluster cells adjacent to cells from the proliferative T-bet⁺ cDC2 clusters 6 and 8. (E). Top panel: plots showing probability of each cell being within 20 nearest neighbors of randomly sampled shortest paths from the Siglec-H⁺ DC to the indicated end points. Middle panel: plots showing the proportion of cells belonging to Siglec-H⁺ DC, T-bet⁺ cDC2, or T-bet^– cDC2 from 20 nearest neighbors of randomly sampled shortest paths. Bottom: plots showing diffusion distance step sizes for each step along the indicated shortest paths (bottom panel). Colors illustrate cluster membership. (F). Graph showing AUC (x axis) for genes differentially expressed between Siglec-H⁺ DC cluster (cluster 11) and all other cDC2 clusters. EMD on the y axis. Dashed lines represents μ_EMD ± 3σ_EMD. (G). Gating strategy for FACS-isolation of MHCII⁺ ILC3s: Lin = CD3, CD19, CD49b, Siglec-F. (H). Heatmap reports scaled expression of 3550 differentially expressed genes (log₂FC > 1, FDR < 0.01) between ILC3s and Rorγt fm cDC2s. Selected genes listed to the right. (I). Representative flow cytometric analysis of phenotypes of splenic progeny from Tbx21^RFP-cre CD45.2⁺Ly6C⁻CD64^–MHCII⁺CD11c⁺Siglec-H⁺ pre-DCs adoptively transferred into sub-lethally irradiated CD45.1 recipient mice 7 days earlier (data from one experiment with n = 3). J. Sort purified T-bet⁺ or T-bet^– cDC2 were cultured for 24hrs in the presence of LPS, CpG, TNF-α or IFN−γ. Representative overlay histogram showing the expression of RFP(T-bet) at 24hrs. Data representative of 2 (TNF-α) or 4 (all other cytokines/TLR agonists) independent experiments, n = 2-3.

Figure 5

Environmental Cues Drive Distinct DC2 Differentiation Pathways within the Spleen (A) Derivation of splenic CD11b⁺ DC subsets from bone marrow (BM) progenitors. Shown are splenic cDC2s 7 days post transfer. Data representative of MDP and CDP recipients (n = 3) or pre-DC recipients (n = 6) from 3 independent experiments. (B) Analysis of YFP in BM progenitors from Rorc^CreRosa26^lsl-YFP mice (top) or RFP expression in cells from Tbx21^RFP-cre mice (bottom). (C) t-SNE embedding of diffusion map of T-bet⁺, Tbet^– cDC2, and Siglec-H⁺ pre-DC clusters identified in Figure 1E. (D) Palantir branch probabilities from Siglec-H⁺ pre-DC to T-bet⁺ cDC2 (top) or Tbet^– cDC2 (bottom) terminal states. Gene expression trends in pseudo-time along corresponding trajectories. (E) Expression of genes identified as varying significantly along the trajectory from the pre-DC cluster to T-bet⁺ or Tbet^– cDC2. (F) Expression of Siglec-H in Lin^–CD90^–Ly6C^–CD64^–CD11c⁺MHCII⁺ cells. Far right: overlay of Siglec-H⁺ DCs against all MHCII⁺CD11c⁺ DCs. (G) Representative plot showing percentage of YFP⁺ cells in splenic cDC1 and cDC2 populations in Rorc^CreRosa26^lsl-YFP mice. (H) Heatmap showing expression of genes differentially expressed between bulk T-bet⁺ and Tbet^– cDC2s across T-bet⁺, Tbet^–, and Rorγt fate-mapped (fm) cDC2s. (I) Frequency of T-bet⁺ cDC2s within the spleen, MLN, large intestine lamina propria (LI), and small intestine lamina propria (SI) at indicated time point post birth. n = 3 per time point, error bars represent mean ± SEM. For the day 7 analysis SI and LI tissues were pooled from 2 mice per sample. (J) Generation of T-bet⁺ cDC2s in mice treated with a broad-spectrum antibiotic cocktail (AVKM) or H₂O (control). Each symbol represents one mouse. Error bars represent mean ± SEM. See also Figure S4 and Table S5.

Figure 6

T-bet⁺ and T-bet^– cDC2s Are Phenotypically and Functionally Distinct (A) Heatmap of select TLRs, chemokines, cytokines, and their receptors for genes differentially expressed between bulk cDC1s, T-bet⁺, and T-bet^– cDC2s. (B) Cytokines detected in culture supernatant 18 h after stimulation with R848, analyzed using a multiplexed cytokine assay. Data are mean ± SEM from triplicate culture wells. (C) Cell surface expression of CD86, PDL1, and MHC class II on indicated DC subset 16 h after intraperitoneal (i.p.) immunization with LPS. Left: representative overlay histogram; right: composite bar graphs of median fluorescence intensity (MFI) (n = 3). (D) Sorted naive OTII CD4⁺ T cells were cultured with indicated DC subset and OVA peptide under non-polarizing (Th0) or polarizing conditions for 4 days. Bar graphs of summary data for intracellular cytokine production following restimulation or expression of Thy1.1(Foxp3) in unstimulated cells. Data shown as mean ± SEM; p values were calculated using two-way ANOVA (Th0) or one way ANOVA. See also Figure S5.

Figure S5

T-bet⁺ and T-bet^– cDC2s Are Functionally Distinct, Related to Figure 6 (A). Cytokines detected in culture supernatant 18 hours after stimulation with CpG using a multiplexed cytokine assay. (B). Proliferation of naive OTII CD4⁺ T cells 5 days after co-culture with OVA peptide and either T-bet⁺ cDC2s, T-bet^– cDC2s or cDC1s. Data, shown as mean ± SEM, are representative of 2 independent experiments (n = 3).

Figure 7

Conservation of DC Subsets across Species (A) t-SNE map of human peripheral blood DCs from Villani et al. (2017), colored by cell type or log-transformed expression of labeled genes (middle and right panel). (B) Expression of CD1c and CLEC10A within peripheral blood cDC2s. Bar graph shows summary frequencies for four individual donors. (C) Distribution of expression changes between T-bet⁺ and T-bet^– cDC2s. Genes up- or downregulated in human “DC2” versus “DC3” cDC2s (Villani et al., 2017) are shown. (D) t-SNE map of 4,465 human splenic DCs. Colors indicate unsupervised clustering by Phenograph (left) or classification based on expression of canonical markers (right). Each dot represents an individual cell. (E) Expression of canonical myeloid genes across the transcriptionally defined DC clusters from (D). (F) Heatmap reporting scaled, imputed expression for the top differentially expressed genes for each cluster identified in (D). (G) Pearson correlation between human and mouse spleen scRNA-seq DC clusters. Mouse spleen DC clusters from Figure 1E. (H) Heatmap reporting scaled, imputed expression for the 103 top varying TF genes across cDC2s in each species. Rows are ordered by the TF gene cluster determined using Phenograph and columns are ordered by the annotated cDC2 cluster. (I) Heatmap reporting overlap of top varying TF genes for mouse and human TF gene clusters identified in (H). Each colored count indicates the number of TFs that belong to both the mouse (rows) and human TF gene cluster (columns) indicated. TF gene clusters are further annotated with the cDC2 subset concordant with their expression profile in (H). (J) t-SNE map of human splenic DCs colored by imputed expression of canonical human DC genes (top) or cDC2 lineage-defining TF and marker genes (middle and bottom) identified in Figures 2A and 3D. (K) Diffusion component analysis of human spleen DC clusters identified in (D) illustrating gradients from “AS” DCs to cDC2A and cDC2B clusters. See also Figure S6 and Tables S5, S6, and S7.

Figure S6

Human DC Heterogeneity, Related to Figure 7 (A). Violin plots showing expression distribution of mouse DC subset marker genes across human peripheral blood DC and monocyte clusters identified in Villani et al. (2017). (B). Representative flow cytometric analysis of mouse peripheral blood cDC2s showing absence of T-bet (RFP)⁺ cDC2s. (C). Gating strategy for FACS-isolation of human spleen DCs for scRNA-seq. DCs were defined as live, LIN(CD3,CD56,CD19)⁻CD14^–CD11C⁺HLA-DR⁺. (D). Representative flow cytometry analysis of human spleen cDC2s gated as Lin(CD3,CD56,CD19)^–CD14^–CD11c⁺HLA-DR⁺CD123^–XCR1^–CLEC4A⁺ cells. Left panel: cell surface expression of CD1c and CLEC10A by cDC2s. Right panel: overlay of CLEC10A⁺ and CLEC10A^– cDC2s distinguished by differential expression of CLEC4A and FcεR1a. Summary bar graphs show frequency of CD1C⁺CLEC10A⁺ and CD1C⁺CLEC10A^– cDC2s as a percentage of cDC2s (n = 4 individuals). (E). t-SNE embedding of 9,315 FACS-isolated CD45⁺ immune cells from two melanoma tumors. Colors indicate unsupervised clustering by Phenograph (left panel) or classification based on expression of canonical markers and correlations with bulk RNA-seq data (right panel). Each dot represents an individual cell. (F). Pearson correlations between cluster centroids in (F) and bulk RNA-seq data from purified immune populations (Jeffrey et al., 2006, Novershtern et al., 2011) (G). t-SNE map of 2,122 myeloid cells identified in (F). Colors indicate patient sample (left) or unsupervised clustering by Phenograph (right panel). Each dot represents an individual cell. (H). Heatmap of normalized, log transformed and MAGIC imputed expression of top 20 differentially expressed genes, defined by the highest earth mover’s distance (EMD), per Phenograph cluster in E. The colored bar at the top of the heatmap shows assignment of cells to clusters labeled in F, right panel. (I). t-SNE map of human melanoma myeloid cells (H) colored by imputed expression of labeled genes.

Figure S7

Identification of cDC2 Heterogeneity, Related to Figures 1 and 2 and STAR Methods (A). Clustering robustness measured by adjusted Rand Index (RI) for ranges of principal components and k supplied to Phenograph. In the top panel, with cDC1 (gray) and T-bet+ cDC2 (black) clusters consolidated, RI values indicate robust partitioning between the remaining clusters. In bottom row, when all other cDC2 clusters are consolidated (black), RI values indicate that partitioning is not reliable within cDC1 and T-bet+ cDC2. (B). Heatmap showing the median Rˆ2 value across all genes for MAGIC imputed values calculated with varying timestep (left panel) and number of neighbors (right panel).