A revised airway epithelial hierarchy includes CFTR-expressing ionocytes

Abstract

The airways of the lung are the primary sites of disease in asthma and cystic fibrosis. Here we study the cellular composition and hierarchy of the mouse tracheal epithelium by single-cell RNA-sequencing (scRNA-seq) and in vivo lineage tracing. We identify a rare cell type, the Foxi1⁺ pulmonary ionocyte; functional variations in club cells based on their location; a distinct cell type in high turnover squamous epithelial structures that we term 'hillocks'; and disease-relevant subsets of tuft and goblet cells. We developed 'pulse-seq', combining scRNA-seq and lineage tracing, to show that tuft, neuroendocrine and ionocyte cells are continually and directly replenished by basal progenitor cells. Ionocytes are the major source of transcripts of the cystic fibrosis transmembrane conductance regulator in both mouse (Cftr) and human (CFTR). Knockout of Foxi1 in mouse ionocytes causes loss of Cftr expression and disrupts airway fluid and mucus physiology, phenotypes that are characteristic of cystic fibrosis. By associating cell-type-specific expression programs with key disease genes, we establish a new cellular narrative for airways disease.

Extended Data Figure 1 |. Identifying tracheal epithelial cell types in 3’ scRNA-seq

a. Quality metrics for the initial droplet-based 3’ scRNA-seq data. Distributions (y axis) of the number of reads per cell (x-axis, left), the number of the genes detected with non-zero transcript counts per cell (x-axis, center), and the fraction of reads mapping to the mm10 transcriptome per cell (x-axis, right). Dashed and blue lines: median value and kernel density estimate, respectively. b. Cell type clusters are composed of cells from multiple biological replicates. Fraction of cells in each cluster that originate from a given biological replicate (color legend, bottom right, n=6 mice); post-hoc annotation and number of cells are indicated above each pie chart. All biological replicates contribute to all clusters (except for WT mouse 1 which did not contain any of the very rare ionocytes: 0.39% of all epithelial cells), and no significant batch effect was observed. c. Reproducibility between biological replicates. Average gene expression values (log₂(TPM+1), x and y axes) across all cells of two representative 3’ scRNA-seq replicate experiments (Pearson correlation coefficient, top left), blue shading: gene (point) density. d. Post-hoc cluster interpretation based on the expression of known cell type markers. tSNE of 7,193 scRNA-seq profiles (points), colored by cluster assignment (Methods, top left) or by the expression (log₂(TPM+1), color bar) of a single marker genes or the mean expression of several marker genes for a particular cell type. e. Cell type clusters. Pearson correlation coefficients (r, color bar) between every pair of 7,193 cells (rows and columns) ordered by cluster assignment (color bar). Inset (right): zoom of 288 cells from the rare types. f. Gene signatures. Relative expression level (row-wise Z-score of log₂(TPM+1) expression values, color bar) of cell type-specific genes (rows) in each epithelial cell (columns). Large clusters (basal, club) are down-sampled to 500 cells.

Extended Data Figure 2 |

301 scRNA-seq profiles (points) colored by region of origin (top left panel), cluster assignment (top second panel, Methods), or, for the remaining plots, the expression level (log₂(TPM+1), color bar) of a single marker genes or the mean expression of several marker genes for a particular cell type. All clusters are populated by cells from both proximal and distal epithelium except rare NE cells, which were only detected in proximal experiments (top left panel).

Extended Data Figure 3 |. High-confidence consensus cell type markers, and cell type-specific expression of asthma-associated genes

a. Cell type clusters in full-length plate-based scRNA-seq data. Cell-cell Pearson correlation coefficient (r, color bar), between all 301 cells (individual rows and columns) ordered by cluster assignment (color bar, as in Extended Data Fig. 2d). Right: zoomed in view of 17 cells (black border on left) from the rare types. b. High confidence consensus markers. Relative expression level (row-wise Z-score of mean log₂(TPM+1), color bar at bottom) of consensus marker genes (rows, FDR <0.01 in both 3’-droplet and full-length plate-based scRNA-seq datasets, likelihood-ratio test) for each cell type (flanking color bar) across 7,193 cells in the 3’ droplet data (columns, left) and the 301 cells in the plate-based dataset (columns, right). Top 15 markers shown, complete sets are in Extended Data Fig. 1f and Supplementary Table 3. c. Cluster-specific transcription factors (TFs) in 3’ scRNA-seq data. Mean relative expression (row-wise Z-score of mean log₂(TPM+1), color bar) of the top TFs (rows) that are enriched (FDR < 0.01, likelihood-ratio test, two-sided) in cells (columns) of each cluster. d-f. Cell type-specific expression of genes associated with asthma by GWAS. d. Relative expression (Z-score of mean log₂(TPM+1), color bar bottom right) of genes (rows) that are associated with asthma in GWAS and enriched (FDR < 0.01, likelihood-ratio test) for cell type (columns) specific expression in our 3’ scRNA-seq data. e. For each gene from (d) shown is the significance (-log₁₀(FDR), Fisher’s combined p-value, likelihood-ratio test, y axis) and effect size (point size, mean log₂(fold-change)) of cell type specific expression in the relevant cell (color legend) and its genetic association strength from GWAS (x axis). f. Distribution of expression levels (y axis, log₂(TPM+1)) in the cells in each cluster (x axis, color legend) for two asthma GWAS genes: Cdhr3 (left; specific to ciliated cells) and Rgs13 (right; specific to tuft cells). FDRs: LRT.

Extended Data Figure 4 |. Krt13⁺ progenitors express a unique set of markers distinct from mature club cells

a. Proximal vs. distal specific club cell expression. Relative expression level (row-wise Z-score, color bar) for genes (rows) enriched in proximal and distal tracheal club cells (FDR<0.05, likelihood-ratio test) in the full-length scRNA-seq data. b. Distal epithelia differentiate into mucous metaplasia. Goblet cell quantification (ln(Muc5ac⁺/ GFP⁺ ciliated cells, y-axis) in Foxj1-GFP mice (n=6, dots) in each of four conditions in (Fig. 2a) (x-axis). p values: Tukey’s HSD test, black bars: mean, error bars: 95% CI. c. Krt8 does not distinguish pseudostratified club cell development from hillock-associated club cell development. Diffusion map embedding of 6,905 cells (as in Fig. 2b) colored either by their Krt13⁺ hillock membership (left: green), or by expression (log₂(TPM+1), color bar) of specific genes (all other panels). d. Immunostaining of hillock strata. Left: Krt13⁺ (green) and Trp63⁺ (magenta) basal (solid outline) and suprabasal (dashed outline) cells. Right: Krt13⁺ (green) and Scgb1a1⁺ (magenta, solid outline) luminal cells, n=3 mice. e,f. Krt13⁺ hillock cells are highly proliferative. e. Co-stain of EdU (magenta) and Krt13 (green), n=4 mice. f. Fraction of EdU⁺ epithelial cells (%, y-axis) in hillock (mean: 7.7%, 95% CI [4.8%, 10.5%]) and non-hillock (mean: 2.4%, 95% CI [1.8%, 3.1%]) areas (x axis). p values: LRT, n=4 mice, black bar: mean, error bars: 95% CI. g. Fraction of Krt13⁺ hillock cells that are club cell lineage labeled (%, y axis) decreases from day 5 (10.2%, 95% CI [0.07, 0.16]) to day 80 (5.2%, 95% CI [0.03, 0.08]). Error bars: 95% confidence interval, n=3 mice (dots). p values: LRT. h. Differential expression (x axis, log₂(fold-change)) and associated significance (y axis, log₁₀(FDR)) for each gene (dot) that is differentially expressed in Krt13⁺ cells (identified using clustering in diffusion map space, Methods) compared to all cells (FDR<0.05, LRT). Color code: cell type with highest expression (green: genes whose highest expression is in Krt13⁺ hillock cells). Dots show all the genes differentially expressed (FDR<0.05) between Krt13⁺ hillock cells and other cells. Those genes with absolute effect sizes greater than log₂ fold-change > 1 are marked with large points, while others are identified as small points (grey). i. Enriched pathways in Krt13⁺ hillock cells. Representative MSigDB gene sets (rows) that are significantly enriched (x axis and color bar, -log₁₀(FDR), hypergeometric test) in Krt13⁺ hillock cells.

Extended Data Figure 5 |. Genes associated with cell fate transitions

Relative mean expression (loess-smoothed row-wise Z-score of mean log₂(TPM+1), color bar at bottom) of significantly (p < 0.001, permutation test) varying genes (a-d) and TFs (e-h) (rows) across subsets of 6,905 (columns) basal, club and ciliated cells. Cells are pseudotemporally ordered (x axis, all plots) using diffusion maps (Fig. 2b and Extended Data Fig. 4c). Each cell was assigned to a cell fate transition if it was within d < 0.1 of an edge of the convex hull of all points (where d is the Euclidean distance in diffusion-space) is assigned to that edge (Methods).

Extended Data Figure 6 |. Lineage tracing using Pulse-Seq

a. Pulse-Seq experimental design schematic. mT: membrane-tdTomato, mG: membrane-EGFP. b. Post-hoc cluster annotation by known cell type markers. tSNE of 66,265 scRNA-seq profiles (points) from Pulse-Seq, colored by the expression (log₂(TPM+1), color bar) of single marker genes for a particular cell type or cell-cycle score (bottom right) c. Pulse-Seq lineage-labeled fraction of various cell populations over time. Linear quantile regression fits (trendline, Methods) to the fraction of lineage-labeled cells of each type (n=3 mice per time point, dots, y-axis) as a function of the number of days post tamoxifen-induced labeling (x-axis). β: estimated regression coefficient, interpreted as daily rate of new lineage-labeled cells, p: p value for the significance of the relationship, Wald test (Methods). As expected, goblet and ciliated cells are labeled more slowly than club cells (Fig. 3d). d. Labeled fraction of basal cells is unchanged during Pulse-Seq time course, as expected. Estimated fraction (%, y-axis, Methods) of cells of each type that are positive for the fluorescent lineage label (by FACS) in each of n=3 mice (points) per time-point (x axis). p values: LRT, error bars: 95% CI. e. Proportion of basal cell lineage-labeled tuft cells at day 0 (0%, n=2 mice, dots) and day 30 (22.9%, 95% CI [0.17, 0.30], bars: estimated proportions, n=3 mice). Error bars: 95% CI, p values: LRT. f-h. Conventional Scgb1a1 (CC10) lineage trace of rare epithelial types shows minimal contribution to rare cell lineages. Fraction of Scb1a1 labeled (club cell trace) cells (y axis, %) of Gnat3⁺ tuft cells (f) at day 0 (n=3 mice, 0.6%, 95% CI [0.00, 0.04]) and day 30 (n=2 mice, 6.3%, 95% CI [0.04, 0.11]), Foxi1-GFP⁺ ionocytes at day 30 (n=2 mice, 2.9%, 95% CI [0.01, 0.11]) (g), and Chga⁺ neuroendocrine (NE) cells at day 0 (n=2 mice, 2.5%, 95% CI [0.01, 0.08]) and day 30 (n=2 mice, 2.6%, 95% CI [0.01, 0.08]) (h) after club cell lineage labeling. p values: LRT. Error bars: 95% confidence interval.

Extended Data Figure 7 |. Club cell heterogeneity and lineage tracing hillock-associated club cells using Pulse-Seq.

a,b. Principal components are associated with basal to club differentiation (PC-1), proximodistal heterogeneity (PC-2), and hillock gene modules (PC-2). a. PC-1 (x-axis) vs. PC-2 (y-axis) for a PCA of 17,700 scRNA-seq profiles of club cells (points) in the Pulse-Seq dataset, colored by signature scores (color legends, Methods) for basal (left), proximal club cells (center left), distal club cells (center right), the Krt13⁺/Krt4⁺ hillock (right), or their cluster assignment (inset, right). b. Bar plots show the extent (normalized enrichment score, y-axis, Methods) and significance of association of PC-1 (left) and PC-2 (right) for gene sets associated with different airway epithelial types (x-axis), or gene modules associated with proximodistal heterogeneity (Extended Data Fig. 4a). Heatmaps shows the relative expression level (row-wise Z-score of log₂(TPM+1) expression values, color bar) of the 20 genes (rows) with the highest and lowest loadings on PC-1 (left) and PC-2 (right) in each club cell (columns, down-sampled to 1,000 cells for visualization only). p values: permutation test (Methods). c. Pulse-Seq lineage tracing of hillock-associated cells. Estimated fraction (%, y-axis, Methods) of cells of each type that are positive for the fluorescent lineage label (by FACS) from n=3 mice (points) per time-point (x axis). p values: LRT. (Methods), error bars: 95% CI. d. Hillock-associated club cells are produced at a greater rate than all club cells. Estimated rate (%, y-axis) based on the slope of quantile regression fits (Methods) to the fraction of lineage-labeled cells of each type (x-axis). p values: rank test (Methods), error bars: 95% CI. e,f. Club cells initially labeled by Pulse-Seq are associated with basal to club cell differentiation. e. Distribution of basal signature scores (y axis) for individual club cells (points) from each Pulse-Seq time point and lineage label status (x axis). p value: Mann-Whitney U-test. Violin plots show the Gaussian kernel probability densities of the data, large white point shows the mean. mT: membrane-tdTomato, mG: membrane-EGFP. f. PC-1 (x-axis) vs. PC-2 (y-axis) for a PCA of 17,700 scRNA-seq profiles of club cells (points), as in (a), highlighting club cells that are lineage-labeled at the initial time point (legend). g. Schematic of the more rapid turnover of basal to club cells inside (top) and outside (bottom) hillocks.

Extended Data Figure 8 |. Heterogeneity of rare tracheal epithelial cell types

a. Cell type-enriched GPCRs. Relative expression (Z-score of mean log₂(TPM+1), color bar) of the GPCRs (columns) that are most enriched (FDR < 0.001, LRT) in the cells of each tracheal epithelial cell type (rows) based on full-length scRNA-seq data. b. Tuft cell-specific expression of Type I and Type II taste receptors. Expression level (mean log₂(TPM+1), color bar) of tuft-cell enriched (FDR<0.05, LRT) taste receptor genes (columns) in each tracheal epithelial cell type (rows, labeled as in a) based on full-length scRNA-seq data. c. Tuft cell-specific expression of the Type-2 immunity-associated alarmins Il25 and Tslp. Expression level (y-axis, log₂(TPM+1)), of Il-25 (left) and Tslp (right) in each cell type (x axis). FDR: LRT. Violin plots show the Gaussian kernel probability densities of the data. d. Morphological features of tuft cells. Immunofluorescence staining of the tuft-cell marker Gnat3 (yellow) along with DAPI (blue). Arrowhead: “tuft”, arrows: cytoplasmic extension. e,f. Tuft-1 and tuft-2 sub-clusters. e. tSNE visualization of 892 tuft cells (points) colored either by their cluster assignment (left, color legend), or by the expression level (log₂(TPM+1), color bar, remaining panels) of marker genes for mature tuft cells (Trpm5), tuft-1 (Gng13), tuft-2 (Alox5ap) subsets. f. Distribution of expression levels (y-axis, log₂(TPM+1)) of the top markers for each subset (x-axis). Violin plots show the Gaussian kernel probability densities of the data, large white point shows the mean. FDR: LRT, n=15 mice. g. Tuft-1 and tuft-2 subtypes are each generated from basal cell parents. Estimated fraction (%, y-axis, Methods) of cells of each type that are positive for the basal-cell lineage label (by FACS) from n=3 mice (points) per time-point (x-axis) in the Pulse-Seq experiment. p values: LRT, error bars: 95% CI. h. Differential expression of tuft cell associated transcription factors between tuft subtypes. Labeled genes are differently expressed in the tuft cell subsets (FDR < 0.01, likelihood-ratio test). i,j. Mature and immature subsets are identified using marker gene expression. The distribution of expression of scores (y-axis, using top 20 marker genes, Supplementary Table 1, Methods) for tuft (i) goblet (j), basal and club cells (label on top) in each cell subset (x axis) (basal and club cells down-sampled to 1,000 cells). p values: Mann-Whitney U-test. k,l. Gene signatures for goblet-1 and goblet-2 subsets. The distribution (k) and relative expression level (l, row-wise Z-score, color bar) of marker genes that distinguish (log₂ fold-change > 0.1, FDR< 0.001, likelihood-ratio test) cells in the goblet-1 and goblet-2 sub-clusters (color bar, top and left) from the combined 3’ scRNA-seq datasets. m. Immunofluorescence staining of the goblet-1 marker Tff2 (magenta), the known goblet cell marker Muc5ac (green), and DAPI (blue). Solid white line: boundary of a goblet-1 cell.

Extended Data Figure 9 |

1) mouse. EGFP appropriately marks Foxi1 antibody-positive cells (left panel, solid white line). EGFP⁺ cells express canonical airway markers Ttf1 (Nkx2–1) and Sox2 (solid white lines). EGFP(Foxi1)⁺ cells do not label with basal (Trp63), club (Scgb1a1), ciliated (Foxj1), tuft (Gnat3), neuroendocrine (NE) (Chga), or goblet (Tff2) cell markers (dashed white lines). b. Ionocytes are sparsely distributed in the surface epithelium. Representative whole mount confocal image of ionocytes EGFP(Foxi1) and ciliated cells (AcTub). c. Expression level of ionocyte markers (rows, FDR<0.05 LRT, full-length scRNA-seq dataset) in each airway epithelial cell type (columns). d. EGFP(Foxi1)⁺ ionocytes extend cytoplasmic appendages (arrows). e-g. Immunofluorescence labeling of GFP(Foxi1)⁺ cells in airway regions. Submucosal gland (SMG, e), nasal respiratory epithelium (f) and olfactory neuroepithelium (g). Dotted line separates surface epithelium (SA) from SMG.

Extended Data Figure 10 |. Functional characterization of ionocytes

a. Ascl3-KO moderately decreases ionocyte TFs and Cftr in ALI cultured epithelia. Expression quantification (ΔΔCT, y-axis) of ionocyte (Cftr: −0.82 ΔΔCT, 95% CI [± 0.20], Foxi1: −0.75 ΔΔCT, 95% CI [± 0.28], Ascl3: −10.28 ΔΔCT, 95% CI [± 1.85]) and basal (Trp63), club (Scgb1a1), or ciliated (Foxj1) markers (x-axis) in hetero-and homozygous KO (color legend) are normalized to wild type littermates. The mean of independent probes (p1 and p2) was used for Cftr. n=10 and 5 hetero-and homozygous KO, respectively and n=4 wild type mice. p values: Holm-Sidak test, Methods, error bars: 95% CI. b. Altered airway surface liquid (ASL) reflectance intensity in Foxi1-KO ALI culture compared to WT. Representative μOCT image of ASL. Red bar: airway surface liquid depth (including the periciliary and mucus layers). Scale bar (white): 10μm. c,d. Ionocyte depletion or disruption does not affect ASL depth (c) as determined by μOCT, nor pH (d) in cultured epithelia derived from homozygous Foxi1-KO (n=9, dots) vs. wild type littermates (x-axis, n=9 mice). p values: Mann-Whitney U-test. e,f. Increased ΔI_eq in Foxi1-KO epithelia. ΔI_eq (y axis) in ALI cultures of wild type (WT), heterozygous (HET) and Foxi1-KO mice (n=5 WT, n=4 HET, n=6 KO, dots) that were characterized for their forskolin-inducible equivalent currents (e, I_eq) and for currents sensitive to CFTR_inh-172 (f). The inhibitor-sensitive ΔI_eq values reported may underestimate the true inhibitor-sensitive current, since the inhibitor response failed to reach a steady plateau for some samples during the time scale of the experiment. g-i. Foxi1 transcriptional activation (Foxi1-TA) in ferret increases Cftr expression and chloride transport g. qRT-PCR expression quantification (ΔΔCT, y-axis) of ionocyte markers (x-axis) in ferret Foxi1-TA ALI (n=4 ferrets) normalized to mock transfection (Cftr: −1.39 ΔΔCT, 95% CI [± 0.44], Foxi1: −5.37 ΔΔCT, 95% CI [± 0.91], Ascl3: −0.87 ΔΔCT, 95% CI [± 0.27], Atp6v0d2: −1.18 ΔΔCT, 95% CI [± 0.58] and Atp6v1e1: −0.070 ΔΔCT, 95% CI [± 0.11] Methods), p values: t-test, bars denote means, error bars: 95% CI. h,i. Foxi1 activation in ferret cell cultures results in a CFTR inhibitor-sensitive short-circuit current (ΔI_sc). Representative trace (h) and quantification (i) of short-circuit current (I_sc, y-axis) tracings from Foxi1-TA ferret ALI after sgRNA reverse transfection (n=6, light blue) vs. mock transfection (n=5, black). j. Evolutionarily conserved ionocyte signatures. Difference in fraction of cells in which transcript is detected (x-axis) and log₂ fold-change (y-axis) between human ionocytes and all other bronchial epithelial cells. Labeled genes are differentially expressed (log₂ fold-change>0.25 and FDR<10⁻¹⁰, Mann-Whitney U-test). Red: consensus ionocyte markers between mouse and human (log₂ fold-change>0.25, FDR<10⁻⁵, LRT).

Figure 1 |. A single-cell expression atlas of mouse tracheal epithelial cells

a. Overview. b. t-distributed stochastic neighbor embedding (tSNE) of 7,193 3’ scRNA-seq profiles, colored by cluster assignment (Methods) and annotated post-hoc (legend). Circled: ionocytes.

Figure 2 |. Club cell differentiation varies by location

a. Distal epithelia differentiate into mucous metaplasia. AcTub (ciliated) and Muc5ac (goblet) cells in cultured epithelia from proximal (top) or distal (bottom) trachea stimulated with recombinant IL-13 (right) vs. control (left). b-c. Differentiation trajectories. Diffusion map embedding (b, Methods) of 6,905 basal (blue), club (green), and ciliated (red) cells colored by cluster assignment (top) or expression (log₂(TPM+1), color bar) of Krt13 (bottom). c. Number of individual cells associated with each trajectory (Methods). d-e. Krt13⁺ cells occur in hillock structures. d. Whole-mount stain of Krt13 (magenta) and AcTub (green), n=3 mice. e. Schematic of squamous hillocks within pseudostratified ciliated epithelium.

Figure 3 |. Tracking differentiation dynamics with Pulse-Seq

a,b. Pulse-Seq tracks the lineage labeling of all cell types. tSNE visualization of 66,265 cells colored by cluster assignment (a, color legend), or lineage label (b, mT: membrane-tdTomato, mG: membrane-EGFP). c. mG⁺ lineage-labeled fractions of each tracheal epithelial cell type (%, y-axis, points: individual mice, bars: estimated proportions, Methods), n=3 mice per time-point (x-axis). Error bars: 95% CI, p values: LRT, PI: pulmonary ionocyte. d. Estimated daily rate of new lineage labeled cells (%, y-axis, Methods, Extended Data Fig. 6c) for each type (x-axis), n=9 mice. Error bars: 95% CI, p values: rank test (Methods). e. Validation in situ. Representative images of unlabeled (dashed outline) and basal lineage-labeled (solid outline) Gnat3⁺ tuft cells. f. Cell types, lineage, and cellular dynamics inferred using Pulse-Seq.

Figure 4 |. Tuft and goblet cell subtypes display unique functional gene expression programs

Tuft-1 and tuft-2 sub-clusters. a. Relative expression (RE, row-wise Z-score of log₂(TPM+1); color scale) of genes (rows) differentially expressed (FDR <0.25, LRT) in tuft cells (columns) of each sub-type (top). b. Immunofluorescence validation of pan-tuft marker Trpm5 (blue) and tuft-1 (Gng13⁺, green, top) or tuft-2 (Alox5ap⁺, magenta, bottom) markers (solid outlines) in vivo with DAPI (grey), n=3 mice, replicates=4. c. Distinct expression programs in tuft-1 and tuft-2 cells. Differential expression in tuft cell subtypes for all genes (left), taste genes (center), and leukotriene synthesis genes (right). Labeled genes are differentially expressed (FDR<0.01, LRT), k=892 cells; n=15 mice. d. Immunofluorescence validation of goblet-1 (Tff2, magenta) and goblet-2 (Lipf, green) cells (solid outlines) with DAPI (blue), n=3 mice, replicates=4.

Figure 5 |. The pulmonary ionocyte is a novel mouse and human airway epithelial cell type that specifically expresses CFTR

a. Mouse pulmonary ionocyte markers. Expression level of ionocyte markers (rows, FDR<0.05 LRT, 3’ scRNA-seq dataset) in each airway epithelial cell type (columns). b. Immunofluorescence co-labeling of EGFP(Foxi1⁺) ionocytes (solid outline) with Atp6v0d2 (left) and Cftr (right). c. tSNE plot of 66,265 Pulse-Seq cells and ionocyte subset (black box, inset) colored by expression of ionocyte markers Foxi1 (left) and Cftr (right). d. qRT-PCR confirms ionocyte enrichment of Cftr. Expression (ΔΔCT, y-axis, Supplemental Table 12) of ionocyte (Cftr, Foxi1) and ciliated cell (Foxj1) markers (x-axis) in ionocytes and ciliated cells (legend) isolated from Foxi1-(n=4, dots) and Foxj1-GFP mice (n=3), respectively. Samples normalized to EpCAM⁺ populations from wild-type mice (n=6) Error bars: 95% CI, p values: t-test. e. Foxi1-KO decreases expression of ionocyte TFs and Cftr in ALI cultured epithelia. Expression (ΔΔCT, y-axis, Supplementary Table 12) of ionocyte markers (x-axis) in heterozygous (n=4) and homozygous KO (n=6, color legend), normalized to wild-type littermates (n=8). Error bars: 95% CI, p values: Holm-Sidak test, Methods. f. Foxi1-KO disrupts mucosal homeostasis in ALI cultured epithelia. Effective viscosity (cP, left) and ciliary beat frequency (Hz, right) assayed with μOCT in homozygous Foxi1-KO (n=9, dots) vs. wild-type littermates (x-axis, n=9 mice), bars: means. p values: Mann-Whitney U-test. g,h. Human pulmonary ionocytes are the major source of Cftr in human bronchial epithelium. g. Human ionocytes detected by FISH of FOXI1 and CFTR in bronchi (Methods), n=3 bronchi. h. tSNE of 78,217 3’ droplet scRNA-seq profiles (points) from bronchial epithelium (n=1 patient), colored by their expression of FOXI1 (left) and CFTR (right).

Figure 6 |. Lineage hierarchy of the airway epithelium

Specific cells are associated with novel cell-type markers, pathways, and diseases.