Ectopic cervical thymi and no thymic involution until midlife in naked mole-rats

Immunosenescence is a hallmark of aging and manifests as increased susceptibility to infection, autoimmunity, and cancer in the elderly. One component of immunosenescence is thymic involution, age-associated shrinkage of the thymus, observed in all vertebrates studied to date. The naked mole-rat (Heterocephalus glaber) has become an attractive animal model in aging research due to its extreme longevity and resistance to disease. Here we show that naked mole rats display no thymic involution up to 11 years of age. Furthermore, we found large ectopic cervical thymi in addition to the canonical thoracic thymus, both being identical in their cell composition. The developmental landscape in naked mole-rat thymi revealed overt differences from the murine T cell compartment, most notably a decrease of CD4+/CD8+ double-positive cells and lower abundance of cytotoxic effector T cells. Our observations suggest that naked mole rats display a delayed immunosenescence. Therapeutic interventions aimed at reversing thymic aging remain limited, underscoring the importance of understanding the cellular and molecular mechanisms behind a sustained immune function in the naked mole rat.


Introduction
Age-related loss of cellularity and weight of the thymus, termed thymic involution, occurs in all mammals (Shanley et al., 2009). Thymic involution contributes to immunosenescence, leading to increased incidence of cancers, autoimmunity, and opportunistic infections in the elderly (Dixit, 2012, Lynch et al., 2009. Thymic involution remains an evolutionary mystery since it occurs in most vertebrates despite its negative effects. In vertebrates aging of the thymus gland precedes age-associated phenotypic changes in other organs, with thymic function gradually decreasing from the first year of human life (Shanley et al., 2009). The process of human thymic involution is incompletely understood, in part due to a lack of animal models showing delayed immunosenescence (Tong et al., 2020).
Naked mole-rats (NMR) display exceptional longevity, with a maximum lifespan of 32 years (Buffenstein, 2008, without increased mortality due to aging (Ruby et al., 2018). In comparison, a similarly sized house mouse has a maximum lifespan of 4 years (de Magalhaes et al., 2005, Turturro et al., 1999. NMRs display negligible senescence (Buffenstein, 2008), characterized by very slow changes in physiological parameters with age, and the lack of an age-related increase in mortality rate (Finch, 1990). Accordingly, NMRs do not exhibit age related changes in basal metabolism, body composition, or bone mineral density (O'Connor et al., 2002, Buffenstein & Ruby, 2021. NMRs show a very low incidence of cancer (Buffenstein, 2008, Delaney et al., 2013, Delaney et al., 2016. NMRs produce abundant high molecular weight hyaluronic acid (HMW-HA) responsible for their resistance to solid tumors (Tian et al., 2013a, Tian et al., 2015 and feature a serum metabolome resembling calorically restricted mice (Lewis et al., 2018, Puppione et al., 2021. Thus, NMRs present a unique model to study the mechanisms of healthy longevity (Edrey et al., 2011. However, the information about the immune system of NMRs is limited. Very recently, whole spleen single cell RNA-Sequencing revealed the absence of canonical Natural Killer Cells (NKC) in NMRs (Hilton et al., 2019). It was also reported that blind mole-rats, an unrelated clade of subterranean rodents that convergently evolved extreme longevity, sustain sustained T cell repertoire diversity in old age (Spalax spp.) (Izraelson et al., 2018). At present there is, however, no data regarding the T cell compartment in NMRs.
Here we set out to characterize the thymus and T-Lymphopoiesis in long-lived NMRs in comparison with short-lived mice. Our results revealed unexpected presence of additional cervical thymi in the NMR, and the absence of thymic involution until 11 years of age, which is the oldest age of animals in our research colony.

Results
Thymic involution and decline of T cell function predisposes to opportunistic diseases and presents a major risk factor in the elderly (Dixit, 2012, Torroba & Zapata, 2003. Since NMRs are characterized by negligible senescence, we set out to examine whether they display thymic involution. To enable the analysis of the NMR hematopoietic cells, we previously developed a flow cytometry (FACS) staining with cross-reactive monoclonal antibodies (moAbs) and functionally characterized purified hematopoietic stem and progenitor (HSPC) fractions (Emmrich et al., 2019). Through CITE-Seq we showed that naked mole-rat T cells (TC) express CD3, LAT and LCK are immunophenotypically Thy1.1 int /CD34 -/CD11b - (Figure S1a, b).
Thymus glands extracted from the thoracic cavity of young NMRs and mice both featured similar sized lymphocytes and largely absent myeloid cell fractions ( Figure S2a). When we extracted lymph nodes (LN) from 43 naked mole-rat necks we noted two types of cervical nodes (Figure 1a): the first type corresponding to LN, and the second resembling a thymus, which we termed cervical thymi. Thoracic and cervical thymi contained a CD34 + fraction while LNs did not, and shared histological features such as many medullary sinuses, no germinal centers and more basophilic cytoplasms than LNs (Figure 1b-c). Conceivably, both thymi types stain positive for cytokeratin as seen for mouse and human thymi, while LNs are cytokeratinwith a much lower TC:BC ratio (Figure 1d, Figure S2b). Embryonic thymus development proceeds through detachment of the primordium, which contains thymus-and parathyroid-anlagen, from the 3rd pharyngeal pouches, followed by separation of the anlagen and migration of thymic lobes into the chest cavity (Gordon & Manley, 2011). We compared the macroscopic anatomy with mice and found undersized thoracic thymus lobes in NMR neonates (Supplementary S2c). To our surprise, the cervical thymus was clearly visible in the ventral pharyngeal region, and appeared markedly larger than their thoracic thymus ( Figure S2d).

Conserved thymopoiesis between mouse and naked mole-rat
Next, we performed whole thymus CITE-Seq from two 3 month-old and two 12 month-old mice versus two 3 year-old and two 11 year-old NMRs ( Figure S3a). The classical thymusspecific FACS pattern of CD4 + /CD8 + double-positive (DP) thymocytes in mice ( Figure S3b) was resolved into 13 cell states of 20019 single cell transcriptomes (Figure 2a). We captured sorted naked mole-rat CD34 + thoracic thymus cells ( Figure S3c) and two LNs from the younger animals and integrated those with naked mole-rat thoracic and cervical thymi. The integrated dataset comprised 70928 cells across 11 droplet libraries from 3 tissues and one FACS-purified sample (Figure 2b). We found 15 distinct communities which were annotated based on canonical marker gene expression for respective TC subsets (Figure 2c). The earliest developmental stage in the thymus are CD4 -/CD8double-negative (DN) early T cell progenitors (ETP), while the expression of CD44 and CD25 further subdivides the DN state. A DN2/3 cluster was the earliest distinct developmental time point detected, in mice showing exclusive CD25 hi expression, a Kit high-tolow-expression gradient and a fraction of CD44 lo cells as measured by CITE-signals ( Figure S3d).
In both species expression of PTCRA, NOTCH1 and its target HES1 was conserved in DN2/3 cells (Figure S3e- initiates TCR-α gene rearrangements and upregulates expression of CD4 and CD8 to yield DPs, which usually progress through an immature cycling CD8 + intermediate single-positive population (ISP) (MacDonald et al., 1988). We found an equivalent cell state DN4/ISP in both species, which was marked by overexpression of cell cycle genes, retention of PTCRA and rising expression of CD4/CD8. Remarkably, sorted CD34 + naked mole-rat thymocytes were strongly enriched for DN2/3 and DN4/ISP clusters (21%, 35%; Figure S3g), while CD34 mRNA and CITE-signal were specific to DN2/3. This trait is shared with humans, who maintain CD34 + primitive ETPs (Terstappen et al., 1992).
Whole thymus is usually comprised of 80-90% DPs ( Figure S3b). Louvain clustering identified several DP communities in each species, which we labelled early (eDP) and late (lDP) according to their CD4/CD8 transcript levels. In NMRs an additional DP cluster was formed, which was detected to 6% in LNs and thus named DP.ext (extrathymic). Total DP cell frequencies were 85% in mice and 64% in naked mole-rats (Figure 2d). Expression of Rag1, anti-apoptotic Dek and Themis, indispensable for proper positive and negative TC selection in mice (Johnson et al., 2009), was conserved in eDPs, whereas lDPs activated Tox, Helios and Gata3 transcription ( Figure S3d-e). Strikingly, cell type distributions between the two NMR thymus types were identical ( Figure S3h), supporting evidence was obtained from quantitation of Thy1.1 int /CD34 + and TC FACS populations, showing no difference between thymus origins (Figure 1c).
Conversely, BCs are prevalent in LNs (Figure S3h), and CD125 + BC frequencies are higher in LNs than in BM, spleen or thymi ( Figure S4a-b), suggesting canonical B-lymphopoietic functions in peripheral LNs.
In summary, we found an additional pair of functional thymi in naked mole-rats. Early and intermittent steps of T-lineage development appear to be conserved, however, there is a stark decrease of naked mole-rat DP proportions compared to mice.

Cryptic γδ T-lymphocytes with a killer cell signature in naked mole-rats
Commitment towards the γδTC lineage is instructed by T cell receptor (TCR) signal strength, whereby PTCRA indicates weak and TCRγδ strong signals at the DN1 stage (Munoz-Ruiz et al., 2017). We found the constant region of the TCRγ chain (TRGC2) marking marrow ETPs (Emmrich et al., 2019), hence the most primitive thymic partition across mouse and naked mole-rat showed conserved TRGC2 overexpression (Figure 3a). Mature murine γδTCs are enriched for TRGC2. In naked mole-rats, however, several mature subsets contain fractions of TRGC2-expressing cells. JAML has been shown to induce γδTC activation (Witherden et al., 2010), conversely the mouse thymic γδTC cluster showed highest JAML expression ( Figure S3f).
However, we did not detect a separate γδTC population in the naked mole-rat dataset, although it had 3.5-fold more cells due to the additional cervical thymi and LNs, thereby increasing clustering resolution. Intriguingly, JAML is one of the top DN cell markers in naked mole-rats and overexpressed in CD8-TCs ( Figure S3f). A CD4 + cluster comprising cytotoxic T-lymphocytes (CTL), specific for GZMA and NKG7, had a cell fraction positive for TRGC2 and JAML.
Interestingly, TRGC2/NKG7 co-expressing cells were also found in CD8-TCs (Figure S4c-d). We therefore compared the specific cluster markers of mouse γδTC with naked mole-rat CD4-CTL (Figure 3b), revealing a strong correlation in a shared subset of genes encompassing S100A4/10/11, RORA, IL7R and CD44 (Supplementary Table 2). Conclusively, putative naked mole-rat γδTCs appeared overtly as CTLs, similar to a human Vγ9Vδ2 + /CD45RA + /CD27 − effector memory γδTC in LNs (Caccamo et al., 2005). By contrast an abundance of murine γδTC subsets with immunomodulatory functions through mainly secreting either IFN-γ of Il-17a, but no CTLs, has been described (Pang et al., 2012). We added further evidence to this by integration of thoracic thymi scRNA-Seq from mouse and naked mole-rats, encompassing 48061 cells across 10871 common genes (Figure 3c). Differential abundance quantitation confirmed the betweenspecies bias in thymic cell composition, with all mature TC clusters significantly more frequent in naked mole-rats (Figure 3d). Strikingly, the integrated NK/CTL cluster mapped to mouse NKCs and γδTCs, which did not co-cluster in the mouse-only dataset, and cross-mapped to naked molerat CD4-CTLs containing the putative killer cell γδTCs (Figure 3e-f). Moreover, by mapping the constant TCR chain region transcripts from all TCR loci using the most recent NMR genome (Zhou et al., 2020), we performed absolute copy number qPCR in sorted PB-TCs and saw significantly less TRGC1/2 expression in naked mole-rats ( Figure S4e). Therefore, NMRs have a cryptic γδTC population with a CTL expression signature in the thymus, resembling human effector memory γδTCs with killer cell function. Unlike their αβTC counterparts that require peripheral activation for effector cell differentiation, γδTCs can be 'developmentally programmed' in the thymus to generate discrete effector subsets with distinctive molecular signatures (Munoz-Ruiz et al., 2017). Our data indicates that NMR γδTCs are overtly programed towards CTLs, potentially compensating the lack of NKCs and a diminished CD8-TC subset.

No signs of thymic immunosenescence in middle-aged naked mole-rats
In mammals, thymus tissue weight decreases with age, largely by decomposition of the thymic niche required for proper TC development, wherein thymic epithelial cells (TECs) and other stromal components are replaced by adipocytes (Palmer, 2013). Mice showed a continuous decline in thymus cellularity, which is proportional to tissue weight loss, from 3 months up to 2 years of age (19-fold, 3m-median 73•10 6 vs 24m 3.9•10 6 cells; Figure 4a), as reported earlier (Sempowski et al., 2002). To our surprise we saw a significant increase in NMR thymus cellularity between 3 and 11 years of age, for both thoracic and cervical thymi. Interestingly, cervical thymi were found to contain 10-fold more thymocytes than thoracic counterparts, which was also evident from visibly smaller thoracic lobe sizes as seen above in neonates ( Figure S4f). However, while murine thymus volumes strongly diminished between 3 and 12 month of age, there was no shrinkage in NMR thymi in older specimen despite the 8 year age difference. We further observed cervical LNs from either species to retain size and cellularity over aging, but while LN cell numbers were consistent even between species (Figure S4g), NMR LNs had a much larger volume than those of mice.
We next probed several correlates of thymic aging and its concomitant loss of T-lineage potential. High autoimmune regulator gene (AIRE) activity of cortical and medullary TECs orchestrates autoantigen presentation during negative selection, and its thymic expression pattern is proportional to thymic involution (Perniola, 2018). As expected, mouse AIRE levels decreased with age as measured by absolute qPCR from whole thymi, whereas NMR AIRE was consistently expressed throughout age in naked mole-rats (Figure 4b). FOXN1 is the major ontogenic thymus marker by regulating TEC development and function from prenatal stage until after birth. In rodents referred to as the 'nude locus', genetic disruption of Foxn1 causes athymia and hairlessness in mice and rats, and a number of TC immunodeficiency syndromes have been linked to the human ortholog (Romano et al., 2013). In mice, Foxn1 levels sharply decline in the postnatal thymus, on the contrary NMR FOXN1 remains at the same expression level as detected in neonates, a clear neotenic feature (Emmrich et al., 2019) (Figure 4c). It is important to note that both AIRE and FOXN1 expression in NMR neonates was lower by 1-2 orders of magnitude than in mouse neonates, but whether this difference pertains to less TECs or less AIRE/FOXN1 in similar number of TECs between species warrants further investigation. We detected a TEC population in NMR thymi by scRNA-Seq with pronounced AIRE expression, whereas the mouse dataset did not feature a distinct TEC cluster ( Figure S3d-e). Another hallmark of immunosenescence is the steady decline in frequency and functionality of DN early thymic progenitors in mice (Min et al., 2004). We showed that the most primitive NMR thymic progenitor compartment is CD34 + , and thus quantified this fraction in thymi across an 11yr timespan.
Remarkably, neither thoracic nor cervical thymi CD34 + populations diminished during this period (Figure 4d), albeit both organs showed a slight trend towards reduction, which likely becomes significant in animals aged >20yrs. The pronounced effect of thymic aging on the thymocyte pool is the dramatic decrease of naïve TCs with an accompanied accumulation of effector and memory subsets in the aged thymus. Although we did not further subset clusters due to lack of established naked mole-rat markers for memory/effector vs naïve, the Simpson index as a measure of species diversity in ecosystems can be used to score fluctuations in cell type frequencies over age (Ferrall-Fairbanks & Altrock, 2021), where represents infinite diversity and 0, no diversity. In the naked mole-rat, the Simpson index remained close to 0.8 between age groups spanning 8 years, regardless thymus type (Figure 4e), while in mice it rose from 0.61 to 0.73 between age groups covering 9 months. Per cell type differential abundance quantitation showed a 5-fold increase of BCs (p=0.006) and 45-fold increase of TCs (p=0.011) in older mice, with 5 further cell types changing >1.5-fold ( Figure S4h). In contrast, no significant up-or down-regulated cell types were found and only 3 cell types were changed by >1.5-fold in naked mole-rats ( Figure S4i). We further found that CD4 and CD8A mRNA levels were as consistent across age as FOXN1 and AIRE in naked mole-rat thymi (Figure S4j-k). These results show that in the naked mole-rat thymocyte pool composition and transcriptional patterns are maintained for over a decade, whereas in mice changes on the cellular level coincide with onset of functional thymic regress as early as within 1 year of life.

Discussion
Here we provide the first characterization of the naked mole rat thymus. We discovered that naked mole rats have an additional pair of cervical thymi. This is an unexpected finding as mammals, including humans and mice, as a rule, have only one bilateral thymus. Cervical thymi can occasionally be detected in mice, but their frequency is rare and they have unilateral appearance (Dooley et al., 2006). Similarly, rare ectopic cervical human thymi had been reported in children (Ahsan et al., 2010). In contrast, cervical thymi are a principal component of NMR ontogenesis. Interestingly, among vertebrates, chickens have seven, sharks five, and amphibians three thymi (Boehm & Bleul, 2007). It is tempting to speculate that the presence of additional thymi in the naked mole rat may contribute to prolonged maintenance of immune function during their lifespan.
The ectopic thymus reflects a failed migration of thymic tissue from the third pharyngeal pouch endoderm during organogenesis, which can be found at any level of the pathway of normal thymic descent, from the angle of the mandible to the superior mediastinum (Saggese et al., 2002). It is possible that the naked mole-rat thymic anlage splits during migration and one remains in the throat. Alternatively, their parathyroid glands may have been repurposed to cervical thymi, which is less likely due to presence of a conserved PTH ortholog in naked mole-rat genome assemblies. Another explanation could be alterations as seen in ephrinB2 mutants (Gordon & Manley, 2011), wherein the thymus remains in the anterior pharyngeal region.
We provide evidence for a delay of thymic involution in naked mole-rats beyond the 1 st decade of their lifespan. Age-associated marker expression and thymic cell composition remained at the level of neonates. The absence of thymic involution up to midlife is unprecedented in mammals. This would translate into similar or even slightly heightened thymic weights and cell counts for humans in their 30's. Thymic involution decreases output of naïve T cells and reduces the ability to mount protective responses against new antigens. In naked mole-rats we did not see thymic involution in animals >10 years old, while markers for thymic function and development, AIRE and FOXN1, were maintained at neonatal levels. Furthermore, the reduction of ETPs accompanying age-related lymphoid decline did not manifest in naked mole-rats, arguing that their intrinsic myeloid bias in the marrow does not predispose HSPCs towards less lymphoid commitment (Emmrich et al., 2019). However, naked mole-rats are not immortal and do show frailty in old age (Edrey et al., 2011). Therefore, an eventual decline in thymic cellularity and immune function is to be anticipated, albeit delayed as opposed to the lifelong steady decline in humans and mice.
Neoteny refers to retention of juvenile phenotypes in adult organisms, hence considering humans neotenic apes (Bufill et al., 2011). NMRs feature an array of neotenic traits (Skulachev et al., 2017), including aspects of their hematopoietic system (Emmrich et al., 2019). Here we found developmental FOXN1 and age-associated AIRE mRNA levels with little to no changes between neonate and 11 year-old adult animals. Similarly, thymic ETPs remain at neonate frequencies in NMRs. Thymic involution occurs in almost all vertebrates (Shanley et al., 2009), hence neotenic retention of a juvenile thymus in mature, aged animals represents a likely function of longevity by maintenance of youthful TC-mediated immune function during adulthood.

Experimental Procedures Animals
Ethical and legal approval was obtained prior to the start of the study by the University of Rochester Committee on Animal Resources (UCAR). All animal experiments were approved and performed in accordance with guidelines instructed by UCAR with protocol numbers 2009-054 (naked mole rat) and 2017-033 (mouse). Naked mole rats were from the University of Rochester colonies, housing conditions as described (Ke et al., 2014). C57BL/6 mice were obtained from NIA.

Primary cell isolation
Marrow from mice and naked mole-rats was extracted from femora, tibiae, humeri, iliaci and vertebrae by crushing. Thymus and lymph nodes were minced over a 70µm strainer and resuspended in FACS buffer. Blood from mice was drawn via retroorbital capillary bleeding, naked mole-rat blood was obtained via heart puncture.

Histology
Imaging and analysis was performed using a using a Nikon Eclipse Ti-S microscope.
Coverslips were applied with DEPEX Mounting media (Electron Microscopy Sciences), except for Alkaline Phosphatase staining where Vectashield Hard Set Mounting Medium for Fluorescence (Vector) was applied. Soft tissues were stored in 10% neutral buffered formalin, processing was done using a Sakura Tissue-Tek VIP 6 automated histoprocessor, paraffin embedding was done using a Sakura Tissue-Tek TEC 5 paraffin embedding center. A Microm HM315 microtome was used to section tissues at a thickness of 5µm, which then were floated onto a slide with a water bath at a temperature between 45°C and 55°C. Sections were deparaffinized and rehydrated to distilled water through xylene and graded ethanol (100% to 70%).

Hematoxylin & Eosin: Sections were stained with Mayers Hematoxylin (Sigma) for 1min
and washed with tap water to remove excess blue coloring. Soft tissue sections were further decolorized with 3 dips in 0.5% acid alcohol and washed in distilled water. The nuclei of sections were blued in 1X PBS for 1 minute and washed again in distilled water. An Alcoholic-Eosin counterstain was applied for 30sec before slides were immediately dehydrated and cleared through 3 changes of 95% ethanol, 2 changes of 100% ethanol, and three changes of Xylene for 1min each.
Cytokeratin: Paraffin sections (4µm thick) of FFPE thymus and lymph node tissues (NMR, mouse, and human control) were stained for cytokeratin (AE1/AE3, Dako GA05361-2) on a Dako Omnis autostainer with pressure cooker antigen retrieval (TrisEDTA; pH 9). A section of normal human thymus resected for routine clinical care at the university of Rochester medical center was stained for H&E and cytokeratin for morphological comparison. The human, mouse, and NMR samples were stained in the same run on the same machine.

Flow Cytometry
Flow cytometry analysis was performed at the URMC Flow Core on a LSR II or LSRFortessa (both BD), or on our labs CytoFlex S (Beckman Coulter). Kaluza 2.1 (Beckman Coulter) was used for data analysis. Staining and measurement were done using standard protocols. Red blood cell lysis was done by resuspending marrow pellets in 4ml, spleen pellets in 1ml and up to 500µl blood in 20ml of RBC lysis buffer, prepared by dissolving 4.1g NH4Cl and 0.5g KHCO3into 500ml double-distilled H2O and adding 200µl 0.5M EDTA. Marrow and spleen were incubated for 2min on ice, blood was lysed for 30min at room temperature. Cells were Sorting was performed at the URMC Flow Core on a FACSAria (BD) using a 85μm nozzle, staining was done as described. Human HSCs were sorted for population RNA-Seq as LIN -/CD34 + /CD38 Lo /CD45RA -/CD90 Dim ( Figure S5A). Naked mole-rat HSPC populations were sorted as described with a lineage cocktail comprised of CD11b, CD18, CD90 and CD125 (NMR LIN).

Quantitative PCR
Mouse and Naked mole-rat sorted TCs and thymic tissue were used for RNA extraction by Trizol (Thermo Fisher). RNA was quantified using a NanoDrop One (Thermo Fisher), and 100ng was used as input for the High Capacity cDNA Reverse Transcription Kit (Thermo Fisher).
RT reaction was performed according to instructions and the 20µl reaction diluted to 200µl, of which 5µl were used per qPCR reaction. We used iTaq Universal SYBR Green Supermix (Bio-Rad) on a CFX Connect® RealTime System (Bio-Rad) with a three-step cycling of 10sec 95°C, 20sec 60°C, 30sec 72°C for 40 cycles. All primers (IDTDNA) were validated to amplify a single amplicon at the above PCR conditions by gel electrophoresis. Gene sequences for primer design by Primer3Plus were retrieved from ENSEMBL, with the exception of the T cell receptor C-region genes for naked mole-rat. Here we used the WBM RNA-Seq from the transcriptome assembly below to map those genes in a recently published naked mole-rat genome (Zhou et al., 2020) using Apollo software and custom scripts. For absolute copy number quantitation, qPCR amplicons were gel-purified using the QIAquick Gel Extraction Kit (Qiagen) and subcloned into the pCR2.1 plasmid using the TOPO-TA cloning Kit (Thermo Fisher). Plasmids were prepared using the QIAprep Spin Miniprep Kit (Qiagen). Sanger sequencing was performed by Genewiz using M13 forward and reverse primers. Standard curves were prepared across a 10-fold dilution range from 20ag to 20pg of plasmid DNA. All amplicon gel images, amplicon plasmids and standard curve data is available upon request.

from either species with
FindIntegrationAnchors with dims = 1:50, anchor.features = 3000, reduction = "cca". For naked mole-rat we integrated 4 thoracic thymus, 4 cervical thymus, 2 lymph node and 1 CD34 + thymocyte libraries (Figure 2b). Scores for G2M and S phases were obtained using Seurat CellCycleScoring as described in the respective Seurat vignette. Clustering was done using Seurat's FindClusters function with resolution = 0.5. Next we used the doublet detection and removal workflow as suggested in the Bioconductor OSCA vignette. Briefly, we run findDoubletCluster from the scDblFinder package, followed by in silico simulation of doublets from the single-cell expression profiles (Dahlin et al., 2018) using computeDoubletDensity from BiocSingular package, and excluded any cluster which was identified in both methods. The DEGs for each cluster were detected by FindAllMarkers function with arguments test.use = "MAST", logfc.threshold = log(2), min.pct = 0.25, return.thresh = 0.05. CITE feature/antibody marker detection was done as described for transcript cluster markers with the exception of test.use = "wilcox". The heatmap of the top overexpressed markers for each cluster in Figure S7A shows that CCA integration worked efficiently. However, cells from sorted CD34 + thymocytes of the integrated clusters PC, and to a lesser extent eCD4-SP, feature a marker signature of combined DN2/3 and DN4/ISP clusters. We attribute this primarily to the difference in Chromium Single Cell 3' Reagent Kits (10X Genomics), which were v3 chemistry for all integrated lymphoid tissues, except for the sorted CD34 + thymocyte library captured with v2 chemistry. We performed manual curation of the cell type annotation based on canonical markers from the literature. Mapping of naked mole-rat γδTCs to the CD4-CTL cluster was done by intersecting 620 mouse γδTC markers with 129 naked mole-rat CD4-CTL markers, yielding 45 genes in both clusters, which we regressed by a linear model from the stats package (Figure 3b). DA testing was performed as described (Lun et al., 2017). For the CCA-integrated naked mole-rat lymphoid dataset we show cell type abundances between age groups across both thoracic and cervical thymi ( Figure S4i). Regardless of DA testing across both or separate testing of either thoracic or cervical thymi, no cell type was significantly (FDR < 0.05) changed. SCTransform was used to integrate scaled, clustered and annotated mouse and naked mole-rat unfractionated thoracic thymus datasets: SelectIntegrationFeatures with nfeatures = 3000, FindIntegrationAnchors with normalization.method = "SCT". Cell cycle scoring, clustering and marker detection was performed as described above. Simpson's index as calculated by D S = 1 − ∑ n i (n −1) N(N−1) using the vegan package was determined as diversity of cell types across libraries (Figure 4e).

Quantification and Statistical Analysis
Data are presented as the mean ± SD. Statistical tests performed can be found in the

Conflict of Interest Statement
The Authors declare that they have no conflict of interest.

Data Availability Statement
The data that supports the findings of this study are available in the supplementary material and Supplementary Tables of this article. The single-cell RNA-Sequencing Raw Data is available at figshare, see doi link and reference number ###. Further information and requests for resources and reagents should be directed to and will be replied by the Lead Contact, Stephan Emmrich (Stephan-Emmrich@gmx.net). Zhou X, Dou Q, Fan G, Zhang Q, Sanderford M, Kaya A, Johnson J, Karlsson EK, Tian X, Mikhalchenko A, Kumar S, Seluanov A, Zhang ZD, Gorbunova V, Liu X, Gladyshev VN (2020) Beaver and Naked Mole Rat Genomes Reveal Common Paths to Longevity. Cell Rep. 32 (4), 107949.

Figures and Legends
shows the top 25 cell type specific mRNA markers for mouse thymus randomly downsampled to ≤ 500 cells, canonical cell type markers from the literature are indicated. Bottom heatmap shows the cell type specific CITE features for mouse thymus randomly downsampled to ≤ 100 cells.
Fold-change cut-off 2, p-value threshold 0.05. The CD137-TC cluster is partially composed of CITE-CD4cells, while TNFRSF9 (CD137) mRNA was 5.3-fold upregulated (expressed in 45.7% CD137-TC vs 1.1% in all other clusters combined; adjusted p<10 -304 ). e, Heatmap of the top 25 cell type specific mRNA markers for the naked mole-rat lymphoid dataset from Figure 6E randomly downsampled to ≤ 500 cells, canonical cell type markers from the literature are indicated. Bottom traces depict differentially expressed CITE features as log-UMI counts. Top color bar labels tissue origin (see legend bottom right), 2 nd color bar labels cell type annotation.    Figure S4. ScRNA-Seq continued and qPCR from mouse and naked mole-rat thymi. a, FACS gating of CD125 + BCs in naked mole-rat BM [left], spleen [middle] and LN [right]. Note that the lymph node CD125 + population has two separate CD125 dim and CD125 bright fractions. This was seen in 1/5 animals with no apparent disease condition noticed elsewise, other animal LNs showed CD125 + staining pattern comparable to BM. b, Frequencies of CD125 + in BM (n=8), Spleen (n=7) and lymph nodes (LN, n=10); p-value determined by Dunnett's One-way ANOVA.