Identification of Meibomian gland stem cell populations and mechanisms of aging

Abstract

Meibomian glands secrete lipid-rich meibum, which prevents tear evaporation. Aging-related Meibomian gland shrinkage may result in part from stem cell exhaustion and is associated with evaporative dry eye disease, a common condition lacking effective treatment. The identities and niche of Meibomian gland stem cells and the signals controlling their activity are poorly defined. Using snRNA-seq, in vivo lineage tracing, ex vivo live imaging, and genetic studies in mice, we identify markers for stem cell populations that maintain distinct regions of the gland and uncover Hedgehog (Hh) signaling as a key regulator of stem cell proliferation. Consistent with this, we show that human Meibomian gland carcinoma exhibits increased Hh signaling. Aged glands display decreased Hh and EGF signaling, deficient innervation, and loss of collagen I in niche fibroblasts, indicating that alterations in both glandular epithelial cells and their surrounding microenvironment contribute to age-related degeneration. These findings suggest new approaches to treat aging-associated Meibomian gland loss.

Fig. 1. snRNA-seq identifies tarsal plate and MG cell subpopulations.

a UMAP plot of snRNA-seq data from pooled tarsal plate samples from four 8-week-old and four 21-month-old male C57BL/6 J mice, respectively, with two replicates for each age. b Schematic diagram of MG structure indicating the spatial locations of MG subpopulations and the stromal environment. c, d RNA Velocity (c) and pseudotime analyses (d) predictions of the differentiation trajectories of MG cellular subpopulations. RNA Velocity (c) predicts that ductular cells can differentiate towards either duct or acinus. See also Supplementary Fig. S1.

Fig. 2. Identification of distinct MG stem cell populations.

a Scheme for lineage tracing. b Lrig1⁺ cells in MG ductal basal layer (yellow arrows) and acinar basal layer (white arrow) self-renew and generate progeny that contribute to MG duct (light blue arrows) and acinus (pink arrow). c Lgr6⁺ cells in MG duct (yellow arrows) and acinar basal layer (white arrow) self-renew and contribute to MG duct (light blue arrow) and acinus (pink arrow). d Axin2⁺ cells in the MG duct (yellow arrows) and acinar basal layer (white arrow) self-renew and replenish MG duct (light blue arrow) and acinus (pink arrow). e At 2 days, Gli2 marks acinar basal cells (white arrowheads), ductule cells (orange arrowhead), and meibocytes (white arrowheads); at 90 days, lineage-traced cells are present in MG duct (yellow and green blue arrowheads), ductule (orange arrowhead) and acinus (white and pink arrowheads). f At 2 days, Slc1a3 marks meibocytes (yellow arrowhead) and acinar basal cells (white arrowhead); at 90 days, lineage tracing labels the acinar basal layer (white arrowhead) and meibocytes (yellow arrowheads). Scale bars: (b–e), 50μm; (f), 25μm. White dashed lines in (b–f) outline MG acini and ducts. N = 3 mice (2 males and 1 female) of each genotype were analyzed in lineage tracing experiments for each line and stage. See also Supplementary Figs. S2–S4; Supplementary Movies S1 and S2.

Fig. 3. Hh signaling is required for normal levels of MG basal cell proliferation.

a–c Whole-mount fluorescence images of MG #6 and/or MG #7 in (a) KRT14-Cre^ERT2 Rosa26^mTmG control mice 18 or 29 weeks after tamoxifen treatment at 8 weeks of age (T18 Weeks and T29 Weeks) (n = 4, 2 male and 2 female); b, c KRT14-Cre^ERT2 Smo^fl/fl Rosa26^mTmG mice of T18 Weeks (b) or T29 Weeks (c) (n = 4, 3 male and 1 female). Yellow asterisks indicate hair follicles; CD, MG central duct. (d-g) IHC for PPARγ (d, e) and FASN (f, g) showing similar expression levels in control (d, f) and Smo-deficient (e, g) MGs at 34 weeks after tamoxifen treatment. h–k IHC for Ki-67 showing reduced acinar and ductal basal cell proliferation in Smo-deficient MGs at 34 weeks after tamoxifen treatment (i, k) compared with control (h, j, green arrowheads). l, m Quantitation of the % of Ki-67⁺ cells in the acinar (l) and ductal (m) basal layer of control and Smo-deficient MGs. Littermate pairs were compared. Statistical significance was calculated using a paired two-tailed Student’s t-test. n = 4 KRT14-Cre^ERT2 control mice (3 male and 1 female) and n = 4 Smo-deficient (KRT14-Cre^ERT2 Smo^fl/fl) mice (3 male and 1 female) were analyzed in (d–m). At least 70 acinar basal cells and 80 ductal basal cells were analyzed per mouse in (l, m). Source Data for (l, m) are provided as a Source Data file. Scale bars: 50 μm. See also Supplementary Fig. S5.

Fig. 4. Increased Hh signaling promotes proliferation and expansion of MG basal cells.

a–f Progressive expansion of acinar basal cell clusters (black arrows) in GLI2ΔN-expressing MG after 2 (a, b), 4 (c, d), and 7 (e, f) days of doxycycline treatment at 8 weeks of age. Meibocytes were reduced by 7 days (f, green arrows). g–l Progressive expansion of GLI2ΔN⁺/KRT5⁺ cells (yellow arrows) and reduction of PLIN2⁺ meibocytes (white arrows) in Gli2ΔN^Krt5rtTA MG acini after 2 (g, h), 4 (i, j) and 7 (k, l) days of doxycycline treatment. N = 3 control mice lacking Krt5-rtTA or tetO-GLI2ΔN (2 males and 1 female) and n = 3 Gli2ΔN^Krt5-rtTA mice (2 males and 1 female) were analyzed at each time point in (a–l); representative samples are shown. m–p Hyperproliferation of GLI2ΔN⁺ cells in Gli2ΔN^Krt5-rtTA MG acini (m, n, yellow arrow) and ducts (o, p, white arrow). q, r Quantification of acinar basal cell (q) and ductal basal cell (r) proliferation. Samples from n = 5 control mice (3 males and 2 females) lacking Krt5-rtTA or tetO-GLI2ΔN and samples from n = 5 Gli2ΔN^Krt5-rtTA mice (3 males and 2 females) were analyzed in (m–r). Statistical significance in (q, r) was calculated with unpaired two-tailed Student’s t-test. Data are presented as mean +/− SEM. At least 100 acinar cells and 110 ductal cells were analyzed from each animal. Source Data for (q, r) are provided as a Source Data file. s, t Overgrowth of Ptch1-deficient acinar basal cells (yellow arrows) in KRT14-Cre^ERT2 Ptch1^fl/fl mice 12 weeks after tamoxifen treatment at 8 weeks of age. u, v Expansion of KRT5⁺ cells in Ptch1-deficient acini (red arrows). Independent biological samples from n = 3 male KRT14-Cre^ERT2 Ptch1^fl/fl mice, n = 2 age-matched male controls and n = 2 age-matched female controls were analyzed in (s–v). Representative images are shown. Scale bars: 50 μm.

Fig. 5. Forced GLI2ΔN expression expands MG stem cell populations and impedes meibocyte differentiation.

a Scheme for bulk RNA-seq of laser-captured MG samples. b GO enrichment analysis of bulk RNA-seq data from GLI2ΔN^Krt5rtTA and littermate control MGs 4 days after induction at 8 weeks of age showing the top 10 enriched pathways. c Volcano plot showing genes upregulated or downregulated at 4 days. Independent biological samples from n = 3 GLI2ΔN^Krt5rtTA mice (2 male and 1 female) and independent biological samples from n = 3 littermate control mice lacking Krt5-rtTA or tetO-GLI2ΔN (2 males and 1 female) were analyzed; differentially expressed genes were defined as padj<0.001 and Log₂FC < -0.5 or Log₂FC > 0.5. FDR calculation was performed by DESeq2 v1.20.0 with the Benjamini-Hochberg procedure. d–g RNAscope showing upregulation of Gli1 and Ccnd1 in GLI2ΔN^Krt5rtTA MGs. h, i IF data showing decreased PPARγ and PLIN2 expression in GLI2ΔN^Krt5rtTA acini. (j-m) RNAscope showing expanded Lrig1 and Lgr6 expression in GLI2ΔN^Krt5rtTA MGs. Independent samples from 3 Gli2ΔN^Krt5-rtTA mice (2 males and 1 female) and 3 littermate controls of genotypes tetO-GLI2ΔN or Krt5-rtTA (2 males and 1 female) were used for RNAscope and IF; all mice were doxycycline-treated for 4 days. (n-s) GLI2ΔN^Krt5-rtTA Rosa26^mTmG mice carrying inducible Cre alleles driven by Lrig1 (n, o), Lgr6 (p, q) or Axin2 (r, s) promoters were tamoxifen-induced at P42 to induce Cre activity and placed on oral doxycycline at P72 to induce GLI2ΔN expression. mGFP expression (red signal) and GLI2 expression (green signal) were analyzed by IF at P74 (n, p, r) or P82 (o, q, s). GLI2ΔN-expressing cells in the acinar basal layer positive for Lrig1 (n), Lgr6 (p), or Axin2 (r) (yellow arrows) give rise to clones that contribute to MG overgrowth (o, q, s, white arrows). n = 3 (1 male and 2 females) samples were analyzed per line per time point. Scale bars: (d–m), 50 μm; (n–s), 25 μm.

Fig. 6. GLI2-, LRIG1- and LGR6-expressing undifferentiated cells are expanded in human MGC.

a, b In normal human MG GLI2 protein localizes to acinar basal cells (a, yellow arrows) and differentiating meibocytes (a, white arrows) but is absent from fully differentiated meibocytes (a, light blue arrows); in MGC samples GLI2 is broadly expressed (b). c, d In normal human MG, GLI1 mRNA is weakly expressed in some acinar basal cells (c, yellow arrows); GLI1⁺ cells are widely present in MGCs (d). e–h LRIG1⁺ and LGR6⁺ cells are primarily present in the acinar basal layer of normal human MGs (e, g, yellow arrows); MGC displays expansion of LRIG1⁺ and LGR6⁺ cells (f, h). i, j Ki-67⁺ cells localize to the normal human acinar basal layer (i, white arrows) and are expanded in MGC (j, yellow arrows). k, l PLIN2-/KRT14⁺ cells are restricted to the acinar basal layer in normal human MG (k, white arrows) but are distributed more broadly in MGC tissue (l, yellow arrows). 7 normal human MG and 10 human MGC samples were analyzed; representative data are shown. Scale bars: 50 μm.

Fig. 7. Aged MGs exhibit fewer cells with Hh activity and elevated association of GLI2 with acetyl-lysine.

a Violin plots derived from snRNA-seq data show relatively reduced percentages of cells expressing Gli2 and Ptch1 in the ductule and acinar basal layer, respectively, in aged compared to young MGs. b, c Ptch1 is expressed in acinar basal cells (b, white arrow), some meibocytes (b, yellow arrow), and surrounding stromal cells (b, blue arrow) at 8 weeks; fewer Ptch1-expressing cells are present in aged MG (c). (d, e) Fewer Gli2-expressing cells are present in aged MG particularly in the ductule (d, e, white arrows). f, g RNAscope using a pan-Hh probe which detects both Ihh and Shh mRNAs reveals similar expression of mRNA for Hh ligands in the acini of young (f) and aged (g) MGs. h, i PLA for GLI2 and acetyl-lysine shows that close association (<40 nm) of GLI2 with acetyl-lysine is elevated in subsets of cells within aged (i) compared with young (h) MG acini (i, white arrows) and stroma (i, blue arrows). 3 male C57BL/6 J mice were analyzed at each age in (b-i); representative data are shown. j, k Association of GLI2 and acetyl-lysine is elevated in the acinar basal layer (k, white arrows) and meibocytes (k, yellow arrows) but not stromal cells (j, k, blue arrows) in Krt5-rtTA tetO-Cre Hdac1^fl/fl Hdac2^fl/fl (Hdac1/2^DcKO) mice doxycycline induced from P48 and analyzed at P58 compared to control mice of genotype Hdac1^fl/fl Hdac2^fl/fl and lacking Krt5-rtTA or tetO-Cre. (j). l, m Ptch1 expression in control acinar basal cells (l, white arrows) and differentiating acinar cells (l, yellow arrow) is reduced in HDAC1/2-deficient MG acini (m); stromal Ptch1 expression is similar in mutants and controls (l, m, blue arrows). 3 male mice of each genotype were analyzed in (j-m); representative data are shown. n Violin plots of snRNA-seq data indicate that there is a relatively lower percentage of Hdac1-expressing acinar basal cells (cluster #16) but a higher percentage of Hdac2-expressing acinar basal cells and ductular cells (cluster #17) in aged compared to young MG. o–r IF for HDAC1 (o, p) and HDAC2 (q, r) shows similar expression in young (o, q) and aged (p, r) MGs. 3 male C57BL/6 J mice of each age were analyzed in (o–r); representative data are shown. White dashed lines in (b–m; o–r) outline MG acini and/or ducts. Scale bars: 25 μm. See also Supplementary Fig. S6.

Fig. 8. HBEGF signaling is disrupted in aged MGs.

a, b CellChat analysis predicts decreased EGF signaling in aged MGs. c CellChat analysis predicts HBEGF signaling between dermal cells and MG ductular cells (green) at 8 weeks and its absence at 21 months. d Violin plots of snRNA-seq data show decreased percentages of Hbegf-expressing cells in the indicated cell populations in aged MG. e RNAscope shows reduced Hbegf expression levels per cell and fewer acinar basal cells (yellow arrows) and surrounding dermal cells (white arrows) expressing Hbegf in aged MG. f IHC shows that p-ERK1/2 levels are decreased, but total ERK1/2 levels are similar in the acini of aged compared with young MG. White dashed lines in (e) and black dashed lines in (f) outline MG acini. Samples from n = 3 8-week and n = 3 21-month-old male C57BL/6 J mice were used for RNAscope and IHC. Scale bars represent 25 μm.

Fig. 9. Aged MGs exhibit reduced peripheral innervation and collagen I expression.

a GO analysis of genes with statistically significantly reduced expression in aged compared with young acinar basal cells (cluster #16). b Violin plots for neuronal guidance genes in acinar basal cells. c, d)Slit3 mRNA expression (yellow arrows) is reduced in aged MG acinar basal cells. e, f GO analysis of genes with reduced expression in aged cluster #23 dermal cells (e) and violin plots for neuronal guidance genes (f). g, h Reduced numbers of PGP9.5 nerve fibers (red) adjacent to the KRT5⁺ acinar basal layer (green) in aged mice. i Quantification of PGP9.5⁺ nerve fibers within 5 μM of each KRT5⁺ acinar basal cell in young and aged MGs. Independent biological samples from n = 4 8-week-old and n = 4 21-month-old male C57BL/6 J mice were used for quantification. At least 185 acinar basal cells were analyzed from each mouse. Significance in (i) was calculated by an unpaired two-tailed Student’s t-test. Data are represented as mean +/− SEM. Source Data are provided as a Source Data file. j, k GO analysis of genes with reduced expression in aged cluster #15 dermal fibroblasts (j) and violin plots for type I collagen genes in dermal fibroblasts (k). l, m Reduced expression of Col1a1 (red) in aged (m) versus young (l) dermis surrounding MG acini. For (a, e, j) FDR calculation was performed by clusterProfiler (v4.4.4) and EdgeR (v3.38.4) with the Benjamini-Hochberg procedure. Dot plots in (a, e, j) were generated with ggplot2. Violin plots in (b, f, k) were generated with the VlnPlot command from Seurat package. White dashed lines in (c, d, l, m) outline MG acini and/or ducts. Independent biological samples from n = 3 8-week and n = 3 21-month C57BL6/J mice were used for RNAscope, IF and IHC. Scale bars represent 25 μm. See also Supplementary Fig. S7.

Fig. 10. Model for the effects of aging in the MG.

Hh signaling is reduced in aging due to inhibited GLI2 transcriptional activity, resulting in lower levels of stem cell proliferation. In parallel, peripheral innervation and collagen fibril density are reduced in aged compared with young MGs.