Prediction of cell states and key transcription factors of the human cornea through integrated single-cell omics analyses

Abstract

The cornea, a transparent tissue composed of multiple layers, allows light to enter the eye. Several single-cell RNA-seq (scRNA-seq) analyses have been performed to explore the cell states and to understand the cellular composition of the human cornea. However, inconsistences in cell state annotations between these studies complicate the application of these findings in corneal studies. To address this, we integrated scRNA-seq data from four published studies and created a human corneal cell state meta-atlas. This meta-atlas was subsequently evaluated in two applications. First, we developed a machine learning pipeline cPredictor, using the human corneal cell state meta-atlas as input, to annotate corneal cell states. We demonstrated the accuracy of cPredictor and its ability to identify novel marker genes and rare cell states in the human cornea. Furthermore, cPredictor revealed the differences of the cell states between pluripotent stem cell-derived corneal organoids and the human cornea. Second, we integrated the scRNA-seq-based cell state meta-atlas with chromatin accessibility data, conducting motif-focused and gene regulatory network analyses. These approaches identified distinct transcription factors (TFs) driving cell states of the human cornea. The novel marker genes and TFs were validated by immunohistochemistry. Overall, this study offers a reliable and accessible reference for profiling corneal cell states, which facilitates future research in cornea development, disease, and regeneration.

Keywords: corneal biology; gene regulatory networks; machine learning; scATAC-seq; scRNA-seq.

Fig. 1.

Integration, identification, and characterization of cell states in the human cornea. A) Uniform manifold approximation and projection (UMAP) of integrated and annotated cell states. Colors indicate overall cell type, and subcolors define the respective cell states. The bar plots on the right depict dataset contributions for each cell state. B) Violin plot of the literature annotated marker gene expression across integrated cell states. The height of the violin plot indicates the expression level of the markers, and the color indicates the median score for all cells in that specific cell state. C) GO-term enrichment of selected cell states on all cluster markers identified with Seurat. Names of cell states (A–C): LSC-1, limbal stem cells 1; LSC-2, limbal stem cells 2; LSE, limbal suprabasal epithelium; CE, central epithelium; Cj, conjunctiva; qSK, quiescent stromal keratocytes; SK, stromal keratocytes; TSK, transitioning stromal keratocytes; CF, corneal fibroblasts; CEC, corneal endothelial cells; B/L EC, blood and lymph endothelial cells; Mel, melanocytes; IC, immune cells; nm-cSC, nonmyelinating corneal Schwann cell; MC, mural cell.

Fig. 2.

Immunofluorescence of novel markers in the human cornea. A), B) Costaining of CPVL and SLC6A6 (green) with p63 (ΔNp63) (red) in limbus and central cornea. Staining of corneal stroma marker keratocan (C), POU3F3 protein (D), and fibulin-1 (E). F) Staining of TNNC1 in corneal endothelium. Arrowheads depict cells in the stroma with nuclear expression of POU3F3. Dotted lines indicate borders between stroma and epithelium (in white) and between stroma and endothelium (in yellow). Cell nuclei were stained with DAPI. Scale bar represents 50 μm.

Fig. 3.

A machine learning–based prediction pipeline, cPredictor, for human corneal cell states in scRNA-seq datasets. Prediction certainty plots of corneal cell states on adult corneal cells (A), 4-month-old organoids from Maiti (54) (B), and 4-month-old organoids from Swarup (55) (C). Left, the x-axis shows the cumulative kernel densities and the y-axis depicts the model confidence (SVM certainty score). The numbers and percentages of cells corresponding to low (<0.3), medium (>0.3 and <0.7), and high (>0.7) certainty scores in each dataset are depicted next to the plots. ^/V indicates cell states most similar, and */# indicates cell states least similar to corneal cell states from the meta-atlas. Right, corneal markers among the top 10 explainable AI genes driving model's decisions for each of the predicted corneal cell states determined by their positive SHAP values are shown.

Fig. 4.

Motif enrichment and prediction of TF binding in cell states of the human cornea: A) heatmaps of the motif scores of the top 10 enriched TFs (upper panel) and of TF expression levels (lower panel) for each cell state. The bottom panel shows the type of associated motif. B) Examples of the consensus motif and TF footprints from 5A. Names of cell states: LSC-1, limbal stem cells 1; LSC-2, limbal stem cells 2; LSE, limbal suprabasal epithelium; CE, central epithelium; Cj, conjunctiva; qSK, quiescent stromal keratocytes; CF, corneal fibroblasts.

Fig. 5.

Prediction of key TFs controlling corneal cell states using GRN analysis of scANANSE. A) An in vitro single-cell human embryonic stem cell (hESC) network was used as a reference against networks of all individual cell states in the corneal cell state meta-atlas. Black boxes enclosing peaks in open chromatin depict putative enhancers where TFs can bind. B) Heatmap of influence scores (>0.8) of TFs shared between limbal/corneal epithelial and stromal cell states, predicted by scANANSE. C) Heatmap of influence scores (>0.8) of TFs specific for limbal/corneal epithelial cell states, predicted by scANANSE. D) Heatmap of influence scores of the top 40 TFs specific for stromal cell states, predicted by scANANSE. Names of cell states: LSC-1, limbal stem cells 1; LSC-2, limbal stem cells 2; LSE, limbal suprabasal epithelium; CE, central epithelium; Cj, conjunctiva; qSK, quiescent stromal keratocytes; CF, corneal fibroblasts.