scRNA-seq reveals mTEC heterogeneity in rat thymus. (A) Uniform Manifold Approximation and Projection (UMAP) embedding of single-cell transcriptomic data from Aire-deficient and control TECs. (B) Dot plot for expression of marker genes. Color represents scaled normalized mean expression of marker genes in each TEC population, and size indicates the fraction of cells expressing marker genes. (C) Module scores for TRA and Aire-dependent (Aire-dep) gene sets computed using Seurat’s AddModuleScore confirm the identity of mTEC II subpopulations. Aire-dependent genes were defined by analyzing bulk transcriptome data from purified Aire-deficient and control mTEC^hi cells. Tissue-restricted genes were inferred as outlined in Materials and methods. (D) Upper panel: number of genes expressed by mTEC subpopulations as a function of UMIs considered. Rarefaction curves were interpolated and extrapolated using the iNEXT R package. Lower panel: transcriptomic diversity of TEC populations, quantified as a Hill number (qD) of order 0 (****P < 1e-4, P values are calculated using Wald-type Z-test based on confidence interval-derived standard errors with mTEC II.2 as a reference group). (E) Hierarchical clustering of TEC population based on gene expression correlation distance supports the functional proximity of mTEC II subpopulations. Red asterisks at each node indicate distinct clusters (P < 0.05) based on bootstrap resampling using pvclust R package. (F and G) Pseudotime trajectory curve describing differentiation of the mTEC lineage overlaid on UMAP showing mTEC subpopulations (F) and pseudotemporal ordering (G). (H) Distribution density plots of mTEC subpopulations ordered along pseudotime showing the successive stages of differentiation. (I) Log-normalized expression of selected genes in Aire-deficient and control rats showing diminished expression of Aire-regulated genes in post-Aire mTEC populations. Throughout the figure, n = 4 for KO and n = 3 for HE animals, and the data are obtained from a single experiment.