Spatiotemporal single-cell analysis elucidates the cellular and molecular dynamics of human cornea aging

Abstract

Background: The human cornea is a transparent and uniquely ordered optical-biological system. Precise coordination of its cellular mechanisms is essential to maintain its transparency and functionality. However, the spatial, cellular and molecular architecture of the human cornea and its intercellular interactions during aging have not been elucidated.

Methods: We performed single-cell RNA sequencing (scRNA-seq) and single-cell SpaTial Enhanced REsolution Omics-sequencing (scStereo-seq) analysis in corneal tissue from eight eyes of donors aged 33-88 years to elucidate the spatiotemporal cellular and molecular dynamics of human cornea aging. Immunofluorescence staining and Western blotting were performed to validate the findings.

Results: Spatiotemporal single-cell analysis revealed the complex cellular landscape, spatial organization and intercellular communication within the human cornea. The subpopulations of major cell types of the cornea were elucidated with precise spatial positions. In particular, we identified novel subpopulations, mapped the spatial positioning of limbal stem cells within the limbal niche, and delineated the interactions between major cell types. We observed that three basal cell subsets migrate centripetally from the peripheral to the central cornea with age, suggesting the "spatiotemporal centripetal pattern" as a novel paradigm for the age-related migration of corneal epithelial cells. Furthermore, we elucidated the age-related, region-specific molecular and functional characteristics of the corneal endothelium, demonstrating differential metabolic capacities and functional properties between the peripheral and central regions.

Conclusions: As the first comprehensive spatiotemporal atlas, our work provides a valuable resource for understanding tissue homeostasis in the human cornea and advances research on corneal pathology, transplantation, senescence and regenerative medicine in the context of corneal aging.

Keywords: Aging; Cornea; Spatial transcriptomics; scRNA-seq; scStereo-seq.

Fig. 1

Spatiotemporal transcriptomic and cell type atlas of the human cornea. a Schematic illustration of the human cornea used for scRNA-seq (10X Genomics), spatial transcriptomics (scStereo-seq) and immunofluorescence. b UMAP visualization of the scRNA-seq data of the 19 primary cell clusters. c Dot plot illustrating the expression of marker genes for each cell cluster. The color indicates the normalized expression level, while the dot size indicates the proportion of expressing cells. d Stacked plot showing the proportion of cell types in each sample. Immune cells include T cells and Langerhans cells. e UMAP visualization of spatial transcriptomics data showing the seven primary clusters. f Spatial visualization of the seven primary clusters in each slide (top), H&E staining of each human cornea (bottom). g Dot plot illustrating the expression of marker genes for each anatomical region. The color indicates the normalized expression level, while the dot size indicates the proportion of expressing cells. h Schematic illustration of cell type mapping from scRNA-seq data to spatial transcriptomics data (left). Spatial visualization of all the cell clusters in each slide (right)

Fig. 2

Spatiotemporal changes in limbal stem cells (LSCs) during corneal aging. a UMAP visualization of two LSC subsets, active limbal stem cells (ALSCs) and quiescent limbal stem cells (QLSCs). b UMAP visualization showing the expression of ATF3 and IFITM3 in LSCs. c Violin plots showing the expression of ATF3 and IFITM3 in two LSC subsets. d GSEA enrichment plots showing the upregulated GO terms of ALSCs (top) and QLSCs (bottom). e UMAP visualization of three ALSC subsets. f Dot plot illustrating the expression of marker genes for each LSC subset. The color indicates the normalized expression level, while the dot size indicates the fraction of expressing cells. g Bar plot showing the results of the GO enrichment analysis of the marker genes of each LSC subset. h Line chart showing the proportion of different LSC subsets in each age stage. i Violin plot showing the stem cell differentiation score for each LSC subset. j Stacked plot showing the proportion of cell phases in each LSC subset. k Trajectory analysis of LSC subsets colored by pseudotime and subsets. l Spatial visualization of each LSC subset in each slide. m and n Immunofluorescence staining of KRT15 (red), CLDN4 (light blue), UBE2S (yellow), and SLC6A6 (green); nuclear DAPI staining is shown in blue

Fig. 3

Spatiotemporal changes in quiescent limbal stem cells (QLSCs) during corneal aging. a Heatmap showing the expression of transcription factors (TFs) for each LSC subset. b Heatmap showing the results of the GO enrichment analysis of HIF1A and its downstream target genes. c UMAP visualization showing the expression of HIF1A in LSCs. d Violin plots showing the expression of IFITM3 and the signature score of wound healing and leukocyte chemotaxis in HIF1A− QLSCs and HIF1A+ QLSCs. e GSEA enrichment plot showing the upregulated GO terms of HIF1A+ QLSCs. f Chord diagram showing the interaction between HIF1A+ QLSCs, HIF1A− QLSCs and vascular endothelial cells. g Lollipop plot showing differences in interaction strength between each cell type and between HIF1A+ QLSCs and HIF1A− QLSCs. h Spatial visualization of HIF1A− QLSCs and HIF1A+ QLSCs, with insets showing enlarged areas of the selected regions. i Immunofluorescence staining of KRT15 (green), HIF1A (red), and CCL14 (lake blue), with nuclei stained with DAPI (blue)

Fig. 4

Diversity and spatial location of corneal basal cells during aging. a UMAP visualization of five corneal epithelial cell (CEC) clusters. b UMAP visualization showing the expression of GJA1, KRT24 and KRT3 in CECs. c Violin plots showing the expression of GJA1, KRT24 and KRT3 in each CEC subset. d Trajectory analysis of CEC subsets colored by pseudotime and subsets. e Dot plot illustrating the expression of marker genes for three basal subsets (left). The color indicates the normalized expression level, while the dot size indicates the fraction of expressing cells. Bar plot displaying the results of the GO enrichment analysis of the marker genes (right). f Spatial visualization of each basal subset in each slide, with insets showing enlarged areas of the selected regions. g Immunofluorescence staining of NTRK2 (purple) and DAPI (blue) in corneal sections. h Schematic representation of the corneal regions. i Western blot showing NTRK2 protein levels in young (8 weeks old) and aged (14 months old) mouse corneas. j Histogram showing the comparative expression levels of NTRK2 in the corneas of young and aged mice. The error bars indicate the standard error of the mean (p* < 0.05). k Violin plot showing the epithelial cell differentiation score for each basal subset. Chord diagram showing the interaction between three basal subsets, wing cells, and squamous cells. l GSEA enrichment plot showing the enriched GO terms associated with genes whose expression was upregulated in NTRK2+THBS1+ basal cells compared to that in NTRK2+ THBS1− basal cells. m Spatial visualization of NTRK2+THBS1 + basal cells, NTRK2+THBS1− basal cells, and other CECs in each slide, with insets showing enlarged areas of the selected regions

Fig. 5

Characterization of stromal cells during aging. a UMAP visualization of five stromal cell subsets. b UMAP visualization showing the expression of DCN and S100B in stromal cells. c Dot plot illustrating the expression of marker genes for four stromal cell subsets and corneal stromal stem cells (CSSCs). The color indicates the normalized expression level, while the dot size indicates the proportion of expressing cells. d Heatmap showing the expression of differentially expressed genes (DEGs) for three stromal cell subsets (left). Bar plot showing the results of GO enrichment analysis of the DEGs of three stromal cell subsets (right). e Spatial visualization of each stromal cell subset in each slide. f Heatmap showing the expression of collagen, extracellular matrix (ECM) proteoglycans and elastic fiber-encoded genes for the four stromal cell subsets in the spatial transcriptomics data. g Heatmap showing the expression of collagen, ECM proteoglycans and elastic fiber-encoded genes for stromal-1 (left) and stromal-2 (right) across all age groups in the spatial transcriptomics data

Fig. 6

Spatial, molecular and functional heterogeneity of corneal endothelial cells. a Spatial visualization delineating the peripheral and central corneal endothelium. b Spatial visualization of endothelial cells with different percentages of mitochondrial gene expression. c Violin plots showing the percentages of mitochondrial gene expression in the peripheral and central corneal endothelium. d Scatter plot showing log2-fold changes in gene expression in the peripheral and corneal endothelium. e Bar plot showing the results of the GO enrichment analysis of the DEGs associated with the peripheral and central corneal endothelium. f GSEA enrichment plot displaying the significantly upregulated metabolic pathways of the peripheral corneal endothelium. g Lollipop chart showing the functional scores of the peripheral and central corneal endothelium across different age groups. h Corneal endothelial cells costained with CA3 (red) and the Na–K-ATPase subunit (green), with DAPI in blue (left). A statistical image of the Na–K-ATPase immunofluorescence results is shown on the right (p*** < 0.001). i Immunofluorescence staining of mouse corneal whole mounts. NDUFV2 (I) (red) and SDHB (II) (lake blue) are expressed at higher levels in the peripheral corneal endothelium than in the central region. An immunofluorescence statistical diagram of NDUFV2 (I) and SDHB (II) is shown on the right (p*** < 0.001)