Introduction

Aqueous humor (AH) dynamics and pressure homeostasis are important for proper ocular inflation and vision. Elevated intraocular pressure (IOP) is a key risk factor for glaucoma, a serious blinding disorder that affects 80 million people worldwide (Quigley and Broman, 2006). IOP is controlled by continually adjusting resistance to drainage (outflow) of AH through the trabecular meshwork (TM) and Schlemm’s canal (SC), the primary route of AH outflow from the eye. AH percolates through the TM before exiting the eye through SC and then entering the venous circulation. As the most abundant constituent of the conventional outflow pathway, TM cells are central in regulating AH drainage and IOP. Despite significant advances over the past few decades, a much deeper understanding of the molecular mechanisms governing TM cell health and function, as well as the pathological alterations responsible for increased outflow resistance, elevated IOP, and increased glaucoma risk are still needed.

The TM is located at the inner aspect of the wall of the eye at the iridocorneal angle just anterior to where the iris and cornea meet (i.e. the limbus). In humans, the TM is morphologically separated into distinct zones, including the uveal, corneoscleral, and juxtacanalicular trabecular meshwork (JCT), which lies adjacent to the SC. AH first drains through the uveal and corneoscleral meshwork whose TM cells have endothelial-like properties (Buffault et al., 2020; Keller and Peters, 2022; Stamer and Clark, 2017). In these regions, the TM consists of beams and plates made of various extracellular matrix (ECM) components including fibrillar collagens (Abu-Hassan et al., 2014; Acott and Kelley, 2008; Fuchshofer et al., 2006; Vittal et al., 2005; Yue, 1996). The beams and plates are covered with monolayers of TM cells. AH drains through a tortuous series of intertrabecular spaces, which run between the beams, ensuring interaction with TM cells. Next, AH flows into the less structured JCT region before entering the SC. The JCT consists of TM cells that are fibroblast-like extending their processes onto neighboring JCT and SC cells. JCT TM cells are embedded within a diffuse ECM ground substance (including glycosaminoglycans, Johnson, 2006; Keller and Acott, 2013; Stamer and Clark, 2017; Vranka et al., 2015) and have contractile smooth muscle-like properties (Stamer and Clark, 2017; Tian et al., 2009). Importantly, cells in the JCT region impose critical resistance to AH outflow with a major site of resistance to outflow being located in the tissue where the JCT and SC meet (Acott and Kelley, 2008; Johnson, 2006; Overby et al., 2009; Tamm, 2009).

Regulation of AH outflow resistance and IOP homeostasis by TM cells is complicated. AH outflow is modulated by: 1) The physical properties (shape/volume and contractility) of the TM cells. For example, some TM cells have smooth muscle-like properties that adjust outflow resistance and thus IOP by controlling cellular contractility (Rao et al., 2005; Syriani et al., 2009; Thieme et al., 2000; Tian and Kaufman, 2005; Zhang and Rao, 2005, Rao et al., 2001; Wang et al., 2018; Yu et al., 2008). 2) The effects of different TM cell types on both the composition and abundance of ECM in the JCT due to ECM synthesis/ metabolism and phagocytic activity. 3) The paracrine or physical effects of TM cells on SC and possibly the vasculature just distal to SC (Balasubramanian et al., 2024; Kizhatil et al., 2014; Stamer and Acott, 2012; Thomson et al., 2021).

Due to their importance, a more detailed molecular characterization of TM cells is still required. The degree to which specific pathways and modulatory roles are compartmentalized or shared by specific TM cell subtypes and how specific subtypes influence each other and the SC are not known. In fact, the number of molecularly distinct subtypes is currently unknown as is the role of the > 120 genes associated with elevated IOP within these subtypes (Gharahkhani et al., 2021). Although the TM has fewer beams and a narrower drainage path in mice, developmental, functional and anatomical similarities between human and mouse TM make mice a valuable model for studying the molecular characteristics, heterogeneity, subtypes and distribution of TM cells (Chen et al., 2016; Smith et al., 2001).

Single-cell resolution transcriptomic profiling is revolutionizing the molecular characterization of cell types and cellular heterogeneity (single-cell RNA sequencing, scRNA-seq; and single-nucleus, snRNA-seq). Initial single-cell-level studies have characterized cell types of the ocular anterior segment, including TM cells (Balasubramanian et al., 2024; Patel et al., 2020; Thomson et al., 2021; van Zyl et al., 2020; van Zyl et al., 2022), but further depth, validation, and more uniform molecular consensus/ naming of cell types are needed. This is especially true for TM cells. Currently, reports differ in the naming and definition of TM cell subtypes and offer limited validation by immunohistochemistry (IHC) or immunofluorescence (IF) on tissue sections or flat mounts.

Here we more deeply characterized the mouse TM and validated findings by immunostaining tissue sections and flat mounts and by in situ hybridization (ISH). We independently isolated and sequenced limbal cells of mice at 2 institutions and characterized the same 3 TM cell subtypes at both. Additionally, single-nucleus assay for transposase-accessible chromatin with sequencing (snATAC-seq) identified LIM homeobox transcription factor 1 beta (LMX1B) as one of the most active transcription factors in TM cells. Variants in LMX1B cause IOP elevation and glaucoma (Choquet et al., 2018; Gao et al., 2018; Gharahkhani et al., 2021; Khawaja et al., 2018; MacGregor et al., 2018; Shiga et al., 2018). By analyzing Lmx1b^V265D/+ mutant mice (develop elevated IOP and glaucoma, Cross et al., 2014; Tolman et al., 2021), we show that mitochondrial pathways are primarily disturbed in the TM3 cell subtype, which most strongly expresses Lmx1b. Treatment with nicotinamide lessened IOP elevation in Lmx1b^V265D/+ mutants supporting testing of therapeutics that boost mitochondrial health and function in human patients.

Results

Three molecularly distinct TM cell subtypes

Using scRNA-seq, we sequenced 17,914 cells from dissected limbal tissue (contains drainage structures) of 2-month-old mice. The dataset included 13,251 cells from strain C57BL/6J (B6) and 4,663 cells from strain 129/Sj (129). The data for both strains were integrated (Figure S1). Computational analysis revealed six distinct clusters of cells that were TM-containing, epithelia, pigmented epithelia (iris and ciliary body), endothelial, immune and neuron (Figure 1A). The identity of each of these cell clusters was based on various well-characterized marker genes (Figure S2A-B, Balasubramanian et al., 2024; Thomson et al., 2021; van Zyl et al., 2020). There were no major differences in the distributions of strain B6 and strain 129 cells within the relevant cell-type clusters (Figure S3A-C). Therefore, to improve statistical power, downstream analyses used the integrated dual-strain dataset.

Figure 1:

Supplemental Figure 1:

Supplemental Figure 2:

Supplemental Figure 3:

Supplemental Figure 4:

Supplemental Figure 5:

Supplemental Figure 6:

Cells of cluster 1 expressed various TM genes, including Acta2 (α-SMA), Pitx2, Tfap2b, and Myoc (Figure S2A-B, Akula et al., 2020; Balasubramanian et al., 2024; de Kater et al., 1992; Ostojic et al., 2008; van Zyl et al., 2020). Unbiased sub-clustering of Cluster 1 revealed 10 distinct sub-clusters (Figure 1B), which included three TM cell clusters based on signature gene expression (Figure 1D-E). We named these clusters TM1, TM2, and TM3. Hierarchical clustering of these cell types indicated that TM1 and TM2 are more molecularly similar to each other than to TM3 (Figure 1C). The remaining seven clusters were identified as cells from other ocular tissue including the iris stroma, sclera, corneal keratocytes, corneal endothelium, ciliary muscle, pericytes, and Schwann cells (Figure S2C, Liu et al., 2003; Lopez et al., 2009; Monje et al., 2018; Nitzan et al., 2013; Patel and Parker, 2015; Stratton et al., 2017; Toyono et al., 2015; van Zyl et al., 2022).

Next, we integrated an additional limbal tissue scRNA-seq dataset (14,912 cells, Duke University) from 3 months old C57BL/6J mice with our initial 2 months old C57BL/6J dataset. Due to batch effects most likely based on differences in cell isolation techniques and environment (Columbia University vs. Duke University, see Methods), TM cells from these datasets occupied adjacent, partially overlapping, but not identical UMAP space (Figure S4A). Individual analysis of the Duke dataset also identified 3 TM cell subtypes with the same marker gene expression and enriched molecular pathways compared to the Columbia University data (Figure S4B-E, Table S1). In datasets from previous mouse limbal scRNA-seq studies, there was also a high degree of overlap of marker gene expression with our 3 TM cell clusters (Figure S5A-B, Table S2, Thomson et al., 2021; Ujiie et al., 2023; van Zyl et al., 2020). However, some TM subtype clusters reported in previous studies expressed genes that overlapped with both of our TM2 and scleral cell clusters. (Figure S5B-C, Thomson et al., 2021; van Zyl et al., 2020). To resolve this discrepancy, we performed IF for two of the discordantly annotated markers previously reported to be expressed in TM cells but present in our scleral cluster (CD34 and LY6C1). This immunolabeling showed that these markers are expressed in scleral but not TM cells (Figure S6A-B). In addition, in previous studies, cells with markers matching TM3 (e.g. Lypd1) were named uveal cells without IF confirmation (Thomson et al., 2021; Ujiie et al., 2023; van Zyl et al., 2020). IF revealed that our TM3 cells are neither located in the uvea nor abundant in the uveal adjacent TM but primarily reside in the anterior TM closer to the cornea (see below).

TM subtypes differentially occupy TM zones

To further explore the sub-anatomical localization of the TM cell subtypes, we used a combination of IF and ISH to identify subtype markers. Each TM cell subtype had a higher percentage of cells expressing, and an increased average expression of, its respective subtype markers (all P < 1E-100) compared to other TM cells, Figure S7A). To analyze the distribution of subtypes, we first divided the TM into zones along its anterior/posterior and inner/outer axes (Figure S7B-E). The total area of TM occupied by each individual subtype marker was assessed in each zone on more than 150 tissue sections. This detected biases for the subtypes to be differentially localized in specific TM zones (Figure 2A-F). For instance, TM1 cells were significantly more frequent in the posterior and outer TM zones, while TM3 cells were overrepresented in the anterior and inner zones (all P < 0.01). The different marker molecules that were assessed for each TM cell subtype gave highly consistent results regarding zonal occupancy (Figure S8A-B, Figure S9). In addition, we confirmed the localization biases in 3D whole mounts of the mouse limbus (Kizhatil et al., 2014, for TM1: MYOC and TM3: α-SMA markers; Figure 2G-H). Collectively, these data indicate clear localization biases that differ for the TM1 and TM3 subtypes.

Figure 2:

Supplemental Figure 7:

Supplemental Figure 8:

Supplemental Figure 9.

To better understand the location of TM cell subtypes, we analyzed marker expression patterns in 50-60 exceptionally high-quality sections for each marker of each TM cell subtype (see Methods). For this analysis, the TM was divided into eight regions (Figure S7F). Not only did this refine the localization bias of TM1 and TM3 cells, but it discovered that TM2 cells are most concentrated in the mid-posterior two-thirds of the TM (mid regarding inner-outer axis, all P < 0.01, Figure 2I-K). Despite their biased distributions, some cells of each TM subtype were detected in most TM regions.

Molecular comparison of TM cell subtypes

To determine the molecular functions of TM cells, we first conducted gene ontology (GO) analysis using genes that were differentially expressed genes (DEGs) in all TM cells as a group compared to other sequenced cells (Table S3). Additionally, we used gene set enrichment analysis (GSEA) to assess whether these GO pathways are relatively enriched or underrepresented in TM cells. Compared to other cells, TM cells were enriched for pathways associated with extracellular matrix (ECM) structure and function, including signaling associated with collagens, proteases, integrins, and glycosaminoglycans (Figure 3A, Figure S10A, Table S4). TM cells were also enriched for growth factor signaling. Conversely, underrepresented pathways in TM cells were associated with desmosomes, peroxisomes, and certain cytoskeletal and plasma membrane elements. These pathways help define TM cells relative to other cells in the region.

Figure 3:

Supplemental Figure 10:

Next, we compared TM cell subtypes to each other (Figure 3B-D). TM1 cells were most enriched for ECM pathways, particularly structural components such as collagens, integrins, and other ECM elements (Figure 3B, Figure S10B). The expression of several ECM genes that are reported to impact IOP are heightened in TM1 cells, including Type VIII collagens, fibronectin, and Ltbp2 (Ali et al., 2009; Desronvil et al., 2010; Roberts et al., 2020). Both TM1 and TM2 cells have enriched expression of glycosaminoglycan (GAG) genes (Figure 3B-C, Figure S10B-C). GAGs in the JCT and intertrabecular spaces form a gel-like consistency and changes in GAG abundance correlate with altered AH outflow resistance (Knepper et al., 1981; Knepper et al., 1996; Knepper and McLone, 1985). TM2 cells are uniquely enriched for the laminin complex, a key component of basement membranes (Aumailley, 2013). Pathways related to clearing debris and monitoring the extracellular environment including phagocytic vesicles, antigen presentation, and lysosomal function are also enhanced in TM2 cells. In contrast, TM3 cells are underrepresented for ECM-related pathways compared to other TM cells but are enriched for pathways related to actin binding and mitochondrial metabolism (Figure 3D, Figure S10D). The actomyosin system plays a crucial role in maintaining TM cell contractility and in regulating outflow facility (Rao et al., 2005; Syriani et al., 2009; Thieme et al., 2000; Tian and Kaufman, 2005; Zhang and Rao, 2005). The enrichment of mitochondrial metabolism pathways in TM3 cells emphasizes their energetic needs and their expected heightened susceptibility to genetic and environmental contexts that compromise metabolic functions. Importantly, genetic variation in mitochondrial/ metabolic pathways contributes risk for IOP elevation and glaucoma (Aboobakar et al., 2023; Khawaja et al., 2016; Zhang et al., 2023).

We next compared growth factor signaling and ligand-target interactions between TM cells subtypes and other cells. TM1 cells were enriched for transforming growth factor beta (TGF-β) pathway signaling (Figure 3B). Receptor ligand activity was enriched in TM2 cells while the membrane raft pathway (a site of signal transduction) was enriched in TM3 cells (Figure 3C-D). Next, we used LRLoop (Xin et al., 2022), which was developed based on NicheNet (Browaeys et al., 2020), to predict significant ligand-target interactions between TM cells and other local endothelial cells that are relevant to AH drainage (Figure 4A-D). TGF-β ligands (Tgfb1, 2, & 3) were predominantly expressed in TM1 cells. By contrast, TM2 cells expressed molecules critical for vascular function, such as Angpt1, Vegfa, and Edn3 (Figure S11A-D, Genovesi et al., 2022; Saharinen et al., 2017; Shibuya, 2011). These TM-secreted molecules have been directly associated with SC maintenance and IOP (Comes and Borras, 2009; Fujimoto et al., 2016; Reina-Torres et al., 2017; Thomson et al., 2020).

Figure 4:

Supplemental Figure 11:

Comparison to human TM cells

We next worked to relate our findings to human data. We compared marker genes for each of our TM cell subtypes to sequenced human TM cells from two previous studies (Patel et al., 2020; van Zyl et al., 2020). These previous studies differed in the number of human TM cell subtypes that they reported and in the identities of their clusters. Van Zyl et al. identified three TM cell subtypes that they called JCT, beam A, and beam B (van Zyl et al., 2020), whereas Patel et al. identified two subtypes that they called fibroblast-like and myofibroblast-like (Patel et al., 2020). TM1 marker genes were most enriched in Van Zyl et al’s human JCT cells (Figure S12A). Similarly, Dcn was expressed in mouse TM1 cells and was localized to fibroblast-like TM cells adjacent to SC by Patel et al. TM1 cell nuclei were typically shorter with greater sphericity than the elongated endothelial-like nuclei of TM3 cells (Figure S13A-D). This nuclear morphology was more consistent with a JCT than beam cell identity. Together, these data support a JCT identity for TM1 cells. TM2 cells were more closely related to human beam A and beam B cells than to the JCT cells reported by Van Zyl et al (Figure S12B). TM2 cells also shared expression of the RSPO2 marker gene with the human myofibroblast-like cells reported by Patel et al. Together, this suggests that TM2 cells are myofibroblast-like beam cells. TM3 cells also shared marker genes with Patel et al’s myofibroblast-like cells (TAGLN and ACTA2), suggesting TM3 cells are also beam cells. There was relatively equal TM3 marker gene enrichment across all human TM cell subtypes of Van Zyl et al (Figure S12C). Nuclei immediately adjacent to TM2 and TM3 cell marker staining had an elongated, endothelial-like shape (Figure S13A-C). This nuclear morphology was also consistent with the location of TM2 and TM3 cells on beams and analogous to the endothelial-like morphology of human beam cells (Stamer and Clark, 2017).

Supplemental Figure 12:

Supplemental Figure 13:

Next, we evaluated the expression of genes implicated in elevated IOP and glaucoma by genome-wide association studies (GWAS). TM1 and TM2 cells had similar expression levels of GWAS genes associated with IOP elevation and primary open-angle glaucoma (POAG) respectively (Figure S14A). TM3 cells had slightly less expression of these GWAS genes compared to TM1 and TM2. To better understand the genes driving GWAS gene enrichment in different TM subtypes, we studied pathways that are significantly enriched for GWAS genes (Figure S14B-C, see Table S6 for complete list). TM1 cells were the most enriched for genes involved in ECM organization, cell-substrate adhesion, and growth factor stimulus (Figure S14D). All 3 TM subtypes had similar expression of GWAS genes associated with hormone stimulus, Vegf signaling, and Rho kinase signaling. TM3 cells expressed lipid and carbohydrate mitochondrial metabolism pathway genes at higher expression levels and in more cells than did TM1 and TM2 cells (Figure S14E). The same was true for actin binding-associated GWAS genes in TM3 cells versus TM1 and TM2. Genetic variation in lipid and carbohydrate mitochondrial metabolism pathways is associated with IOP elevation and glaucoma (Khawaja et al., 2016).

Supplemental Figure 14:

Regulation of TM cell gene expression

To better understand the overlap and differences in transcriptional control of TM cell gene expression, we performed single-cell resolution multiome sequencing. Using limbal strip tissue isolated from a separate cohort of mice than those used for scRNA-seq, we simultaneously performed both single nucleus RNA sequencing (snRNA-seq) and snATAC-seq. We thus captured both gene expression and open chromatin data throughout the genome in the same cell. Using the snRNA-seq data, we clustered cells into various cell types. The generated clusters and their marker genes closely agreed with our scRNA-seq analyses (Figure 5A-B, Figure S15A). Notably, there was a significant correlation between gene expression (RNA-based clustering) and chromatin accessibility (ATAC-based clustering) with adjusted Rand index of 0.20 (p < 0.001, permutation test, Figure S15B-D). This included an increase in chromatin accessibility at the promoters of marker genes in their respective TM cell subtypes (Figure 5C-E). To discover the key transcriptional regulators within TM cells, we examined the correlation between the expression of each transcription factor (TF) and their chromatin accessibility at predicted binding sites within all TM cells. Various TFs that either activated or repressed transcriptional activity exhibited good correlations with gene expression across TM cells (Figure 5F). Such activating TFs include Tcf21, Arid3b, and Myc, while repressing TFs include Meox2, Junb, and Sox6. Notably, LMX1B, an important gene in glaucoma GWAS (Choquet et al., 2018; Gao et al., 2018; Gharahkhani et al., 2021; Khawaja et al., 2018; MacGregor et al., 2018; Shiga et al., 2018), is among the most strongly activating TFs in TM cells.

Figure 5:

Supplemental Figure 15:

Metabolism-supporting treatment lessens IOP elevation in Lmx1b^V265D mice

We evaluated the utility of our new TM cell atlas to better understand the functions of Lmx1b in TM cells and how Lmx1b mutations lead to glaucoma. To do this, we analyzed scRNA-seq data generated for limbal cells from mice with a dominant mutation in Lmx1b (Lmx1b^V265D/+). These mutant C57BL/6J mice develop early onset IOP elevation (Cross et al., 2014; Tolman et al., 2021; Zhang et al., 2024). We compared gene expression in V265D mutant TM cells after integration with our wild-type data (all data B6 background and postnatal day 60, Figure 6A, Figure S16A). Results show that Lmx1b is minimally expressed in TM1 and TM2 but is highly expressed in TM3 cells (Figure 6B). Next, we analyzed DEGs and pathway differences between WT and V265D/+ cells across each TM cell subtype. Across all TM subtypes, various genes related to ECM function and growth factor signaling were downregulated in mutants, including glycosaminoglycans and insulin growth factor binding proteins (Figure 6C, Table S6 for complete list of pathways). Pathways enriched in mutant TM cells were associated with ribosome and calcium signaling. We then examined the differential responses to the Lmx1b mutation by analyzing pathways dysregulated in each TM subtype that were not common across all TM cells (Figure 6D-F, Figure S16B). V265D mutant TM1 cells had a downregulation of collagen synthesis/metabolism gene expression relative to WT cells, particularly in the production and turnover of fibrillar collagens. TM2 cells with Lmx1b mutations had upregulation of immune signaling pathways. The highly Lmx1b expressing TM3 cells exhibited perturbed mitochondrial metabolism, with a downregulation of genes in complex I and ATP synthase, and upregulation of genes involved in ubiquitin pathways. These results document altered mitochondrial function and proteostasis in TM3 cells.

Figure 6:

Supplemental Figure 16:

Our findings implicated perturbed metabolism within TM3 cells as responsible for IOP elevation in an Lmx1b glaucoma model. Nicotinamide (NAM) is known to boost nicotinamide adenine dinucleotide (NAD), supporting healthy metabolism/ mitochondrial function (Schondorf et al., 2018; Williams et al., 2017b). Thus, we hypothesized that the metabolic abnormalities in TM3 cells underlie the IOP elevation that causes glaucoma, and that NAM treatment may protect the Lmx1b mutant TM. Beyond LMX1B, this hypothesis is relevant more commonly to POAG treatment because metabolism-relevant genes are implicated by GWAS (and their expression is enriched in TM3 cells, Figure S14E). To test our hypothesis, we treated mice by adding NAM to their drinking water starting at postnatal day 2 and continuing into adulthood. Results showed that NAM supplementation significantly protected against anterior chamber deepening (ACD, a symptom of IOP elevation, Figure 7A-B) and IOP elevation (Figure 7C-D). These results support NAM as a treatment to prevent TM damage and IOP elevation in LMX1B-mediated glaucoma with potential general relevance to POAG.

Figure 7:

Discussion

This study provides a comprehensive multimodal analysis of mouse limbal cells that focuses on the TM. We defined three TM cell subtypes with unique marker panels that can be used for their identification across different strains, laboratories, and processing methods. We validated key marker gene expression by IF and ISH, discovering biases for specific TM subtypes to reside in different parts of the TM. For example, TM1 cells have a location bias for the outer, posterior TM while TM3 cells are more abundant in the inner, anterior TM. Importantly, our pathway analyses suggest that the different TM subtypes have overlapping roles. All 3 TM subtypes are enriched in ECM synthesis and ECM metabolism compared to other limbal subtypes. However, there is also compartmentalized enrichment of some processes to specific subtypes. For example, ECM production and turnover, secreted ligand signaling, debris removal, immune surveillance, actin-associated cell contraction and increased metabolic gene activity (as discussed below) are enriched in specific cell-subtypes. Using ATAC-seq data, we discovered key transcription factors controlling gene expression in TM cells, including LMX1B. Importantly, we implicate LMX1B in metabolic control in TM cells. More specifically, we demonstrate that Lmx1b is most active in TM3 cells, while mutation of Lmx1b induces prominent mitochondrial defects in TM3 cells. Finally, we show that therapeutically targeting TM3 dysfunction with a mitochondrial/ metabolism supporting nutrient (nicotinamide) protects from IOP elevation in the Lmx1b glaucoma model. We also suggest that a TM3-like cell is generally relevant to metabolic pathways in human glaucoma based on enriched expression of human glaucoma GWAS gene in TM3 cells. Together these data provide a wealth of new molecular information about TM cells and suggest a new approach for preventing IOP elevation.

Mouse trabecular meshwork subtype classification

Here we provide evidence for three TM cell subtypes in mice that are robust across datasets from different institutions, using different tissue processing and sequencing protocols. Although others (Thomson et al., 2021; Ujiie et al., 2023) subdivided TM cell subtypes beyond the three defined in the present study, at this point we decided against further subdivision (See Figure S5) for the following reasons. We lacked confidence in the higher resolution clustering, as it was not consistent across datasets, did not produce unique marker gene sets, and could have resulted in over-clustering due to stochastic or stress-induced differences between subpopulations of cells. Moreover, the number of reported mouse TM cell clusters in previous studies were inconsistent. Although three TM subtypes were reported by Van Zyl et al. and Ujiie et al. (Ujiie et al., 2023; van Zyl et al., 2020), the molecular features of cells in these clusters were not the same across these studies. Differences between studies may be due to batch effects, different computational methods and possibly over-clustering. Most importantly, there was very limited validation by IF, IHC, or other methods in previous studies. Only the JCT type of TM cell in these previous studies was assessed by means other than scRNA-seq (IF and ISH, van Zyl et al., 2020). Additionally, Thompson et al. treated cells with Y27632, a ROCK inhibitor, which could alter TM cell transcriptomes. It is also important to point out that both Thompson et al. and Van Zyl et al. used albino mouse strains, which can impact TM cell development and function (Libby et al., 2003). Finally, Ujiie et al. sequenced TM cells at 3-4 weeks, which is prior to full maturation of the TM structure (Smith et al., 2001), while our ongoing developmental studies show substantial molecular maturation of TM cells between P21 and P60.

As Van Zyl et al’s data is deposited in the Broad Institute’s single-cell portal, we further compared their data with ours. Although both Van Zyl et al. and our current paper identified 3 TM subtypes, our work does not fully agree. Some of the cells named beam A and beam Y by Van Zyl et al. express markers present in our TM2 cell subtype, while others express markers that are absent in TM2 cells but present in our scleral cells. Our IF studies confirm that these discordant markers are not expressed in TM cells, but instead show that they are in fact scleral cells. These scleral cells highly express fibroblast markers including Mfap5, Clec3b, and Tnxb, suggesting they are scleral fibroblasts. Cells initially named uveal cells by Van Zyl et al. (and subsequently by Ujiie et al. and Thompson et al.) express markers of our TM3 cells. Consistent with our gene expression data, our IF shows that these Van Zyl markers are expressed by cells that are primarily located in the TM region that is closer to the cornea and not in the uvea or the TM region adjacent to the uvea. Thus, our current study resolves some previous inconsistencies in classifying mouse TM cell subtypes. It lays a firmer foundation for continued understanding of TM cell types and their biology. Moving ahead, sequencing cells at greater depths will enhance TM cell characterization and improve comparisons between datasets.

Differential roles of specific mouse TM cell subtypes and comparison to human

As a group, all TM cells share various molecular properties that distinguish them from other ocular cells. Our analyses determine the expression of IOP, and glaucoma genes identified by GWAS in TM cell subtypes as well as in other limbal cells. Pathways enriched across all TM cell subtypes based on RNA-seq are largely related to ECM (including collagens, integrins, and glycosaminoglycan pathways), cell-cell signaling, and growth factor signaling. This fits with an essential role of the TM in both ECM synthesis/turnover and in signaling to the SC to maintain IOP (Balasubramanian et al., 2024; Fuchshofer and Tamm, 2009; Keller et al., 2009; Stamer and Acott, 2012; Thomson et al., 2021). Additionally, genes associated with intermediate filaments, which function in determining cell shape and structure and can also act as mechanical stress absorbers (Herrmann et al., 2007), are underrepresented in TM cells. Lower levels of intermediate filaments may enable TM cells to alter their shape to respond to changes in IOP.

Mouse TM1 cells closely resemble sequenced human JCT cells. Like human JCT cells, the majority of TM1 cells are located nearest to the SC inner wall, particularly adjacent to the posterior of the mouse SC. Human JCT cells are generally accepted to be embedded in a more diffuse ECM and to not cover the collagenous trabecular beams (Keller and Acott, 2013; Stamer and Clark, 2017; Vranka et al., 2015). In vivo, a proportion of TM1 cells have a shorter, more spherical nuclear morphology, consistent with JCT cell identity. Importantly, the JCT region is critical in determining resistance to AH outflow and in regulating IOP (Acott and Kelley, 2008; Ethier, 2002; Johnson, 2006; Stamer and Acott, 2012). JCT cells have fibroblast-like properties, including the secretion of ECM proteins and degradation enzymes to support continuous ECM remodeling (Acott et al., 1988; Keller et al., 2009; Stamer and Clark, 2017). The subset of TM1 cells that lie distant from SC have more beam cell-like, elongated nuclear morphology. This suggests that either our TM1 cluster contains more than 1 TM cell subtype or that the morphology of a single TM1 cell subtype is influenced by local environment. This needs to be resolved by future studies. Compared to the other TM cell subtypes, TM1 cell pathways are enriched in both ECM structural molecules, such as collagens and integrins, and TGF-β signaling, which increases ECM production (Munger and Sheppard, 2011). Notably, excessive TGF-β signaling and ECM production can lead to fibrosis and elevate IOP (Biernacka et al., 2011; Frangogiannis, 2020; Meng et al., 2016). These data suggest that TM1 cells play an important role in the constant remodeling of the ECM. However, dysregulation of these pathways can be pathogenic.

TM2 cells molecularly resemble previously sequenced human TM beam cells (Patel et al., 2020; van Zyl et al., 2020). This includes overlap of expression of top TM2 signature gene expression with sequenced human beam cells. Consistent with a beam cell identity, TM2 cells are located more towards the mid and inner TM than most TM1 cells. TM2 cells have an elongated nuclear morphology, consistent with the previously established beam cell morphology of cells in the regions they occupy (Overby et al., 2014; Smith et al., 2001; Stamer and Clark, 2017). Beam cells have endothelial-like properties, such as maintaining the patency of the AH outflow tract through the production of antithrombotic molecules. Proteases and molecules classically known as antithrombotic and thrombolytic maintain fluid flow through the TM by preventing build-up of protein aggregates to maintain fluid drainage (Acott and Kelley, 2008; Stamer and Acott, 2012). Such molecules include glycosaminoglycans, like heparin, the serine protease tPA and matrix metalloprotease that are all highly expressed in TM2 cells. Beam cells are also phagocytic and this activity along with turnover of ECM further cleans the TM (called filtering) to enhance fluid drainage (Stamer and Clark, 2017). Human Beam cells also participate in antigen presentation and modulation of inflammation through cytokines and major histocompatibility proteins (Shifera et al., 2010; Stamer and Clark, 2017), consistent with TM2 cells being enriched for these immune pathways/ molecules. Additionally, our data suggest that signaling from TM2 cells to SC endothelial cells is mediated by soluble molecules including Vegfa, Edn3, and Angpt1, which are critical for SC maintenance and function. This links TM2 cells to ocular development and glaucoma. VEGF receptor and ANGPT signaling genes modulate SC development/ function as well as contribute to developmental and primary open-angle glaucoma (Aspelund et al., 2014; Kabra et al., 2017; Khawaja et al., 2018; Kizhatil et al., 2014; MacGregor et al., 2018; Souma et al., 2016; Thomson et al., 2020).

TM3 cells also have a location and elongated nuclear morphology resembling beam cells with enriched localization in the inner, anterior TM. TM3 cell-enriched processes include actin binding, and modulation of the actomyosin system, which is known to modulate outflow resistance and IOP (Rao et al., 2001; Yu et al., 2008). As such, TM3 cells express a large percentage of glaucoma-associated genes involved in actin binding at a high level. Importantly, TM3 cells are enriched for various mitochondrial metabolic and antioxidant pathways, including those associated with IOP elevation and glaucoma (Aboobakar et al., 2023; Khawaja et al., 2016). Inhibiting metabolic pathways alters AH outflow resistance (Reina-Torres et al., 2020). Since maintaining contractile tone is an energy-intensive process (DeWane et al., 2021), this could explain the enrichment of energy metabolism pathways in TM3 cells compared to other TM cells.

TM3 cells have a unique susceptibility to mitochondrial dysfunction

In addition to being enriched in mitochondrial/ energy metabolism processes, TM3 cells express Lmx1b at significantly higher levels than both the other TM cell subtypes and other limbal cells. Importantly, in heterozygous mutant V265D/+ mice, TM3 cells had pronounced gene expression changes that implicate mitochondrial dysfunction and oxidative stress, but that were absent or much lower in other cells including TM1 and TM2. Together, these data suggest that mutation of Lmx1b has a primary effect on metabolism and oxidative stress in TM3 cells that subsequently leads to IOP elevation and then glaucoma. Relevantly, homozygous conditional knockout of both Lmx1b and Lmx1a as well as siRNA knockdown of Lmx1b alone impact mitochondrial functions in neurons and may contribute to Parkinson’s disease (Bergman et al., 2009; Doucet-Beaupre et al., 2016; Jimenez-Moreno et al., 2023). Our data extend these published findings by showing that inheritance of a single dominant mutation in Lmx1b similarly affects metabolism in TM3 cells. In addition to modulating mitochondria, LMX1B was recently identified as a transcription factor that regulates responses to cell stress. This includes promotion of autophagy under stress conditions to enable recycling of cellular resources to maintain cellular functions and enable survival (Jimenez-Moreno et al., 2023). Thus, mutation of Lmx1b may also result in a deficiency of the beneficial autophagic response to stress (Jimenez-Moreno et al., 2023).

Although further experiments are needed, dysfunctional autophagy may further prolong and exacerbate TM cell stress, resulting in more extensive depletion of cellular resources, exacerbated cellular dysfunction, and IOP elevation. Other consequences of the Lmx1b mutation that were evident in the V265D/+ cells may also contribute to IOP elevation, including the overexpression of genes associated with ribosomes and calcium signaling, and depletion of genes involved with ECM synthesis and metabolism. However, these processes were not limited to TM3 cells or even to cell types that express detectable Lmx1b, suggesting that they are secondary damaging processes that are subsequent to the initiating, Lmx1b-induced, metabolic disturbances and cell stress in TM3 cells. Thus, we hypothesize that mitochondrial abnormalities in TM3 cells are of primary importance to IOP elevation in V265D/+ glaucomatous mice.

Nicotinamide protects against IOP elevation in Lmx1b^V265D-mutant mice

Based on our findings in TM3 cells, we tested if mitochondrial and metabolic disturbances are important in IOP elevation by administering nicotinamide (NAM). NAM is a form of vitamin B3 that is well-established to boost NAD⁺ levels, promote mitochondrial health and energy metabolism, relieve oxidative stress and thus promote cellular resilience (Gomes et al., 2013; Imai and Guarente, 2014; Jang et al., 2012; Mouchiroud et al., 2013; Verdin, 2015; Williams et al., 2017b, Jang et al., 2012; Kang and Hwang, 2009; Song et al., 2021). Supporting our hypothesis of primary metabolic dysfunction in TM3 cells, NAM treatment strongly protected V265D/+ mice from IOP elevation. Thus, treatments that support normal metabolism may protect from IOP elevation in individuals with LMX1B variants including both developmental glaucoma and POAG patients (Chen et al., 1998; Choquet et al., 2018; Choquet et al., 2017; Gao et al., 2018; Gharahkhani et al., 2021; Khawaja et al., 2018; MacGregor et al., 2018; Sweeney et al., 2003; Vollrath et al., 1998). As NAM and other NAD⁺ boosting treatments are also known to be directly neuroprotective in glaucoma (De Moraes et al., 2022; Hui et al., 2020; Williams et al., 2018; Williams et al., 2017a; Williams et al., 2017b), our current data suggest that such treatments will have a double benefit for patients. As mitochondrial defects are becoming broadly associated with glaucoma risk and progression (Abu-Amero et al., 2006; Lo Faro et al., 2021; Petriti et al., 2024), these treatments could have wide-ranging potential to complement and augment existing IOP-lowering therapeutics through a completely distinct mechanism.

Our results provide a thorough molecular characterization of mouse TM cells, providing much needed molecular information on subtype specializations in regulating IOP and the roles of specific subtypes in glaucoma. This comprehensive TM atlas provides a new foundation to guide development of new glaucoma treatments. We validate its utility by implicating metabolic changes in a single TM cell subtype in both mouse and human glaucoma, with success of a metabolic treatment targeting this subtype in a glaucoma model.

Methods

Animal husbandry and ethics statement

Experimental animals were either C57BL/6J (Stock# 664) or 129/Sj (unique substrain of 129, Tolman et al., 2021) inbred mice. The Lmx1b^V265Dmutation was discovered in an N-ethyl-N-nitrosourea (ENU) mutagenesis screen (Cross et al., 2014; Thaung et al., 2002; Tolman et al., 2021). These mice were then backcrossed to the C57BL/6J (Stock #000664) for at least 30 generations. All mice were treated in accordance with the Association for Research in Vision and Ophthalmology’s statement on the use of animals in ophthalmic research. The Institutional Animal Care and Use Committee of either Columbia University or Duke University approved all experimental protocols performed at their respective institutions.

NAM administration

Nicotinamide (NAM or niacinamide; Sigma-Aldrich, St. Louis, Missouri) was dissolved in the standard institutional drinking water to a dose of 550 mg/kg/day (low dose) based on the average volume mice consume. Untreated groups received the same drinking water without NAM. Water was changed once per week. Treatment was started at postnatal day 2. Births were checked daily between 9am-12pm to determine the pup’s age.

Preparation of ocular cells

Columbia animals: Eyes were pooled from 4 animals to generate 1 sample (8 eyes total per pool). Two independent samples were processed. Mice were 2 months of age. Sex was balanced in each group.

Duke animals: One sample was generated from a pool of 4 animals (8 eyes total). All mice were male. Mice were 3 months of age.

Dissections: After euthanizing the animals, eyes were enucleated into fresh Dulbecco’s Modified Eagle Medium (DMEM; Gibco) and subsequently dissected in DMEM. All dissection tools and surfaces were treated with RNAse zap (Ambion). For dissections, we removed the posterior eye by cutting between the limbus and sclera. The lens was then removed. Finally, we made a small circular cut in the anterior cornea, removing the middle 1/3^rd of the area. Only limbal tissues were collected and all other tissues were removed. The remaining limbal tissue was minced and prepared for digestion.

Single Cell RNA Sequencing

Single-cell RNA sequencing (Columbia University): Tissues were digested enzymatically using Papain and Deoxyribonuclease I (Worthington Biochemical Corporation, Cat. no. LK003153) for 20 min at 37°C and stopped using Earl’s balanced salt solution (EBSS, Thermo Fisher Scientific, Cat no. 24010-043). Cells were triturated using an 18-gauge needle, centrifuged at 300g at 4 °C, washed with cold DMEM, and filtered using a 100 μm Flowmi Cell Strainer (H13680-0040, Bel-Art. Cells were resuspended in cold DMEM and placed on ice immediately. Cells were counted using Countess II automated cell counter and processed for 10x library preparation immediately (see below). All data will be deposited in the Gene Expression Omnibus (GEO) database and single cell portal.

Single-cell RNA sequencing (Duke University): Single cell dissociation was performed using Collagenase IV (Worthington Biochemical Cat. no. LS004188), Dispase II (Sigma Cat. No. D4693) and DNAse I (Roche Cat no. 50225500) for 60 min and then trypsin/EDTA for 10 min at 37°C. Cells were triturated using a 1 ml pipette tip, centrifuged at 300g at 4°C, washed with cold DMEM containing 10% FBS, and filtered using a 70 μm cell strainer (Fisher). Cells were resuspended in cold DMEM containing 5% FBS and placed on ice immediately and sent to Duke Genome Core Facility. Cells were counted using Countess II automated cell counter (the cell viability was 88%) and then immediately processed immediately (see below). All data will be deposited in the Gene Expression Omnibus (GEO) database and single cell portal.

Multiome (Columbia University): Limbal strip samples for multiome analysis were snap-frozen on dry ice immediately following dissection and stored at −80C until processing. To isolate nuclei, tissue was homogenized in 100 μL of chilled lysis buffer (10 mM Tris-HCl, pH 7.4, 10 mM NaCl, 3 mM MgCl2, 0.1% NP40, 0.1% Tween20, and 0.01% digitonin – from (Corces et al., 2017) supplemented with 1% BSA). Afterwards, samples were incubated for 5 min on ice, gently agitated by pipetting, and incubated for an additional 10 minutes on ice. We then added 1 mL of chilled wash buffer (10 mM Tris-HCl, pH 7.4, 10 mM NaCl, 3 mM MgCl2, 0.1% Tween 20, 1% BSA) to the cells and centrifuged at 500 rcf for 5 min at 4°C. Cells were resuspended in chilled Diluted Nuclei Buffer (PN-2000153, 10x Genomics) and filtered through a 40 μm Flowmi Cell Strainer (H13680-0040, Bel-Art). Cells were resuspended in cold DMEM and placed on ice immediately. An aliquot of nuclei suspension was stained with SYTOX green to count nuclei using Countess II automated cell counter and processed for 10x library preparation immediately (see below). All data will be deposited in the Gene Expression Omnibus (GEO) database and single cell portal.

Single-cell and single-nucleus RNA sequencing

Single-cell and single-nucleus RNA-seq was performed at the single-cell sequencing core in Columbia Genome Center (Columbia University) or Duke Molecular Genomics Core (MGC, Duke University). Single cells or nuclei were loaded into chromium microfluidic chips with v3 chemistry and barcoded with a 10x chromium controller (10x Genomics). RNA from the barcoded cells was reverse-transcribed, and sequencing libraries were constructed with a Chromium Single Cell v3 reagent kit (10x Genomics). Sequencing at both institutions was performed on NovaSeq 6000 (Illumina).

Sequencing data processing

Single cell and single nucleus sequencing

Raw reads were mapped to the mm10 reference genome by 10x Genomics Cell Ranger pipeline. “Seurat” (Version 4.0.1, Hao et al., 2021) was used to conduct all single-cell and single-nucleus sequencing analyses. Briefly, the dataset was filtered to contain cells with at least 200 expressed genes and genes with expression in more than three cells. Cells were also filtered for mitochondrial gene expression (<20% for single-cell and <5% for single-nucleus). The dataset was log-normalized and scaled. Biological incompatibility based on gene expression was used to identify doublets. For the scRNA-seq dataset generated at Columbia scRNA-seq, we performed unsupervised clustering followed by comparison to previous annotations (Balasubramanian et al., 2024; Patel et al., 2020; Thomson et al., 2021; van Zyl et al., 2020; van Zyl et al., 2022). A resolution parameter of 0.03 was used to match these data and previous annotations. Cluster 1 was then subclustered using unsupervised methods with a resolution parameter of 0.3. All cluster identities based on a combination of known marker genes and validation experiments. For additional datasets (Duke scRNA-seq, and Columbia snRNA-seq), parameters were manually adjusted cluster cells based on the Columbia scRNA-seq dataset. Marker genes were identified using Wilcoxon test implemented in Seurat using default parameters.

Multiome

Single-nucleus multiome data were quantified by 10x CellRanger and preprocessed according to integratedSeurat v4 (Hao et al., 2021) and Signac (Stuart et al., 2021). We filtered cells which have <500 and >80000 ATAC reads, have <500 and >25000 RNA reads, and have >10% mitochondrial reads. RNA-seq processing was performed using the same methods as scRNA-seq (above). ATAC-seq processing was performed using ArchR (version 1.0.2) (Granja et al., 2021), where we used default latent semantic indexing (LSI) approach to derive low-dimensional embedding followed by single-cell neighborhood construction and graph-based clustering using the top 50 LSI dimensions excluding the first dimension. We assigned cell types using similar process based on marker peak accesibility and transferred cell type labels to the ATAC data.

Lmx1b^V265D analysis

Data from limbal cells with the Lmx1b^V265D mutation (C57BL/6J background) were integrated with Wild-Type (WT) cells from the C57BL/6J strain only. Multiple approaches were used for normalizing data for integration in Seurat, including log normalization and SC transform. The relative number of cells in each cell cluster between WT and V265D mutant cells changed based on the integration technique. Therefore, we chose not to assess TM or other cell proportions between genotypes. Log normalization was selected for downstream analysis based on the log-transformed integrated WT clusters being more molecularly similar to the WT only clusters compared with the SC transform results.

Statistical data analysis

Hierarchical clustering

TM-containing cells were assessed by hierarchical clustering based on the similarity of the expression of the top 500 marker genes (ranked by P value) of each individual cluster. This analysis filtered out low-expressing genes (expressed in >10% of cells, with a logFC > 0.25 compared to other cells).

GO analysis

Gene ontology and pathway enrichment analysis were performed by comparing the differentially expressed genes in three scenarios:

1. All TM cells versus all other sequenced cells.

2. Each individual TM cell subtype against other TM cell subtypes.

3. WT versus Lmx1b^V265D/+ TM cells.

For comparisons of all TM cells to all sequenced cells, pathway enrichment was compared against a background expression universe of genes expressed in any sequenced cell. For comparisons of individual TM cell subtypes, only genes expressed in trabecular meshwork cells were used as a background (P value cut-off 0.01, Timmons et al., 2015).

Gene set enrichment analysis (GSEA)

was performed using the same set of differentially expressed genes from the above comparisons with the GseGO function. All pathways were assigned an enrichment score (< 0, underrepresented; > 0, enriched). Pathways that were significantly enriched or underrepresented were also analyzed (P value cut-off 0.01). The R package ClusterProfiler was used for analysis (Yu et al., 2012).

Module score

From genome-wide association studies of IOP elevation (Blue Mountains Eye and Wellcome Trust Case Control, 2013; Choquet et al., 2017; Gao et al., 2018; Hysi et al., 2014; Khawaja et al., 2018; MacGregor et al., 2018; Nag et al., 2014; Ozel et al., 2014; Springelkamp et al., 2015; Springelkamp et al., 2017; van Koolwijk et al., 2012) or primary open-angle glaucoma (Bailey et al., 2016; Bonnemaijer et al., 2018; Choquet et al., 2018; Craig et al., 2020; Gharahkhani et al., 2021; Shiga et al., 2018; Verma et al., 2024), we curated a list of genes associated with the disease. Seurat’s AddModuleScore function was applied to assess glaucoma disease risk for each cell type, using default parameters (24 bins for aggregation with 100 control features per analyzed feature).

Predicted ligand-receptor interactions

Predicted ligand-target links between interacting cells were identified using LRLoop (Xin et al., 2022), which was developed based on NicheNet (Browaeys et al., 2020). Briefly, the expression of genes in cell types is linked to a database of signaling and gene regulatory networks curated from prior information, enabling viable predictions of potential interactions between cell types. These interactions were calculated using target molecules from one cell type and ligands from multiple cell types. The top interactions, regardless of the cell type origin of the ligand, are represented in a Circos plot. Ligands were assigned to an individual TM cell subtype/subtypes based on the following criteria:

1. Exclusively assigned to a TM cell subtype if expressed > 1 transcript per ten thousand in a given cell type and expressed on average 4X more in that subtype compared to all other TM cell subtypes.

2. Assigned to multiple TM cell subtypes if expression is > 1 transcript per ten thousand in multiple cell types and average ligand expression is < 4X between those TM cell subtypes.

3. If ligand expression is low (< 1 transcript per ten thousand) across TM cell subtypes, the ligand is assigned to all subtypes that express the ligand.

Comparison to published datasets

To compare TM cell clusters across datasets, the top 10 marker genes (by P value) for TM cell subtypes in our dataset were compared to published marker genes in other TM cell datasets (Patel et al., 2020; Thomson et al., 2021; Ujiie et al., 2023; van Zyl et al., 2020; van Zyl et al., 2022). Comparisons to both human and mouse TM cell data in the Van Zyl et al. dataset were performed using the Broad Institute single-cell browser: https://singlecell.broadinstitute.org/single_cell/study/SCP780.

ATAC

To refine ATAC-seq analysis, we applied ArchR (Granja et al., 2021) to remove any additional doublet cells and iteratively estimated LSI embedding using 500bp peaks tiled across the genome, followed by ArchR clustering and UMAP projection workflows to generate ATAC-seq cell clusters and two-dimensional embeddings (same cluster resolutions as above). Using marker gene expression, clusters were named and related to the single cell cluster identities. Using single-cell RNA-seq data collected earlier with annotated cell types, we projected single-cell ATAC-seq data onto the reference using ArchR’s addGeneIntegrationMatrix function with default parameters. The integrated data with RNA-seq defined clusters were used for peak calling with MACS2, followed by marker peak identification (getMarkerFeatures) and TF enrichment (addMotifAnnotations) using the CisBP database (Weirauch et al., 2014). Finally, TF activity was quantified using chromVAR (Schep et al., 2017) and gene-peak links were inferred using the addPeak2GeneLinks function with default parameters. All significant TFs (VAR > 1) for each TM cell subtype were filtered by average RNA expression within the entire cluster using the Columbia University scRNA-seq data (better read depth than snRNA-seq). All significant TFs with no RNA expression were removed and all TFs with expression > 0.5 counts per 10,000 were prioritized for analysis. For a schematic of the analysis workflow, see Figure 4A. We then performed integrative analysis correlating the activity of each TF with its gene expression profile (both directly measured from multiome data and integrated expression data from reference scRNA-sequencing data). We used ArchR’s correlateMatrices function, which reports Pearson correlation of normalized counts (TF activity matrix and gene expression / gene integration matrix) across 100 randomly subsampled neighborhoods of cells. We designated influential TFs as those with correlation > 0.5, adjusted p-value of correlation < 0.01, and chromVAR deviation above 75% of all TFs.

ATAC tracts: We chose marker genes (TM1: Fn1, TM2: Gpc3, TM3: Acta2) and plotted accessibility tracks. Motif footprinting analysis was performed using ArchR getFootprints function using CisBP-defined motif positions for RNA-defined cell types with default parameters.

Immunoflourescence and in situ hybridization

Sections with antibody staining

Enucleated adult eyes were fixed for 1 hour at 4°C in 4% paraformaldehyde (PFA, Electron Microscopy Science, Hatfield, PA) prepared in phosphate-buffered saline (1X PBS, 137 mM NaCl, 10 mM phosphate, 2.7 mM KCl, pH 7.4). Following fixation, a small window was made through the optic nerve head so that eyes would not shrivel during dehydration. Eyes were dehydrated in a 30% sucrose (Sigma-Aldrich, St. Louis, Missouri) in 1X PBS at 4°C until eyes sunk to the bottom of a 2 ml glass vial. Once sunken, eyes were embedded in optimal cutting temperature embedding medium (Fisher Scientific, Waltham, MA) and flash frozen on dry ice. Eyes were cryosectioned at 10 μm thickness in the sagittal plane. Sections were evenly spaced throughout the peripheral and central ocular regions.

Sections were washed three times with 1× PBS with 0.3% Triton X-100 for 5 min to quench residual fixation. The sections were then incubated with 10% Donkey serum (Sigma-Aldrich, St. Louis, Missouri) and 0.3% Triton X-100 in 1× PBS (blocking buffer) at room temperature for 1 h to block nonspecific binding of antibody and to permeabilize the tissue. The sections were then incubated with primary antibodies (see Table S7) in 200 µl blocking buffer overnight at 4°C. The sections were then washed three times for 5 m with 1× PBS with 0.3% Triton X-100. Primary antibodies were detected with the appropriate species-specific secondary antibody (all Alexa 488, 594, or 647 at 1∶1000 dilution, Life Technologies, Grand Island, NY) diluted in 1× PBS with 0.3% Triton X-100 for 2 h at room temperature. The secondary detection solution also had DAPI at 1:1000 (Thermo Scientific, Waltham, MA) to label nuclei. The immunostained sections were washed three times for 5 m in 1× PBS with 0.3% Triton X-100. Sections were then cover-slipped using Flouromount (Sigma, St. Louis, MO).

Tyramide signal amplification (TSA)

We performed a TSA reaction for a subtype of antibodies (See Table S7). TSA was performed on slides containing serial transverse cryosections of enucleated mouse eyes using the Akoya Biosciences TSA Plus Cyanine 3 detection kit (NEL744001KT). Freshly cut sections were allowed to come to room temperature and were then quenched for 10 minutes in a 1% hydrogen peroxide, 10% methanol in 0.01M phosphate buffered saline solution (1X PBS). The slides were rinsed in 1X PBS and allowed to block for 1 hour at room temperature in a blocking solution of 10% donkey serum, 0.3% Triton-X 100 in 1X PBS. Primaries for TSA were diluted 1:100 in blocking solution and approximately 200 microliters were added to the slides overnight at 4 degrees Celsius. The slides were then rinsed using a 0.3% Triton-X 100 in 1X PBS solution before being blocked in Tris-NaCl-blocking (TNB) buffer for 1 hour at room temperature. Approximately 200 microliters of Horseradish peroxidase conjugated secondaries diluted 1:200 in TNB buffer were applied to the slides for 2 hours at room temperature. The slides were subsequently rinsed in Tris-NaCl-Tween (TNT) buffer. TSA Plus Working Solution was made by diluting the TSA Plus Stock Solution 1:50 in 1X Amplification Diluent. The TSA Plus Working Solution was added at a volume of 200 microliters for 7 minutes at room temperature. The slides were then washed in TNT buffer and approximately 200 microliters of non-amplified primaries for counterstaining diluted 1:200 in 10% donkey serum, 0.3% Triton-X 100 in 1X PBS blocking solution was applied overnight at 4 degrees Celsius. Secondary antibody detection and mounting were performed as above.

in situ hybridization

We used the RNAscope® Multiplex Fluorescent Reagent Kit v2-Mm (Advanced Cell Diagnostics a Bio Techne Brand, Newark, CA). Probes for mouse Lypd1 (#318361-C2), Inmt (#486371-C2), Edn3 (#505841-C1), Acta2 (α-SMA) (#319531-C3), and Myoc (#460981-C1) (also purchased from ACD Bio) were used. Fluorescent dyes Opal 690 (FP1497001KT) and Opal 620 (FP1495001KT) (Akoya Biosciences, Marlborough, MA) were used with each of these probes. Briefly, eyes from 3-month-old C57BL/6J were enucleated and fixed in 4% PFA for 24 h at 4°C. A window was made to the back of the eye cup, lens removed, and the anterior cup placed in 30% sucrose overnight at 4°C. Tissue was cryopreserved in OCT and placed at −80°C. Sections were cut at 12 μm thickness using SuperFrost Plus slides (Fisher Scientific, Hampton, NH).. RNAscope® was performed over two days as per manufacturer’s protocol.

Whole mounts

Enucleated eyes were fixed for 2 hours at 4°C in 4% paraformaldehyde (PFA, Electron Microscopy Science, Hatfield, PA) prepared in phosphate-buffered saline (1× PBS, 137 mM NaCl, 10 mM phosphate, 2.7 mM KCl, pH 7.4). The postnatal eyes were dissected out following fixation. The anterior part of the eye was cut just posterior to the limbus (transition zone between cornea and sclera) as described in Kizhatil et al. Briefly, the iris, lens, ciliary body, and thin strip of retina were carefully removed to obtain the anterior eye cup. The anterior cup includes the cornea, limbus, and a small portion of retinal pigmented epithelium. Four centripetal cuts were made to relax the eye-cup and facilitate eventual mounting onto a slide after immunofluorescence.

Anterior eye cups were rinsed multiple times with 1× PBS. The eye cups were incubated with 3% bovine serum albumin and 1% Triton X-100 in 1× PBS (blocking buffer) at room temperature for 1 h in 2 ml glass vials to block nonspecific binding of antibody and to permeabilize the tissue. The anterior cups were then incubated with primary antibodies (see Table S7) in 200 µl blocking buffer for 2 days, with rocking, at 4°C. The anterior cups were then washed three times over a 3-h period with 1× PBS. The primary antibodies were detected with the appropriate species-specific secondary antibody (all Alexa 488, 594, or 647 at 1∶1,000 dilution, Life Technologies, Grand Island, NY) diluted in blocking buffer, which also had DAPI to label nuclei. The immunostained eye cups were washed three times over a 3-h period in 1× PBS. Eye cups were then whole-mounted on slides in Fluoromount (Sigma, St. Louis, MO).

Microscopy of drainage structures

Sections: Microscopy of sectioned postnatal eye slides was performed using an LSM SP8 confocal microscope (Leica) and 40×1.1 NA water immersion objective.

In situ: Images were taken using a Nikon Eclipse 90i confocal laser-scanning microscope (Melville). Z-stack images of 1 µm thickness, using the 40X objective lens, were taken of each section and converted into maximum projection images that were used for analysis.

Whole mounts: Microscopy was performed using an LSM SP8 confocal microscope (Leica) using a 63×1.4 NA glycerol immersion objective. The Mark and Find mode was used to automate collection of images, generating a folder full of Z stacks (Z = 0.3 µm) at various individual overlapping positions around the limbus. We examined n = 6 eyes total. We imaged a certain distance around the ocular circumference, with an equal representation from each ocular quadrant (1.5-2.5 mm total).

Postprocessing of whole-mount images

3D rendering: We followed the guidelines detailed in Kizhatil et al for visualizing SC and TM in 3D. Individual confocal Z stacks (.lsm files) were processed directly using Imaris 9.5 (Oxford Instruments, Carteret, NJ). The maximum intensity projection setting in Surpass mode of Imaris was used for 3D rendering. Controls and experimental image sets were treated identically. Images were oriented so that structures of interest were visible. For whole-mounted eyes, multiposition Z stacks generated using the Mark and Find setting on the confocal microscope were first imported into Imaris and converted into Imaris files (.ims files). The resulting Imaris files were then stitched to generate a comprehensive Z stack encompassing the limbus using ImarisStitcher 9.5 (BitPlane AG, Zurich, Switzerland). The stitch was exported as an Imaris file and processed appropriately. For sections, images were taken in the Z-plane in focus. When making figures, the snapshot feature of Imaris was used to collect images at high resolution (1024×1024 pixels, 300 dpi).

TM marker location analysis

We defined the TM as the region between the inner wall of Schlemm’s canal (SC, outer TM) and the anterior chamber (inner TM), extending from the anterior edge of the pars plana (posterior TM) to the posterior edge of Descemet’s membrane of the corneal endothelium (anterior TM). The SC was defined using a panel of validated markers that varied to allow for different antibody species and different secondary antibodies to be used in conjunction with the SC markers (See Table S7 for list of antibodies, Kizhatil et al., 2014). Sections were only included in the scoring analysis if they met the following criteria:

1. Absence of obvious sectioning artifacts or abnormalities in the SC or TM region.

2. Distinguishability of all landmarks defining the TM.

3. Acceptably low background signal and absence of staining for non-specific or non-biological entities within the section.

For each section, the TM was segmented and analyzed in two separate ways. The anterior-posterior distance was measured, and the TM was divided in half along this axis at the midway point (Figure S7C-D). Using the surface feature in Imaris, the area of individual marker staining was measured, and the total area in both the anterior and posterior zones was calculated. Because markers are expressed at different levels, we normalized across markers by calculating the percent of individual marker expression in the anterior and posterior halves within each section. These percentage values in anterior and posterior halves were compared both within markers and across all markers of a given TM cell subtype by Student’s t-test. The same analysis was repeated for the inner-outer halves by measuring the inner-outer TM distance at the midpoint of the TM (anterior-posterior axis) and dividing the TM in half at this value. A total of 130-160 sections were examined for markers of each individual TM cell subtype (see Table S8 for breakdown across markers). These sections were from 10-12 individual eyes per cluster. All eyes were from WT mice of C57BL/6J background and between 3-6 months of age.

Additional analysis was performed by selecting a subset of examined sections with staining for select markers (TM1: MYOC antibody, TM2: CRYM antibody, Inmt in situ, TM3: α-SMA antibody). These sections were of exceptionally high quality and had: 1) Anterior-posterior TM length between 125-200 µm. 2) Inner-outer TM distance is between 15-30 µm.

With these selected sections, we subdivided the TM into eight different zones along both inner/outer and anterior/posterior axes (Figure S7F). For each individual section, we used Imaris to measure the anterior-posterior distance and divide the TM into three equal sectors (anterior, central, and posterior. Within each sector, the average depth of the TM was measured in Imaris. This calculation was used to divide the central and posterior sectors into three zones each (inner, intermediate, and outer). Because the anterior TM depth is thinner, we divided the anterior sector into two zones (inner and outer). In total, there were eight distinct zones.

The percent of total marker area was calculated for each of the eight zones for each individual section as described above. Percent of marker area in each zone was compared as shown in Figure 2 by one-way ANOVA followed by Tukey’s honestly significant difference test.

Scleral cell marker analysis

For CD34 and LY6C1 expression, we examined greater than 60 sections from 4 eyes (CD34) and 3 eyes (LY6C1) respectively. The observed locations of these markers within the iridocorneal angle are described.

Whole mount marker location analysis

To analyze the expression of markers within the TM, we segmented the TM based on anatomy. Because various anatomical features including the pars plana were removed during the whole mount process, we used the SC (marked by PECAM) as the major anterior-posterior landmark. The TM was defined as the group of cells inner to the SC and not outside the anterior or posterior boundary of the SC.

Nuclear morphology

Anterior segment whole-mounts were used to examine nuclear morphology in 3D (TM1: MYOC, TM3: α-SMA, all antibodies). The TM was segmented as described above and the ‘Surface’ feature in Imaris used to demarcate TM cell nuclei. Nuclei that were analyzed were immediately adjacent to the immunolabeled TM marker and had to be well spaced from other nuclei. To ensure accuracy we used manual inspection to exclude: 1) Nuclei that were very close to other nuclei or groups of nuclei and 2) nuclei of cells that were labeled with a Schlemm’s canal endothelial cell marker (PECAM) and/ or that projected into the lumen of SC. The ‘Sphericity’ tool was used to determine the degree of sphericity in Imaris. For each cell subtype, we assessed a total 60 nuclei from 6 eyes, assessing the first 10 randomly selected nuclei that met our inclusion criteria in each eye. Groups were compared by Student’s t-test.

Availability of data and materials

All data generated or analyzed during this study are included in this published article (and its supplementary information files). Raw sequencing data is available on GEO and SCP.

Acknowledgements

The authors would like to thank the members the Simon John Laboratory for experimental and technical assistance, Drs. Krish Kizhatil and Chi Zhang for reading and editing the manuscript, the Institute of Comparative Medicine and Devanshi Ragha at Columbia University for animal and veterinary care, the JP Sulzberger Columbia Genome Center and Duke University Molecular Genomics Core for assistance with genomics studies, and the P30EY019007 Columbia Shared Resources grant for histological experiments. Dr. Simon John is the Robert L. Burch III Professor of Ophthalmic Sciences at Columbia University’s Vagelos College of Physicians and Surgeons.

Additional information

Funding

This project was supported by National Eye Institute grant EY032507 (to S.W.M.J). Partial support was provided by EY11721, EY032062, and EY018606 (to SWMJ), Brightfocus grants CG2020004 (to S.W.M.J.), and G2021007 (to R.B.:), startup funds at Columbia University including the Precision Medicine Initiative, and the New York Fund for Innovation in Research and Scientific Talent (NYFIRST; EMPIRE CU19-2660). At Duke, the work was also supported by NEI grants EY028608, EY022359 (to WDS). Supported also came from Vision Core grants P30EY019007 (Columbia University), P30EY005722 (Duke University) and an unrestricted departmental award from Research to Prevent Blindness (Columbia University and Duke University). The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health. SWMJ is the Robert Burch III Professor of Ophthalmic Sciences.

Authors’ contributions

NT – conceived and designed and performed experiments, analyzed data, wrote the manuscript

RB – conceived and designed and performed experiments, analyzed data, edited manuscript

TL – conceived and designed analyses, analyzed data, edited manuscript

VBC –designed analyses, performed experiments

RK –designed analyses, performed experiments

GL –designed analyses, performed experiments

SM – performed experiments

MS – assisted data anal edited manuscript

CM – edited manuscript assisted w lab management and mouse colonies /physiology

DS – conceived and designed experiments, analyzed data, edited manuscript, sought and obtained funding

JQ – conceived and designed analyses, edited manuscript

SJ – conceived and oversaw overall study and designed experiments, analyzed data, wrote and edited manuscript, sought and obtained funding