Maladaptive functional changes in alveolar fibroblasts due to perinatal hyperoxia impair epithelial differentiation

Infants born prematurely worldwide have up to a 50% chance of developing bronchopulmonary dysplasia (BPD), a clinical morbidity characterized by dysregulated lung alveolarization and microvascular development. It is known that PDGFR alpha–positive (PDGFRA+) fibroblasts are critical for alveolarization and that PDGFRA+ fibroblasts are reduced in BPD. A better understanding of fibroblast heterogeneity and functional activation status during pathogenesis is required to develop mesenchymal population–targeted therapies for BPD. In this study, we utilized a neonatal hyperoxia mouse model (90% O2 postnatal days 0–7, PN0–PN7) and performed studies on sorted PDGFRA+ cells during injury and room air recovery. After hyperoxia injury, PDGFRA+ matrix and myofibroblasts decreased and PDGFRA+ lipofibroblasts increased by transcriptional signature and population size. PDGFRA+ matrix and myofibroblasts recovered during repair (PN10). After 7 days of in vivo hyperoxia, PDGFRA+ sorted fibroblasts had reduced contractility in vitro, reflecting loss of myofibroblast commitment. Organoids made with PN7 PDGFRA+ fibroblasts from hyperoxia in mice exhibited reduced alveolar type 1 cell differentiation, suggesting reduced alveolar niche-supporting PDGFRA+ matrix fibroblast function. Pathway analysis predicted reduced WNT signaling in hyperoxia fibroblasts. In alveolar organoids from hyperoxia-exposed fibroblasts, WNT activation by CHIR increased the size and number of alveolar organoids and enhanced alveolar type 2 cell differentiation.

Contraction was calculated with Image J software by measuring the diameter of the pellet, followed by a 2-tailed Student's t-test (Figure 4 L) or ANOVA followed by Tukey's multiple comparison (Figure 6 A).

RNA-seq and bioinformatic analysis
MACS-sorted PDGFRA + (CD140 + ) fibroblasts from either RA or O2 PN4, PN7, and PN10 mouse lungs were prepared for bulk RNA-seq. RNA sequencing was performed by Cincinnati Children's Hospital Medical Center's Gene Expression Core utilizing the Illumina HiSeq2500. RNA-seq FASTQ files were aligned using Bowtie to mouse genome version mm10 (4). Raw gene counts were obtained using Bioconductor's Genomic Alignment, which were subsequently made into normalized FPKM values using Cufflinks (5,6). DeSeq (Bioconductor) was used to calculate differential gene expression from raw expression values. Genes were deemed differentially expressed if they satisfied the following requirements: gene has a fold change >2, binomTest pvalue <.01, and RPKMs >1 in 2 of the 3 replicates in at least one condition being compared. Gene patterns were determined by comparing differential expressed genes from all three time points. Genes that were significantly changed or unchanged at the same time points were grouped together. Heatmaps of genes in particular patterns were z-score normalized and generated using Partek Genomics Suite (http://www.partek.com/pgs). Gene enrichment analysis was carried out using ToppGene's ToppFun, and functional enrichments within each profile were identified and all profiles were compared to each other using Toppcluster (7). Pvalues of functional enrichment hits where -log10 transformed for graphical visualization. Signature genes for Matrix and Myo fibroblast were determined by downloading and comparing PN7 and PN10 markers genes for respective cell types from LGEA https://research.cchmc.org/pbge/lunggens/mainportal.html (8). Only signatures genes identified at both time points for a particular cell type were assessed.
Top 50 Lipofibroblast signature genes were obtained from a recently published mouse lung scRNA-seq study (9). Fold changes of any Matix, Myo or Lipo fibroblast marker genes significantly altered for a particular cell type, at any time point, were visualized in a heatmap generated by Pheatmap (https://cran.rproject.org/web/packages/pheatmap/index.html). Wnt related genes changes and predictive network was generated by Qiagen Ingenuity Pathway Analysis (IPF) using genes significantly altered at D4, D7 and/or Day10 (10).

Immunofluorescence
Immunofluorescent staining was performed on 5-µm slides sectioned from paraffinembedded lung tissue blocks. Slides were deparaffinized in xylene, rehydrated in a series of graded ethanol and washed in 1X PBS. When required, antigen retrieval in 10 mM citrate buffer (pH 6.0) was performed. Non-specific antigens were blocked in 4% normal donkey serum in PBS with 0.1% Triton X-100 (PBST) for 2 hours. Slides were incubated in primary antibodies (Supplemental Table 1) diluted in blocking buffer overnight at 4°C.
After washing in PBST, slides were incubated in fluorescent secondary antibodies (1:200) and DAPI (1µg/ml) diluted in blocking buffer for 1 hour at room temperature. Slides were subsequently washed in PBST and mounted in Prolong Gold (Thermo Fisher Scientific).
Z-stack images were captured on a Nikon AR1 inverted confocal microscope and further analyzed using Nikon Elements software. Antibody staining was quantified on Nikon Elements using a previously established protocol (11).

In vitro r-spondin treatment of IMR90 cells
Human IMR90 fibroblasts (ATCC ® CCL-186 ™ ) were cultured in growth medium (DMEM Hams F12; 10% FBS; 0.1% penicillin and streptavidin; 0.1% gentamycin and amphotericin) until passage 3 in a 100mm tissue culture plate. Cells were seeded on a 24-well plate, and upon reaching 70% confluency, treated with r-spondin conditioned media at a dilution of 1:500. Cells were harvested at 12, 24, and 36 hrs after r-spondin treatment and processed for RNA isolation and RT-qPCR gene expression analysis using the same methods as previously described.

Morphometrics
Morphometric analysis with a FIJI-macro was used to quantify alveolar structure including volume density of alveolar septa (Vvsep), mean linear intercept of airspaces (Lm), mean transsectional wall length (Lmw), and surface area density of airspaces (Svair) on PN4, PN7 and PN10 room air and hyperoxia exposed mouse lungs. Lungs were inflation-fixed and processed as previously described and 7-μm thick sections were stained with hematoxylin and eosin for analysis. FBs were treated for 24hrs continuously with CHIR after plating. Primary fibroblasts in collagen contraction assays were treated for 24hrs continuously with CHIR or PDGF-AA after seeding on top of collagen. Organoids were allowed to grow for 7 days before treatment with either CHIR or PDGF-AA added to the organoid media, continuously for 14 days (with media changes every two days).
Supplemental Table 1 List of primary antibodies used in this study Supplemental Table 2 List of Taqman primers used in this study  Supplemental Figure 6. Basal cells, proximal, and distal progenitors unaffected in the epithelium in organoids made with hyperoxia-exposed primary PDGFRA + fibroblasts A. Immunofluorescence images of organoids made from PN7 room air-and hyperoxia-exposed PDGFRA + fibroblasts and adult epithelial cells. Stained with P63, SOX2, CCSP, ECAD, and SOX9. Scale bar=100 μm B. Percentage of P63 + DAPI + nuclei over total DAPI + nuclei. C. Percentage of SOX2 + DAPI + nuclei over total DAPI + nuclei. D.
Percentage of SOX9 + DAPI + nuclei over total DAPI + nuclei. In panels (B-D), a 2-tailed Student's t test was used. E. ARL13B immunofluorescence in organoids quantified using Nikon elements as antibody area over DAPI area. N=3-10 organoid transwells (replicates) used, 3 slides per transwell. One-way ANOVA followed by Tukey's multiple comparison was used to determine significance. F. Quantification of brightfield images of organoids using Nikon Elements software by size, organoids >500um diameter. One-way ANOVA followed by Tukey's multiple comparison was used to determine significance between three or more groups.