Cellular blueprint of healthy and diseased human epiglottis and subglottis-a study of the Canadian Airways Research (CARE) group

Abstract

Background: The larynx consists of the supraglottis, glottis, and subglottis and each differ in tissue composition, lymphatic drainage, ability to counter infections, and response to injuries. However, the cellular mechanisms driving laryngeal homoeostasis remain largely unexplored. As a result, understanding disease pathogenesis within the larynx including idiopathic subglottic stenosis (iSGS) and intubation-related traumatic stenosis has been challenging. Here, we sought to characterise the cellular processes governing laryngeal health and disease.

Methods: As part of the prospective Canadian Airways Research (CARE) iSGS study, we characterised 122,004 high-quality transcriptomes using single nucleus RNA-sequencing to profile 11 human epiglottis and 17 human subglottis biopsies across three different conditions: control, iSGS, and intubation-related traumatic stenosis to define cell populations and pathways associated with disease. We validated our results using cohort-level bulk transcriptomics using 114 human epiglottis and 121 human subglottis.

Findings: We defined the single-cell taxonomy of the human subglottis and epiglottis using single-nucleus sequencing in both healthy and disease states. Mechanistically, we discovered the presence of unique epithelial and fibroblast progenitor subsets within the control subglottis but not within the anatomically adjacent epiglottis. The uncontrolled proliferation of these cellular subsets exhibited skewed sex hormone signalling and orchestrated a fibro-inflammatory cascade. We leveraged cohort-level bulk transcriptomics to define hallmarks of iSGS associated with disease covariates and introduced the first biomarker associated with recurrent relapse. Longitudinal sampling demonstrated that the subglottic microenvironment in patients with iSGS is changing dynamically with and without therapeutic intervention.

Interpretation: Together, our data refines our understanding of laryngeal biology, nominates candidate compounds for iSGS treatment, and serves as a transformative platform for future clinical investigations to further precision laryngology.

Funding: This study was funded by a grant from the American Laryngology Association (#1082), an Academic Medical Organisation of Southwestern Ontario innovation fund grant (INN21-016), grant support from the Departments of Otolaryngology-Head and Neck Surgery at University of Toronto, Canada and Western University, Canada. ACN was supported by the Wolfe Surgical Research Professorship in the Biology of Head and Neck Cancers Fund. PYFZ was supported by a Vanier Canada Graduate Scholarship and PSI Foundation fellowship.

Keywords: Biomarkers; Genomics; Idiopathic subglottic stenosis; Larynx; Single-nucleus RNA-sequencing; Transcriptomics; Traumatic subglottic stenosis.

Fig. 1

A single cell atlas of the healthy and diseased epiglottis and subglottis. a) UMAP visualisation of cell states and annotated with disease status, tissue type, sample type, coarse cell type (level 1) labels, and finely annotated (level 4) cell type labels. b) Distribution of coarse cell type annotations by patient (cell type level 1). c) Proportion of profiled cells by sample type (disease status and tissue type).

Fig. 2

Epithelial heterogeneity and dysfunction in subglottic stenoses. a) UMAP visualisation of 25 epithelial cell states and annotated with disease status, tissue type, sample type, level 3 cell type labels, and finely annotated (level 4) cell type labels. b) Sample type tissue prevalence estimated by R_o/e score. FDR derived from statistical compositional difference testing using sccomp package, corrected with Benjamini-Hochberg procedure. Statistical comparisons of iSGS epiglottis, trauma SGS epiglottis, and control subglottis made against control epiglottis, while iSGS subglottis and trauma SGS subglottis were compared against control subglottis. c) Abundance of each estimated of epithelial cell type in the epiglottis (top) and subglottis (bottom) from bulk RNA-seq using CIBERSORT. FDR from Benjamini-Hochberg adjusted Wilcoxon test. ∗: FDR <0.05, ∗∗: FDR <0.01, ∗∗∗: FDR <0.001, ∗∗∗∗: FDR <0.0001. d) Gene ontology enrichment value for androgen response in each epithelial subset. e) Pathway output activity of androgen and oestrogen inferred using PROGENy, normalised between groups. f) Inferred developmental trajectory of epithelial cells from the subglottis using pseudotime using scFate package, colour by their cell type level 3 annotations. g and h) PCA-representation of the developmental trajectory of epithelial cells from the subglottis inferred using RNA velocity and coloured using either the disease status (g) or cell type level 3 annotation (h). i) PCA visualisation of transcript abundance of CD55, MUC5AC, SCGB3A1, MUC5B, ERBB4, and SPDEF, along with PROGENy inferred pathway output activity score of androgen, oestrogen, and WNT-signalling. j–l) Scatter plot showing the specificity scores of regulons in E19 (j), E20 (k), and E22 (l) cell populations. The top 10 regulons are highlighted.

Fig. 3

Distinctive fibroblast populations within the subglottis primed for expansion. a) UMAP visualisation of eight fibroblast cell states. b) Sample type tissue prevalence estimated by R_o/e score. FDR from sccomp scores, corrected with Benjamini-Hochberg procedure. Statistical comparisons of iSGS epiglottis, trauma SGS epiglottis, and control subglottis made against control epiglottis, while iSGS subglottis and trauma SGS subglottis were compared against control subglottis. c) Abundance of each estimated fibroblast cell states in the epiglottis (top) and subglottis (bottom). FDR from Benjamini-Hochberg adjusted Wilcoxon test. ∗∗: FDR <0.01, ∗∗∗: FDR <0.001, ∗∗∗∗: FDR <0.0001. d) Inferred developmental trajectory of RNA velocity of fibroblast found in the control subglottis, trauma SGS subglottis, and iSGS subglottis. Pathways are coloured by activity of TGFβ, VEGF, and myogenesis pathways. e) PCA visualisation of transcript abundance of SLIT2, ABCA10, CREB5, AXIN2, SOX6, and COL1A1. f, g) Scatter plot showing the regulon specificity scores of transcription factor in F6 (F) and F4 (g) fibroblasts. Top 10 regulons are highlighted. h) Connectivity Map analysis of potential therapeutics to reverse gene expression of F4 fibroblasts. Norm connectivity score (CS) and adjusted P-values from ConnectivityMap query module. i) UMAP visualisation of immune cell subset. j) Sample type tissue prevalence estimated by R_o/e score. FDR from compositional analysis performed using the sccomp package, corrected with Benjamini-Hochberg procedure. Statistical comparisons of iSGS epiglottis, trauma SGS epiglottis, and control subglottis made against control epiglottis, while iSGS subglottis and trauma SGS subglottis were compared against control subglottis. k) Abundance of each estimated immune cell states in the epiglottis (top) and subglottis (bottom). FDR from Benjamini-Hochberg adjusted Wilcoxon test. ∗: FDR <0.05.

Fig. 4

Cellular and molecular classification of iSGS. a) Cell type modules identified using consensus clustering of cell type abundance (n = 83). Only cell types with R_o/e above 0.01 in at least one of control, iSGS, or trauma SGS subglottis is included. b) SUBglottiS gEne subtype (SUBSET) microenvironment archetypes identified using clustering of identified relative abundance. c) Association between each SUBSET cluster and total number of interventions (n = 72). Only the first sample for patients with multiple samples were considered. Estimate, 95% confidence interval, and P-values from Poisson regression. d) Summary of patients stratified by SUBSET. e) SUBSET changes with or without treatment demonstrating microenvironment plasticity patients with iSGS with multiple samples.