Single cell transcriptomics in a treatment-segregated cohort exposes a STAT3-regulated therapeutic gap in idiopathic pulmonary fibrosis

Abstract

Idiopathic pulmonary fibrosis (IPF) is a progressive fibrotic pulmonary disease with unknown etiology. Since approved idiopathic pulmonary fibrosis (IPF) drugs only slow disease progression, novel therapeutics are required that improve clinical outcomes. Here, we report a single cell RNA-Seq and regulatory network analysis of the largest IPF cohort assembled to date. Segregating this cohort based on status of treatment with approved antifibrotics (untreated, nintedanib- and pirfenidone-treated), we describe for the first time the transcriptional landscape of untreated IPF across 40 lung cell types, and the elements of this program that are impacted by approved antifibrotics. On average, 60% of the untreated IPF-dysregulated transcriptome is refractory to treatment with these drugs, a transcriptional deficit we refer to as the IPF therapeutic gap. Regulatory network analysis indicated a dominant functional footprint for the transcription factor STAT3 in both untreated IPF and in the IPF therapeutic gap. Validating our analysis in a translational precision cut lung slice platform that recapitulates IPF explants, treatment with a STAT3 inhibitor reduced the IPF therapeutic gap in numerous lung cell types. Finally, we implicated STAT3 as a master transcription factor that regulates a network comprising numerous profibrotic transcription factors in IPF alveolar fibroblasts, a critical fibrotic cell lineage. Our study represents a comprehensive resource for translational lung fibrosis research and establishes a novel strategy for drug discovery in human disease more broadly.

Figure 1.. A treatment status single cell transcriptional Atlas of IPF.

(A) Overview of experimental design. (i) Collection of disease lung explants (n = 74) and donor lungs (n = 64). The IPF cohort was segregated into three subcohorts: untreated IPF (n = 22), pirfenidone-treated IPF (n = 24) and nintedanib-treated IPF (n = 28). (ii) Dissociation into single-cell suspension. (iii) scRNA-seq library preparation and sequencing. (iv) Differential gene expression analysis. (v) Regulatory network and pathway analysis. (vi) Prioritization and validation of STAT3 as a novel clinical IPF target. (B-E) UMAP plots and accompanying bar charts showing effect of approved antifibrotics on cell counts in untreated IPF for selected cell types: (B) PTBSC and AT2 cells; (C) MACA and MACMOND; (D) FBALV and FBPERIB; (E) VEN and CAPG. C, control; UNT, untreated IPF; TRT, approved antifibrotic-treated IPF.

Figure 2.. The transcriptional landscape of untreated IPF at single cell resolution.

(A) Hierarchical clustering of the top significantly enriched GO terms in untreated IPF-induced gene sets. GO enrichment analysis was performed on genes induced in untreated IPF (log₂FC>0.32, FDR<0.1) relative to controls. Color palette represents −log10p enrichment scale from p adj < 0.05 (light orange) to p adj <1E-70 (deep red); blue indicates no significant enrichment. Terms that were significantly enriched (p adj < 0.05) in at least six cell types were subjected to hierarchical clustering as described in the methods section. Panels B-H represent predicted TF regulatory network plots for unique gene sets from the indicated panel A GO term clades across all cell types; in these plots log₂ odds ratio is on the x axis and −log₁₀p is on the y axis, such that TFs that have the strongest and most significant footprints within the gene set of interest are distributed towards the top right of the plot. (B) Regulatory network plot for clade 1. (C) Regulatory network plot for clade 7. (D) Regulatory network plot for clade 8. (E) Regulatory network plot for clade 10b. (F) Regulatory network plot for clade 12. (G) Regulatory network plot for clade 16. (H) Regulatory network plot for clade 25. Full GO data are in table S4; full clade regulatory network data are in table S5.

Figure 3.. Single cell transcriptional biology of approved antifibrotics in IPF.

(A) Stacked bar chart representation of differentially expressed genes derived from IPF Atlas scRNA-Seq data. Bars represent total numbers of untreated IPF vs. control-dysregulated (log2 FC>0.32 or <−0.32, FDR<0.1) genes in each cell type at the indicated cutoffs. Subsets of these genes whose dysregulation is reversed by nintedanib only (green), pirfenidone only (blue), both antifibrotics (yellow) and neither (red) are indicated. Full data are in table S6. (B) Heatmap indicating the extent to which untreated IPF-induced GO terms are reversed by treatment with nintedanib or pirfenidone in each cell type. The 508 untreated IPF-induced enriched GO terms from Fig. 2A are arranged in the identical order on the vertical axis; the 40 Atlas cell types are arranged on the horizontal axis. Within the heatmap, the color palette represents −log10p enrichment scale from p adj < 0.05 (yellow) to p adj <1E-25 (red); blue indicates no significant enrichment. On the horizontal edges, for a given cell type-GO term intersection, UNT (red) indicates enrichment of that GO term in untreated IPF-induced genes; NTD (blue) indicates enrichment of that GO term in nintedanib-IPF-repressed genes; and PFD (blue) indicates enrichment of that GO term in pirfenidone-IPF-repressed genes. Full data are in table S4. Panels C-G represent scatterplots of pirfenidone-IPF v untreated IPF (x axis) and nintedanib-IPF vs. untreated IPF (y axis) log₂ FC values for untreated IPF-induced genes mapped to the following GO terms: (C) cadherin binding in AT2 cells; (D) immune response-activating signal transduction in NK cells; (E) secretory granule lumen in MACIPs; (F) SMAD binding in FBADs; and (G) GTPase regulator activity in CAPA cells.

Figure 4.. STAT3 is a dominant transcription factor in untreated IPF and in the IPF therapeutic gap.

(A) Heat map showing predicted top 50 TFs in untreated IPF-induced pan-Atlas regulatory network. Regulatory network analysis was performed on untreated IPF-induced genes in the 40 cell types. TFs are ranked based on (i) number of significant footprints (and (ii) average footprint size. Full data are in table S9. Color palette represents −log10p enrichment scale from p adj < 0.05 (yellow) to p adj <1E-35 (red); blue indicates no significant enrichment. (B) Predicted STAT3 footprint in untreated IPF-induced genes across the 40 cell types. (C) Regulatory network plot showing predicted TFs driving severe IPF in the LTRC dataset GSE32537. Full data are in table S10. (D) Regulatory network plot for untreated IPF-induced genes in Atlas alveolar fibroblasts showing TFs that are transcriptionally induced in untreated IPF. (E) Volcano plot of untreated IPF vs. control Atlas alveolar fibroblast DEGs highlighting genes encoding TFs labelled yellow in panel D. (F) Master regulatory network plot showing TFs predicted to regulate genes encoding untreated IPF-induced network TFs in alveolar fibroblasts. Full data are in table S14. (G) Same plot as panel D highlighting TFs encoded by genes that were hits in an RNAi screen of lung fibroblast differentiation (Oh et al.).

Figure 5.. Characterization of the IPF-induced therapeutic gap.

(A) Heatmap showing enrichment of prioritized GO terms in the IPF-induced gap (i.e., genes induced in untreated IPF but not repressed by either antifibrotic) in each cell type. GO terms were prioritized as described in the text. Color palette represents −log10p enrichment scale from p adj < 0.05 (yellow) to p adj <1E-25 (red); blue indicates no significant enrichment. Panels B-F indicate induction in untreated IPF vs. control of genes mapped to the following gap-prioritized GO terms: (B) cell-cell junction in AT2 cells; (C) T cell receptor signaling pathway in TCD8s; (D) cell leading edge in MACMONDs; (E) collagen-containing extracellular matrix in MFBs; (F) endothelial cell development in CAPAs. Full data are in table S8. (G) Scatterplot comparing TF rankings in the untreated IPF-induced (x axis) and IPF-induced gap (y axis) pan-Atlas regulatory networks. Full data are in table S11. (H) Regulatory network plot for the IPF-induced gap in Atlas IPF alveolar fibroblasts.

Figure 6.. STAT3 inhibition reduces the IPF therapeutic gap in precision cut lung slices

(A) Bar chart showing enrichment −log10p values for the intersection of PCLS IPF-induced and Atlas untreated IPF-induced gene sets across common cell types. (B) Volcano plot showing enrichment of Atlas untreated IPF-induced alveolar fibroblast genes among PCLS IPF-alveolar fibroblast induced genes (C) Functional footprint of STAT3 among TTI-101-repressed genes in PCLSs. Dotted line indicates p adj < 0.05. (D) Manhattan plot comparing the effect on PCLS alveolar fibroblast IPF-induced genes of treatment with nintedanib, pirfenidone and TTI-101. (E) Volcano plot showing enrichment of PCLS alveolar fibroblast IPF vs. CON-induced genes among TTI-101-repressed PCLS alveolar fibroblast genes. (F) Manhattan plot comparing the effect on PCLS AT2 IPF-induced genes of treatment with nintedanib, pirfenidone and TTI-101. (G) Manhattan plot comparing the effect on PCLS capillary aerocyte IPF-induced genes of treatment with nintedanib, pirfenidone and TTI-101. Panels H-L compare numbers of IPF-induced and nintedanib-, pirfenidone- and TTI-101-repressed genes mapping to the indicated IPF-induced gap-prioritized GO term and PCLS cell type; the right panels represent the TTI-101-IPF vs. IPF volcano plot for the corresponding cell type showing TTI-101 repression of the indicated GO term genes: (H) cell-cell junction in AT2 cells; (I) T cell receptor signaling pathway in TCD8s; (J) cell leading edge in MACMONDs; (K) collagen-containing extracellular matrix in MFBs; (L) endothelial cell development in CAPAs. See fig. S4D-H for corresponding TTI-101 volcano plots.

Figure 7.. STAT3 inhibition disrupts the IPF alveolar fibroblast regulatory network

(A) Regulatory network plot for IPF-induced genes in PCLS alveolar fibroblasts showing enrichment of top-ranked Atlas alveolar fibroblast IPF network TFs. Full data are in table S15. (B) Volcano plot of IPF v control PCLS alveolar fibroblast DEGs showing enrichment of genes encoding the 23 untreated IPF-induced FBLAV network TFs. (C) Scatterplot comparing −log₁₀p footprints of Atlas (x axis) and PCLS (y axis) alveolar fibroblast IPF master regulatory network TFs. Full data for atlas master network are in table S16. (D) Manhattan plot comparing effect of approved antifibrotics and TTI-101 on IPF-induced genes encoding the 23 untreated IPF-induced FBALV network TFs. (E) Volcano plot of TTI-101-IPF vs. IPF PCLS alveolar fibroblast DEGs highlighting TTI-101 repression of the 23 IPF-induced alveolar fibroblast IPF network TF genes. (F) Regulatory network plot for TTI-101-repressed alveolar fibroblast genes showing enrichment of top TFs in regulatory network for IPF-induced genes. Full data are in table S15. (G) Generalized schematic of STAT3 predicted mechanism of action in IPF alveolar fibroblasts supported by the Atlas and PCLS analysis. STAT3 has a dual role to (i) co-ordinate transcriptional induction of genes encoding TFs predicted to form IPF regulatory networks across numerous cell compartments; and (ii) as a functional member of many of those networks. Although we have focused here on STAT3, autoregulation has been reported for many other members of the FBALV master regulatory network, including AP-1 members and NFATC1, and this mechanism likely contributes to replenishment of the master network in IPF.