Long-Term In Vitro Expansion of Epithelial Stem Cells Enabled by Pharmacological Inhibition of PAK1-ROCK-Myosin II and TGF-β Signaling.

Despite substantial self-renewal capability in vivo, epithelial stem and progenitor cells located in various tissues expand for a few passages in vitro in feeder-free condition before they succumb to growth arrest. Here, we describe the EpiX method, which utilizes small molecules that inhibit PAK1-ROCK-Myosin II and TGF-β signaling to achieve over one trillion-fold expansion of human epithelial stem and progenitor cells from skin, airway, mammary, and prostate glands in the absence of feeder cells. Transcriptomic and epigenomic studies show that this condition helps epithelial cells to overcome stresses for continuous proliferation. EpiX-expanded basal epithelial cells differentiate into mature epithelial cells consistent with their tissue origins. Whole-genome sequencing reveals that the cells retain remarkable genome integrity after extensive in vitro expansion without acquiring tumorigenicity. EpiX technology provides a solution to exploit the potential of tissue-resident epithelial stem and progenitor cells for regenerative medicine.

INTRODUCTION

Tissue-resident stem cells ensure homeostasis and tissue repair throughout the lifetime of an individual. In various epithelia, the stem and progenitor cells residing in the basal layer are marked by KRT5 and TP63 and have infinite self-renewal capability in vivo (Blanpain and Fuchs, 2014; Donati and Watt, 2015; Hogan et al., 2014; Rock et al., 2010). However, it has been difficult to extensively expand epithelial cells in vitro in feeder-free condition due to the CDKN2A-dependent stasis (Shay and Wright, 2007). Immortalization using telomerase reverse transcriptase (TERT) or viral genes (SV40T or HPV16 E6/E7) significantly alters epithelial cells behavior, limiting their utility for studying normal biology or as drug-screening models (Miller and Spence, 2017). Lack of suitable long-term expansion methods has hampered epithelial stem cell biology study in vitro and greatly stalled advances in regenerative medicine exploiting their potential. Pluripotent stem cells (PSCs), including induced PSCs, have been the subject of intense research in the hope that they offer physiology-relevant models and solutions for regenerative medicine. However, they face challenges including donor variability, acquired oncogenic mutations, and inefficient differentiation toward mature cell types (Avior et al., 2016; Merkle et al., 2017).

Encouraging progress has been made in developing methods that allow continuous in vitro propagation of epithelial cells. Liu et al. proposed that feeder cells and Rho-kinase (ROCK) inhibitor Y-27632 “conditionally reprogrammed” (CR) epithelial cells to proliferate continuously (Butler et al., 2016; Chapman et al., 2010; Liu et al., 2012; Suprynowicz et al., 2012). The Stingl group used a similar approach to expand mammary repopulating units, an indication of the expansion of mammary epithelial progenitors (Prater et al., 2014). The CR method has garnered interest due to its successful use in expanding patient-derived epithelial cells to identify effective therapy (Crystal et al., 2014; Yuan et al., 2012). Wang et al. (2015) used feeder cells and several small molecules regulating TGF-β, WNT, and NOTCH pathways to expand “ground-state intestinal stem cells.” However, the use of feeder cells complicates the interpretation of signaling events that govern cell proliferation and creates challenges in meeting regulatory expectation for manufacturing cell therapy products (Lipsitz et al., 2016).

The Clevers group has led the way in developing feeder-free 3D organoids for intestinal stem cells (Sato et al., 2009, 2011), which has later expanded to epithelial cells from liver, pancreas, and stomach (Boj et al., 2015; Huch et al., 2013, 2015). Stem cells, progenitors, and differentiated epithelial cells are present in the organoid, making it a good in vitro model for epithelial cell biology. Katsuda et al. (2017) reported the use of small molecules, including Y-27632, A83-01, and CHIR99021, which converted rodent hepatocytes into proliferative bipotent cells; however, it did not work for human hepatocytes.

To develop medium formulations that address aforementioned issues, including safety, reproducibility, and scale-up compatibility, we set off to identify small molecules that support long-term epithelial cell expansion without feeder cells. We found that the combination of TGF-β signaling inhibition, PAK1-ROCK-Myosin II inhibition, and low extracellular [Ca2+] were key components that transformed traditional culture medium to enable long-term propagation of epithelial cells from various tissues. High single-cell cloning efficiency and the ability to differentiate into tissue-specific mature epithelial cell types suggested that stem and progenitor cells were preserved during expansion. Remarkably, the cells retained genome integrity with no tumorigenic mutations after extensive expansion as assessed by multiple approaches including whole-genome sequencing. Gradual changes in DNA methylation landscape were the by-product of long-term culture and had little impact on overall gene expression profile.

RESULTS

TGF-β Signaling Inhibition and ROCK Inhibition Synergistically Support Long-Term Epithelial Cell Expansion in the Absence of Feeder Cells

As epithelial cells quickly cease proliferation when the feeder cells or Y-27632 are omitted from the CR method (Liu et al., 2012), we designed a proliferation assay to pharmacologically screen a focused collection of small molecules modulating diverse biological pathways governing stem cell self-renewal and differentiation (Table S1) to develop feeder-free condition that supports continuous cell proliferation. Human prostate epithelial cells (PrECs), bronchial epithelial cells (HBECs), and foreskin keratinocytes (HFKs) were transduced by lentiviruses expressing nuclear-localized red fluorescent protein (nRFP) to facilitate automatic cell count. The lentivirus transduction did not alter cell proliferation. Late-passage cells were cultured for 7 days in the absence of feeder cells and Y-27632 to screen small molecules. We found that small molecules that inhibit TGF-β signaling (A83-01, RepSOX, GW788388, SB431542) or ROCK (SR3677, Y-27632, Thiazovivin, GSK429286) supported continuous cell proliferation at micromolar concentrations in the F medium (Figure 1A). Small molecules that affect the WNT pathway by inhibiting GSK3 (BIO, CHIR99021, BIO-acetoxime, endo-IWR1), or inhibit Abl kinase (AP24534), or increase intracellular cAMP level (forskolin), or target the NOTCH pathway (DAPT), also supported continuous cell proliferation to certain degrees (Figures 1A and S1). Importantly, we found that TGF-β signaling inhibitor (0.5–2 μM A83-01) and ROCK inhibitor (5–10 μM Y-27632) synergistically promoted epithelial cell proliferation (Figures 1B and S1). PrECs grew for >28 population doublings (PDs) in the F medium plus A83-01 and Y-27632 (F+Y+A) (Figure S1), much longer than in the control PrGM. Still, the cells stopped proliferation after several passages in F+Y+A medium, which was significantly shorter than in the CR condition (Figures 1C and S1) and suggested that additional optimization could further improve their expansion.

Figure 1.

TGF-β Signaling Inhibition, ROCK Inhibition, and Low Calcium Synergistically Support Long-Term Epithelial Cell Proliferation

(A) Small molecules inhibiting the TGF-β signaling or ROCK supported the proliferation of late-passage PrECs/nRFP cells in the absence of feeder cells in the F medium. Data are represented as mean ±SD, n = 4.

(B) Synergy between A83-01 and Y-27632 in the F medium (four replicates per condition).

(C and F) PrECs/nRFP cells proliferated for 10 PDs in the F medium plus Y-27632 and A83-01 (F+Y+A) but continued to proliferate in the CR condition (C). Many cells in F+Y+A exhibited differentiated morphology (F).

(D) ROCK inhibitors synergistically promoted the proliferation of HBECs/nRFP cells in KSFM plus 1 μM A83-01.

(E) Synergy between A83-01 and Y-27632 in KSFM (four replicates per condition).

(G) Morphology of HFKs over successive passages in KSFM (P1 and P5) or the EpiX medium (P2, P11, and P20).

(H) TP63 was ubiquitously expressed in late-passage HFKs (P16) cultured in the EpiX medium.

(I-L) Expansion of HFKs (I), PrECs (J), HBECs (K), and mammary epithelial cells (L) in KSFM or EpiX.

Low Calcium Concentration Enhances the Growth-Promotion Effect of TGF-β Signaling Inhibition and PAK1-ROCK-Myosin II Inhibition

F medium contains ~1.4 mM CaCl2, which is known to promote epithelial cell differentiation (Hennings et al., 1980; Martin et al., 1991), so we conjectured that low [CaCl2] might further promote cell growth. We used keratinocyte-SFM (KSFM) that contains 90 μM CaCl2 to pan the small molecules and observed similar synergy between TGF-β and ROCK inhibitors (Figures 1D, 1E, and S2). HBECs and PrECs achieved >1,000,000-fold expansion in KSFM plus A83-01 and Y-27632 (K+A+Y) versus control (Figure S2). This confirmed that neither TGF-β nor ROCK inhibition alone supports extended epithelial cell proliferation (Chapman et al., 2014; Natarajan et al., 2006). In KSFM, TGF-β inhibition had the biggest impact on cell growth, while in the F medium ROCK inhibition was the primary hit (Figures 1D and 1E). Interestingly, elevating the [Ca2+] in KSFM reduced the growth-promoting effect of A83-01 but boosted the effect of Y-27632 (Figure S3). Additionally, we found that IPA-3 (inhibitor for PAK1) or blebbistatin (inhibitor for Myosin II) synergistically promoted cell proliferation with A83-01 (Figures 1D and S2). HBECs, keratinocytes (neonatal and adult), and PrECs expanded for >1,000,000-fold in KSFM plus A83–01 and IPA-3 or blebbistatin over KSFM control (Figure S2). These findings suggested that, in the absence of feeder cells, inhibiting PAK1-ROCK-Myosin II axis promoted epithelial cell proliferation in low-calcium medium when the TGF-β signaling was also suppressed.

We sieved the small molecules with the K+A+Y medium and found some agents that further promoted cell proliferation, including IBMX, 8-bromo-cAMP, prostaglandin E2, and DPPIV inhibitor (Figure S2). Notably, the synergy between A83-01 and Y-27632 was independent of KSFM, as they also supported long-term epithelial cell expansion in BEGM (Lonza), Bronchia-Life (Lifeline), or LHC-9 (GIBCO) (data not shown). We used KSFM plus 1 μM A83-01, 5 μM Y-27632, and 3 μM isoproterenol and dubbed it epithelial expansion medium (EpiX). Panning the small molecules using the EpiX medium did not reveal any significant hits.

We have established epithelial cell cultures from neonatal foreskin, adult skin, nasal, bronchial, prostate, mammary, and other tissues using the EpiX medium and propagated them for 40–90 PDs (Figures 1I-1L and S4), far beyond the 15–25 PDs obtained in conventional media. We routinely generated >50 million epithelial cells from nasal brushing samples in 3–4 weeks (Figure S4). EpiX medium supported multi-round single-cell cloning with high efficiency, suggesting that epithelial stem and progenitor cells were preserved during the clonal expansion. This allowed us to knock the GFP gene into the AAVS1 locus in HFKs using CRISPR/Cas9 for clonal expansion (Figure S5).

Transcriptome Dynamics of Epithelial Cells Expanded in the EpiX Medium

We examined HFKs and HBECs cultured in the EpiX medium by qRT-PCR and found the cells maintained steady expression of basal cell markers (ITGA6, ITGB4, KRT14, KRT5, TP63; Table S2). Immunofluorescence staining demonstrated that all cells strongly expressed TP63 (Figure 1H). PSC markers (LIN28A, NANOG, OCT4) and other stem cell markers (CD34, AXIN2, LGR5, PROM1) were not activated, indicating that EpiX did not reprogram the cell identity (Table S2). Markers of differentiated cells (IVL, LOR, MUC5AC, SFTPB, SGCB1A1, FOXJ1) were absent in the proliferating cells. Interestingly, qRT-PCR assay detected weak expression of TERT in PrECs and HBECs cultured in EpiX similar to that in the CR method (Figure S6). Genes involved in stress responses and senescence (AKT1, ATM, CDKN2A, GADD45A, GLB1, PLAU, SERPINE1) expressed at much lower levels in cells cultured in the EpiX than in KSFM (Figure S6), suggesting that the EpiX medium alleviated various stresses.

To assess the global impact of EpiX medium on gene expression, we compared a total of 14 HFK transcriptomes under various conditions by RNA sequencing (RNA-seq) (Figure 2A; Table S3), which recapitulated the sustained expression of basal cell markers and absence of PSCs and terminal differentiation markers (Figure S6; Table S2). Globally, two transcriptomic shift trends were illustrated by principal-component analysis (PCA) (Figure 2B), with a major shift associated with media change and a minor shift associated with long-term culture. We identified 1,277 differentially expressed genes in HFKs expanded in EpiX versus in KSFM, comprising 315 upregulated and 962 downregulated genes (Figure 2C; Table S3). Downregulated genes were enriched in cell-cell interaction, interferon signaling, extracellular matrix organization, and cytokine signaling (Figure 2F). As expected, we observed significant downregulation of genes involved in senescence (PLAU, SERPINE1, VIM) and the TGF-β pathway (FN1, SKIL, ITGB6) (Figures 2C, 2D, and S6). Importantly, the downregulated genes were de-repressed when HFKs were withdrawn from the EpiX medium (Figure 2E), suggesting that repression of these genes in EpiX was not a permanent but reversible process.

Figure 2.

Transcriptome Analysis of HFKs Expanded with the EpiX Medium

(A) Experimental scheme of HFK expansion in KSFM or the EpiX medium and time points for samples collection. HFKs underwent three different routes: (1)stasis in KSFM (P2→P5); (2) transient exposure and withdrawal from the EpiX medium (P2→P3→P4); and (3) expansion in the EpiX medium and withdrawal at late passage (P3→P12→P19→P20). Each sample was collected in duplicate.

(B) PCA on HFK transcriptomes in various conditions. The red arrows represented the direction of culture condition changes.

(C) MA plot of differentially expressed genes in two different media. The red dots represented differentially expressed genes.

(D) Temporal expression changes of two top downregulated genes, VIM (solid line) and FN1 (dashed line).

(E) Heatmap of 962 downregulated genes in HFKs expanded with the EpiX medium. These downregulated genes were de-repressed when the HFKs were withdrawn from EpiX medium. The red arrows indicated EpiX medium withdrawal.

(F) Gene ontology (GO) and pathway enrichment analysis of upregulated and downregulated genes by metascape tool.

Epithelial Cells Expanded with EpiX Medium Retain Remarkable Genome Stability

Epithelial cells generally encounter stasis around passage 5–6 in conventional media (Shay and Wright, 2005, 2007), while a few cells may acquire spontaneous suppression of CDKN2A to evade stasis and eventually result in chromosome abnormality and oncogenic mutations (Romanov et al., 2001). We examined HFKs, HBECs, and PrECs cultured for over 40 PDs in the EpiX medium and found they all retained normal diploids (Figure 3A). Whole-genome sequencing (WGS) of HBECs from a cystic fibrosis patient (CF) and HFKs only found 5,083 de novo single-nucleotide variations (SNVs) in the CF sample after 43 PDs, and 4,800 de novo SNVs in the HFK sample after 42 PDs (Figure 3D; Table S4). Average SNV rate was 1.97 × 10−8/bp/generation in the CF and 1.90 × 10−8/bp/generation in the HFKs, close to the 1.5 × 10−8/bp/generation rate observed in germline cells (Rahbari et al., 2016). Less than 0.4% (39/9,844) of those SNVs were common between the CF and HFK samples and they all located in non-coding areas and none had known detrimental effects (Figure 3D; Table S4). Over 99.2% of the 9,844 SNVs located outside coding regions, and the remaining (25 in CF and 24 in HFK) led to heterozygous missense variants (Table S4). We did not detect any mutation in oncogenes (e.g., MYC, RAS) or tumor suppressor genes (e.g., TP53, RB), suggesting that the EpiX medium neither induced nor favored oncogenic mutations. Whole genome-sequencing revealed 60 small insertions and deletions (InDel) in the CF sample and 26 InDel in the HFK sample, all located in introns or intergenic regions (Table S4). Twenty-four copy number variations (CNVs) affecting 1.7 Mb (CF) and one 86-kb CNV (HFK) all located in repeat elements (Table S4; Figure S7). We found one structural variation (SV) in the CF sample and two in the HFK sample (Table S4). None of these CNVs or SVs are functionally linked to tumors. We injected 1 × 107 young HFKs (<20 PDs) differentiated in a high-calcium medium for 7 days (see below), old HFKs (>50 PDs) differentiated for 7 days, or old HFKs (>50 PDs) in active proliferation subcutaneously into nude mice, and found no tumors after 12 weeks (Figure 3B), indicating that EpiX-expanded epithelial cells were not tumorigenic.

Figure 3.

Genome Stability of Epithelial Cells Expanded with the EpiX Medium

(A) HFKs, HBECs, and PrECs expanded with the EpiX medium over 40 PDs retained diploid karyotypes.

(B) EpiX-expanded HFKs did not form tumors after 12 weeks in female nude mice (1 × 107 cells per mouse, n = 6 in each group). Data are represented as mean ±SD, n = 6.

(C) Data analysis workflow for whole-genome sequencing.

(D) Distributions of de novo SNVs identified by whole genome-sequencing in late-passage HFK and CF samples. Venn diagram showed the overlap between the CF and HFK samples.

(E) Genes that were affected by missense de novo SNV (all heterozygous) after 40 PDs expansion in the EpiX medium. Two separate missense SNVs were found in the DDX11 gene in the CF and HFK samples.

Tissue-Restricted Differentiation of Epithelial Cells after Long-Term Expansion in EpiX Medium

Increasing [Ca2+] above 1 mM in the EpiX medium quickly led to the formation of tight junctions stained positive for ZO-1 and Occludin (Figures 4A and 4B). Mini-domes occurred sporadically in confluent HFKs or HBECs cultured in EpiX with 1 mM CaCl2 (Figures 4C and 4D). An intact epithelium sheet could be released by Dispase to create cultured epidermal grafts (Green, 2008) (Figure 4E). HFKs cultured at the air-liquid interface (ALI) matured into stratified epithelium after 2 weeks (Figure 4F). As expected (Doran et al., 1980), EpiX-expanded HFKs injected subcutaneously into immunodeficient mice formed epidermis-like cystic epithelium after 5 weeks, with basal cell layer stained positive for KRT14 and Ki67 (Figures 4H and 4I), spinous cell layers, a granular cell layer with evident keratohyaline granules, and several cornified layers (Figure 4G). Bioluminescence imaging (the HFKs were engineered to express firefly luciferase) confirmed the presence of HFKs 4 months after engraftment in the mice, indicating that EpiX-expanded HFKs retained self-renewal capability in vivo (Figure S8).

Figure 4.

Differentiation of HFKs Expanded with the EpiX Medium

(A) Addition of 1 mM CaCl2 to the EpiX medium induced the HFKs to differentiate in 24 hr.

(B) Immunofluorescence staining of tight junctions (ZO-1 and Occludin) in HFKs cultured in EpiX plus 1 mM CaCl2 for 7 days. ZO-1, zonula occludens-1.

(C and D) Confluent HFKs in EpiX (C) or EpiX plus 1 mM CaCl2 (D). Many domes with liquid accumulated underneath occurred in the culture with EpiX plus 1 mM CaCl2.

(E) HFKs cultured in a T-75 flask in EpiX plus 1.5 mM CaCl2 for 7 days form an intact epithelium sheet, which was released from the flask after 30-min incubation in dispase at 37°C.

(F) HFKs were differentiated at ALI for 14 days and formed a stratified epithelium.

(G) EpiX-expanded HFKs were subcutaneously injected into immune-comprised mice. After 5 weeks, the cells formed cysts which resembled epidermis.

(H) Most cells in the basal layer were stained positive for the proliferation marker Ki67.

(I) The basal and supra-basal layers of the cystic epithelium stained positive for KRT14 (K14).

HBECs cultured in standard BEGM lose the ability to differentiate into mucociliary epithelium after four to five passages (Neuberger et al., 2011). HBECs from several healthy and CF donors grew for 45–60 PDs (12–16 passages) in the EpiX medium (Figure 5A), and readily differentiated into mucociliary epithelium at ALI or bronchospheres in Matrigel. Cilia of the multiciliated cells spontaneously beat (Videos S1 and S2) and stained positive for acetylated tubulin, and mucin-producing goblet cells stained positive by alcian blue or an anti-MUC5AC antibody (Figures 5C and 5D). High expression of differentiated cell markers (FOXJ1A, MUC5AC, SCGB1A1) was revealed by qRT-PCR (Figure 5B). The mucociliary epithelium also maintained physiological functions, showing a significant increase in MUC5AC+ goblet cells upon IL-13 stimulation (Figures 5D). Additionally, ion channel physiology expected for their CFTR genotypes was maintained after >30 PDs (Figure 5E). For HBECs from a healthy donor, CFTR activity was readily stimulated by 10 μM forskolin, and diminished by the inhibitor CFTRinh-172. For CF cells with a rare CFTR variant (ΔF508:Q685TfsX4), the CFTR corrector (3 μM VX-809) with the potentiator (100 nM VX-770) increased transepithelial Cl− transport (Figure 5E). Amiloride-sensitive epithelial sodium channel and calcium-dependent chloride channel activities were also observed. Importantly, the response of CF cells to CFTR modulators did not decline in late-passage cells (Figure 5E).

Figure 5.

Differentiation of HBECs Expanded with the EpiX Medium

(A) HBECs from healthy and CF donors (n = 4) were expanded in the EpiX medium or conventional medium (BEGM). EpiX medium supported million-fold more expansion than BEGM (n = 2).

(B) HBECs from a healthy (UNC42I) or a CF donor were expanded in the EpiX medium for 30 PDs and differentiated at ALI for 21 days. The expression of basal cell markers (TP63, NGFR), multiciliated cell markers (CFTR, FOXJ1A), goblet cell marker (MUC5AC), club cell marker (CC10), or type II cell marker (SFTPD) were checked by qRT-PCR. Gene expression levels were depicted as relative to that of β-actin, which was set at 1. Data are represented as mean ± SD, n = 3.

(C) UNC42I cells (P8) were differentiated at ALI for 21 days. Paraffin sections were stained with H&E, or with an anti-acetylated tubulin antibody, or with alcian blue to show multiciliated cells and goblet cells.

(D) UNC42I cells (P7) were differentiated at ALI and treated with 1 ng/mL IL-13. Multiciliated cells and goblet cells were stained with anti-acetylated tubulin (in green) or anti-MUC5AC (in red) antibodies respectively. IL-13 led to goblet cell hyperplasia and a decrease of multiciliated cells.

(E) Early- (15 PDs) and late-passage (30 PDs) CF cells were differentiated at ALI for Ussing assays, using either the Vertex or the PneumaCult-ALI protocol. The activities of ENaC (ΔAmiloride, 30 μM amiloride), CFTR (ΔFsk peak, 10 μM forskolin; ΔVX-770 peak and ΔVX-770 plateau, 100 nM VX-770 and 3 μM VX-809) and CaCC (ΔUTP, 100 μM UTP) were measured. VX-809, a CFTR trafficking corrector; VX-770, a CFTR potentiator. The responses of mutant CFTR variants to the CFTR corrector (3 μM VX-809) and CFTR potentiator (100 nM VX-770) were similar between early and late passages. Data are represented as mean ± SD, n = 3.

Telomere Length and DNA Methylation Landscape Gradually Change over Long-Term In Vitro Expansion

As shown in Figure 1G, early-passage epithelial cells in the EpiX medium were small with a high nuclear/cytoplasm ratio and bright under microscope. The population doubled in <24 hr (Figures 1I-1L). Large cells accumulated slowly over passages and eventually the majority of the population appeared senescent or differentiated. Comparing the transcriptomes of early- versus late-passage HFKs revealed that many genes involved in the senescence-associated secretory phenotype (Coppé et al., 2010), including IL6, IL8, IGFBP3, IGFBP5, THBS1, SERPINE1, and TGFB1I1 increased significantly in late passage (Figure S9). We examined the telomeres length of EpiX-expanded HFKs by a qPCR method (Cawthon, 2002) and observed an inverse correlation between the T/S ratio and the PDs (Figure 6A). Interestingly, HFKs expanded by the CR method also exhibited telomere erosion at a similar rate (Liu et al., 2012), suggesting that the erosion of telomeres contributed to ultimate growth arrest in both conditions.

Figure 6.

Epigenetic Changes in HFKs Expanded with the EpiX Medium

(A) Telomere length (T/S ratio) gradually decreased in HFKs expanded in the EpiX medium or the CR method. Data are represented as mean ±SD, n = 3.

(B) The number of differentially methylated regions (DMRs) identified between different culture conditions. The numbers in upper right cells were DMRs with DNA methylation gain, and the numbers in lower left cells were DMRs with DNA methylation loss.

(C) PCA of DNA methylomes. DNA methylation levels of DMRs were used as input data. The first principal component (passage) explained most variances among the DNA methylomes.

(D) DNA methylation levels in HFKs gradually changed over successive passages in the EpiX medium. The methylation levels of DMRs at P12 was intermediate between P3 and P19. Also see Figure S10.

(E) The chromatin states of DMRs using 18 chromHMM expanded states in keratinocytes. Lots of DMRs were in regulatory regions such as promoters and enhancers. Also see Figure S10.

(F) Expression changes of genes whose promoters underwent DNA methylation changes. Expression changes were calculated as log2 of fold change between P3 and P19. The red dots indicated the genes whose expression was downregulated as expected from the gain of promoter methylation.

(G) Predicted cumulative population doublings (pcPDs) of HFKs based on the DNA methylation levels at six CpG sites correlated with the actual PDs at different passages.

We integrated methylated DNA immunoprecipitation sequencing (MeDIP-seq) and methylation-sensitive restriction enzyme sequencing (MRE-seq) (Li et al., 2015; Maunakea et al., 2010) methods to further examine genome-wide DNA methylation (DNAme) changes in the HFK samples that were used in the RNA-seq study (Figure 2A; Table S5). Average DNA methylation levels in each passage were around 70%–75% with a slight decrease at higher passages (Figure S10), consistent with previous reports (Wilson and Jones, 1983). Bimodal distribution revealed a lower number of highly methylated CpGs in late passages, yet the global DNA methylation over genic regions was stable across passages (Figure S10). In contrast to the transcriptomic shifts, PCA showed that passage number was the only factor that drove DNA methylation change (Figure 6C). Over 2,400 differentially methylated regions (DMRs) (1,940 DNA methylation gain and 479 DNA methylation loss) were identified along with the increase in passages, with no statistically significant DMRs identified between adjacent passages in different media (Figures 6B and S10). Interestingly, mid-passage (P12) HFKs exhibited intermediate DNA methylation levels in these DMRs, confirming the DNA methylation change was gradual across passages (Figures 6D and S10). The DNA methylation changes were maintained after withdrawal of EpiX medium (Figure S10), suggesting that these changes were driven mainly by the passage number but not media change.

We next investigated the potential impact of DNA methylation changes on gene expression. Taking advantage of the Human Reference Epigenome Map (Kundaje et al., 2015; Lowdon et al., 2014), we annotated the DMRs using histone modification data from reference foreskin keratinocytes (E057 and E058) and found that DMRs occurred more often than expected in promoters or enhancers (Figures 6E and S10; Table S5). The DNA methylation changes were weakly associated with gene expression change, and no significant DNA methylation change was found in or around genes involved in important signaling pathways (TGF-β, WNT) and cell cycle (Figure S10). Of all genes that exhibited promoter DNA methylation change, 18/126 down-regulated genes and 0/6 upregulated genes showed expected gene expression changes (Figures 6F and S10; Table S6). The majority (112/132, 85%) of genes with promoter DNA methylation changes remained at low expression levels (transcripts per million [TPM] < 10; Figure S10), matching previous report that DNA methylation changes accumulate predominantly at inactive gene promoters in long-term human cell cultures (Gordon et al., 2014). We used a model derived from long-term fibroblast culture (Koch et al., 2012) to predict the cumulative PDs (pcPDs) of the HFKs using DNA methylation levels at certain CpGs, and found a good correlation between the predicted and actual PDs (Figure 6G; Table S7). Overall, this suggested that the majority of DNA methylation changes accompanied with EpiX expansion were by-products of long-term culture and not responsible for gene expression changes.

DISCUSSION

As cell behaviors are influenced by diverse external signals, we postulated that screening small molecules could yield hits if we focused on desired outcome (e.g., continuous cell proliferation) even without a priori knowledge on the native stimuli. Through this approach, we found that TGF-β signaling inhibition and ROCK inhibition synergistically supported epithelial cell proliferation in the absence of feeder cells. TGF-β is a well-known cytostatic factor (Siegel and Massagué, 2003), presumably through the activation of p21Cip1 and p15INK4b (Bhowmick et al., 2003; Denicourt and Dowdy, 2003). Noteworthily, Rheinwald proposed that feeder cells protected epithelial cells from TGF-β growth inhibition by efficiently degraded TGF-β (Rollins et al., 1989). Attenuation of TGF-β signaling suppresses premature senescence in a p21Cip1-dependent manner (Lin et al., 2012) or antagonizes ATM-mediated growth-arrest response to genotoxic stress (Kirshner et al., 2006), and TGF-β inhibition rescues hematopoietic stem cell defect (Zhang et al., 2016). Nevertheless, inhibiting TGF-β signaling with small molecule alone cannot immortalize epithelial cells, but can assist TERT in doing so (Natarajan et al., 2006).

Modest expressions of TGF-β were detected in epithelial cells cultured in KSFM by ELISA and RNA-seq (Figures S3 and S6), suggesting that autocrine TGF-β might contribute to cell growth arrest. Although RNA-seq detected little change in TGF-β genes transcription in EpiX versus KSFM (Figure S6), ELISA showed that A83-01 and Y-27632 together significantly reduced active TGF-β protein produced by the epithelial cells, although either had little effect when used alone (Figure S3). The suppression on TGF-β protein may help explain the synergy between A83-01 and Y-27632 on promoting cell proliferation, but how the TGF-β protein is reduced in the absence of transcriptional changes warrants further investigation.

Recently, the Rajagopal group reported that “dual SMAD” inhibition enables long-term expansion of epithelial basal cells from airway, skin, and epididymis (Mou et al., 2016), while they also used Y-27632 and CHIR99021. They concluded that both TGF-β and BMP pathways needed to be suppressed by A83-01 and DMH-1, respectively. DMH-1 was included in our small-molecule collection but did not have a significant impact on its own nor with either TGF-β or ROCK inhibitors, suggesting that BMP pathway inhibition may not be critical at least under in vitro condition.

ROCKs influence a wide variety of signal pathways in eukaryotic cells and play pivotal roles in regulating actin cytoskeleton, cell polarity, microtubule dynamics, membrane transport pathways, and transcription factors activities (Etienne-Manneville and Hall, 2002). ROCK inhibitor has been widely used in mammalian cell cultures presumably for its ability to suppress anoikis (Watanabe et al., 2007). What has not been highlighted before is the synergy between TGF-β signaling inhibition and ROCK inhibition in promoting epithelial cell proliferation as we delineated in this study. We also found that PAK1 inhibitor (IPA-3) or Myosin II inhibitor (blebbistatin) had similar synergistic effect on promoting epithelial cell proliferation when used with A83-01. The Watt laboratory proposed that Hippo effector YAP nuclear translocation and its co-factor WBP2 as essential mediator to promote cell proliferation and drive clonal expansion when cytoskeleton and microtubule dynamics are pharmacologically modulated through the PAK1-ROCK-Myosin II axis (Walko et al., 2017). We found that Y-27632 but not A83-01 led to significant cytoskeleton reorganization and YAP nuclear translocation in HFKs cultured in KSFM with 1.5 mM CaCl2 (Figure S3), which provided an explanation for ROCK inhibitor as the primary hit in high-calcium F medium (Figures 1 and S1).

High extracellular [Ca2+] induces epithelial cell differentiation through promoting intercellular interaction (Martin et al., 1991). When [Ca2+] is below 0.3 mM, NOTCH-1 is constitutively active and allows cell-autonomous signaling in the absence of reciprocal cell-cell interaction (Dalrymple et al., 2005; Rand et al., 2000). Interestingly, the low [Ca2+] in KSFM resulted in cytoskeletal re-organization and YAP nuclear translocation similar to adding Y-27632 in the presence of high [Ca2+] (Figure S3). Extracellular Ca2+ also gates cadherin-cadherin homotypic interaction and influences β-catenin nuclear trafficking and WNT target genes expression (Nusse and Clevers, 2017). Such miscellaneous effects of low extracellular [Ca2+] probably overshadow any activity imparted by the GSK3 inhibitor CHIR99021, which showed little effect in the EpiX medium.

Importantly, our study revealed that EpiX-expanded epithelial cells maintained remarkable genome integrity. The cells retained normal diploid and had extremely low SNV rate comparable to germline cells (Rahbari et al., 2016). No mutations occur in oncogenes or tumor suppressor genes. On the contrary, many PSC lines are reported to acquire genomic variations including TP53 mutations with successive passages (Merkle et al., 2017). Current PSC differentiation protocols produce epithelial progenitors that only propagate for a few passages in convention medium (Firth et al., 2014; Umegaki-Arao et al., 2014). It is tempting to use EpiX medium to expand PSC-derived epithelial progenitors, which could reduce the production cost and assure the genomic quality of end cell products.

Transcriptomic study of HFKs cultured in EpiX versus KSFM revealed significant downregulation of genes involved in senescence, cell-cell interaction, interferon signaling, extracellular matrix organization, and stress responses, together with upregulation of genes involved in various metabolic processes and cell cycle (Figure 2F). Importantly, the impact was not permanent and reversed upon the withdrawal of EpiX medium (Figure 2E), allowing the cells to differentiate along their tissue lineages. Genome-wide DNA methylation study revealed gradual accumulation of DMRs that were mainly associated with successive passages. These DMRs were weakly associated with gene expression change, and occurred mostly in the promoters or enhancers of low-expression genes, confirming similar observation in long-term culture of human cell lines (Gordon et al., 2014). The DNA methylation changes were maintained after withdrawal of EpiX medium (Figure S10), suggesting that these DMRs were mainly driven by long-term in vitro expansion.

Cell proliferation eventually ground to a halt in the EpiX medium, presumably due to telomeres erosion that lead to ATM/TP53-dependent DNA-damage response (Lazzerini-Denchi and Sfeir, 2016). Conversely, this might be a major contributor of the remarkable genome integrity due to prompt elimination of the cells with DNA damages via apoptosis. The accumulation of senescent cells in late-passage population is an important quality concern for ex vivo cell manufacturing, as they may impair healthy cells’ functions via destructive paracrine effect (Campisi, 2005; Rodieretal., 2009). Recent advances in the identification of senolytics that can postpone senescence or eliminate senescent fibroblasts, such as rapamycin (Iglesias-Bartolome et al., 2012), IL-1Ra (Uekawa et al., 2004), ABT-263 (Chang et al., 2016), and FOXO4-DRI (Baar et al., 2017), might be useful to further extend the expansion of epithelial cells.

In summary, the EpiX medium supports over one trillion-fold expansion of epithelial cells from diverse tissues in the absence of feeder cells. The EpiX technology provides unique solutions to unleash the potential of tissue-resident epithelial stem and progenitor cells for cell therapy and regenerative medicine.

STAR★METHODS

Detailed methods are provided in the online version of this paper and include the following:

CONTACT FOR REAGENT AND RESOURCE SHARING

Further information and requests for resources and reagents should be directed to and will be fulfilled by the Lead Contact, Chengkang Zhang (ck.zhang@propagenix.com)

EXPERIMENTAL MODEL AND SUBJECT DETAILS

Primary epithelial cells were purchased from Lonza or ThermoFisher or obtained from academic tissue bank (UNC CF Center). Fresh human foreskin tissue was purchased from Zenbio for keratinocyte isolation. Fresh human nasal brushing samples were collected from CF patients at UNC CF Clinic following approved IRB protocol. The epithelial cells were cultured in the CR condition, commercial media (Lonza BEGM, Lonza PrGM, Lonza MEGM, GIBCO KSFM), or the EpiX medium. The EpiX medium was KSFM supplemented with 1 μM A83–01, 5 μM Y-27632 and 3 μM isoproterenol. The epithelial cells were cultured in collagen-coated vessels (Corning) until they reach 80%–90% confluence, when they were passaged with Trypsin-EDTA and re-seeded into new vessels.

For in vivo animal studies, female nude mice (6-8 weeks old) were used and each group had 6 animals. Body weight and tumor measurements were recorded twice weekly for 12 weeks. The study was conducted in compliance with the IACUC protocols established at Noble Life Sciences.

METHOD DETAILS

Lentivirus transduction and live-imaging cell growth analysis with IncuCyte

Epithelial cells cultured using the CR method were seeded into 24-well plate and cultured without feeder cells. The cells were transduced with IncuCyte® NucLight Red (nRFP) Lentivirus Reagent (Essen BioScience) and stable clones were established by puromycin selection. Late-passage nRFP-expressing cells were seeded into collagen-coated 96-well or 384-well plate (Corning). A collection of small molecules (Tocriscreen Stem Cell Toolbox, Tocris; Table S1) were diluted in culture medium to the desired concentrations and added to the culture. The cells were cultured in IncuCyte ZOOM (Essen BioScience) for live-imaging analysis following the manufacturer’s instructions for 5–7 days. Small molecules with positive effects were purchased individually from Sigma-Aldrich or Sellekchem for validation. Stocks of the chemicals were prepared by dissolving in DMSO to 10 mM and added to culture media to desired final concentrations.

Keratinocytes isolation and expansion

Fresh human skin purchased from ZenBio Inc. was cut into small pieces and placed in Dispase solution (Corning) at 4C overnight. The next day, the epidermis was separated from dermis with forceps and digested in Trypsin-EDTA for 15 minutes. The cell suspension was filtered through a 40 μm strainer (BD Bioscience) and cultured in EpiX medium.

Cultured cell immunofluorescence

Cells cultured on collagen-coated CultureSlides (Corning) were fixed in 4% PFA for 15 min at room temperature, washed 3 times in PBST (PBS+0.2% Triton X-100) (5 min/wash) and incubated with the primary antibodies for 2 hours at room temperature or 4C overnight in PBS +1% normal goat serum. Following incubation, the cells were rinsed 3 times in PBST and incubated with secondary antibodies at room temperature for 1–2 hours. After rinsing 3 times in PBST, the nuclei were stained with ProLong® Gold with DAPI (ThermoFisher) and imaged with a fluorescence microscope (EVOS-FL, ThermoFisher).

Telomere length measurement

Genomic DNA was extracted from the cells using the Quick-DNA Miniprep Plus Kit (Zymo Research) and quantitated on a NanoDrop 2000. Telomere length measurements by quantitative PCR was performed by the Risques Lab at University of Washington.

Quantitative RT-PCR

Total RNA was isolated using the TRIzol Plus RNA Purification Kit (ThermoFisher) and PureLink RNA Mini Kit (ThermoFisher). Total RNAs from human lung, small intestine and mammary gland were purchased from Clontech and used as control. RT2 Profiler PCR Array Human Cellular Senescence (QIAGEN) was used to analyze the expression of genes involved in cellular senescence using twenty nanogram total RNA for each reaction. The primers (forward: 5′-TGACACCTCACCTCACCCAC-3′; reverse: 5′-CACTGTCTTCCGCAAGTTCAC-3′; Taqman probe: 5′-ACCCTGGTCCGAGGTGTCCCTGAG-3′) were used to determine TERT expression with the TaqMan RNA-to-CT 1-Step Kit (ThermoFisher), using one hundred nanogram total RNA per reaction.

Karyotype and Whole genome sequencing

Live cells were seeded in flasks and sent to Cell Line Genetics for G-Band Karyotyping service. Genomic DNA was extracted from the cells using the Quick-DNA Miniprep Plus Kit (Zymo Research) and quantitated on a NanoDrop 2000. Whole genome library preparation and DNA sequencing on Illumina HiSeq X and data analysis were performed by Novogene Inc.

Differentiation of HBEC at ALI

HBEC were seeded onto polyester Transwell membranes (Corning) at a density of 400,000 cells/cm2 in EpiX+1.5 mM CaCl2. After the cells reached confluence on the insert, the medium was replaced with Pneumacult-ALI medium (STEMCELL Technologies) only in the lower chamber to initiate air-liquid interface culture. The medium was changed every 2–3 days for 21–28 days until differentiation was well established. Ciliogenesis was monitored by inverted-phase microscopy. The membranes were fixed in 4% PFA at room temperature for 10 min, followed by washing and permeabilization in PBST, for immunofluorescence staining. Paraffin section, H&E and alcian blue staining were performed by VitroVivo Biotech LLC.

Ion channels activity assays

Ussing assays were performed by ChanTest, a subsidiary of Charles River Laboratories. Frozen cells were sent to ChanTest for plating and differentiation at CRL using the Vertex (Neuberger et al., 2011) or Pneumacult-ALI medium. Short circuit current (ISC) was evaluated after 32 and 36 days in ALI culture.

Differentiation of keratinocytes at ALI and production of cultured epidermal grafts

HFK were seeded onto polycarbonate Transwell membranes (Corning) at a density of 400,000 cells/cm2 in EpiX+1.5 mM CaCl2. After 2–3 days, the medium was removed from top chamber and only added in the lower chamber to initiate ALI differentiation. The medium was changed every 2–3 days for 14 days. To make cultured epidermal graft, HFK were cultured in collage-coated T-75 flask in the EpiX medium until they approached confluence, then the medium was switched to EpiX+1.5 mM CaCl2 and cultured for another 7 days. The whole epidermal sheet was released using Dispase solution at 37C for 30 min.

In vivo tumorigenicity and differentiation

HFK were expanded in EpiX or differentiated for 7 days in EpiX+1.5 mM CaCl2. The cells were harvested with Trypsin-EDTA and resuspended in HBSS for in vivo tumorigenicity assays in female nude mice, performed by Nobel Life Sciences as contract service in compliance with established IACUC protocols. Briefly 1.0×107 cells were administered subcutaneously in a volume of 0.1 mL in Matrigel on the flanks of 6–8 weeks old mice (n = 6 for each group). Following cell injection, body weight and tumor measurements were recorded twice weekly for 12 weeks. Tissue samples were harvested at specified time points and fixed in formalin for paraffin sections.

RNA-Seq library preparation

Half million HFK per replicate were collected and stored at −80°C in 1 mL of RNAlater solution until use. For late-passage HFK cultured in KSFM (P5 and P20), 66,000-150,000 cells per replicate were used as the cells quickly ceased expansion in KSFM. The mRNA was extracted directly from the cells by using Dynabeads mRNA DIRECT Purification Kit according to the manufacturer’s instructions. To avoid any genomic DNA and rRNA contamination, eluted mRNA was treated with DNase (TURBO DNA-free Kit) and bound again to the same Dynabeads used for the original isolation, following the manufacturer’s instructions. RNA-Seq libraries were generated using ScriptSeq RNA-Seq Library Preparation Kit (Illumina) according to the manufacturer’s instructions.

RNA-Seq data processing

RNA-Seq libraries were sequenced on the Illumina NextSeq 500 platform. The paired-end sequenced reads were aligned to the reference genome hg19 and transcriptome (gencode v19) using STAR (Dobin et al., 2013) v.2.5.3a with the following parameters:–outFilterType BySJout–outFilterMultimapNmax 20–outFilterMismatchNmax 999–outFilterMismatchNoverReadLmax 0.04–alignIntronMin 20–alignIntronMax 1000000–alignMatesGapMax 1000000–alignSJoverhangMin 8–alignSJDBoverhangMin 1. The total number of reads overlapping each gene were counted using featureCounts (Liao et al., 2014) (Subread v.1.5.3) with gencode v19 gtf file and the following parameters: -O -s 1–primary -p. Transcripts per million (TPM) were calculated from the read counts obtained from featureCounts.

Differential gene expression analysis

Differential gene expression was performed by using DESeq2 bioconductor package (Love et al., 2014) v.1.14.1, following standard workflow suggested by the package. The matrix of genes and read counts generated by featureCounts were used as count matrix input. The dataset was pre-filtered by removing genes with no or only one count across all samples. For PCA analysis and heatmap, regularized-logarithm transformation (rlog) was performed for count data to stabilize the variance across the mean. Differentially expressed genes were identified by comparing two different media conditions, control KSFM (P2, P4, P5, and P20) and EpiX (P3, P12, and P19). Only genes with log2FC > 1 and adj.p <0.01 were considered significant. The enrichment of gene ontology terms and pathways for DEGs was analyzed by using Metascape (Tripathi et al., 2015).

MeDIP-seq and MRE-seq library preparation

Half million HFK per replicate were collected and stored as a pellet at −80°C until use. For late-passage HFK cultured in KSFM (P5 and P20), 66,000-150,000 cells per replicate were collected due to difficulty of culture expansion. The genomic DNA was extracted by incubating the cells in genomic DNA extraction buffer (50 mM Tris, 1 mM EDTA, 0.5% SDS, 1 mg/mL Proteinase K) followed by phenol-chloroform extraction. MeDIP-seq and MRE-seq libraries were generated as described previously (Li et al., 2015; Maunakea et al., 2010), with minor modifications. For MeDIP-seq, 100 ng of genomic DNA was sonicated to a fragment size of 100-500 bp, end processed and ligated to paired-end adapters. The DNA was then denatured and immunoprecipitated using 100 ng of mouse monoclonal anti-methylcytidine antibody in 400 μL of immunoprecipitation buffer (10 μM sodium phosphate, pH 7.0, 140 mM NaCl and 0.05% Triton X-100) overnight at 4°C. Antibody/DNA complexes were isolated by addition of 0.1 μL of rabbit anti-mouse IgG secondary antibody (2.0 mg/mL, Jackson Immunoresearch) and 20 μL protein A/G agarose beads (Pierce Biotechnology) for 2 h at 4°C. Beads were washed ten times with immunoprecipitation buffer and then DNA was eluted in TE buffer with 0.25% SDS and 0.25 mg/mL of proteinase K for 2 h at 50°C. DNA was then purified with MinElute PCR Purification kit (QIAGEN). DNA was amplified by 17 cycles of PCR with the standard Illumina index primers and size selected (150–500 bp) by Agencourt AMPure XP beads (Beckman Coulter). For MRE-seq, five parallel digests (HpaII, Hin6I, SsiI, BstUI and HpyCH4IV; New England Biolabs) were performed, each with 20 ng of genomic DNA. The digested DNA was size selected using Agencourt AMPure XP beads (Beckman Coulter), end processed and ligated to adapters. Then DNA was amplified by 18 cycles of PCR and size selected (150–500 bp) by Agencourt AMPure XP beads (Beckman Coulter).

DNA methylome data processing

MeDIP-seq and MRE-seq libraries were sequenced on the Illumina NextSeq 500 platform. The sequenced reads were adaptor-trimmed by using cutadapt (Martin, 2011) v.1.9 paired-end mode with the parameters -q 10 -m 20. Trimmed reads were aligned to the hg19 genome assembly using BWA-MEM (Li, 2013) v.0.7.10 with the default parameters. The aligned MeDIP-seq reads were further processed using methylQA (Li et al., 2015) v.0.1.6 medip mode with the default parameters. The aligned MRE-seq reads were processed using methylQA (Li et al., 2015) v.0.1.6 mre mode with the parameter -c 4. Methylation levels at single CpG resolution were estimated by integrating MeDIP-seq and MRE-seq data using methylCRF with default parameters as described previously (Stevens et al., 2013).

Identification of differentially methylated regions

Differentially methylated regions between two culture conditions were identified by using methylMnM package (Zhang et al., 2013) with the default parameters. Briefly, the coverage of MeDIP and MRE sequencing data and genomic CpG information were calculated in each 500-bp genomic bin. The CpGs in the human blacklisted genomic regions (ENCODE Project Consortium, 2012) and the mitochondrial genome were excluded from the analysis. DMRs with a q-value (false discovery rate) <1 × 10−5 were selected for each pairwise comparison (Figure S10A). Highly reproducible set of DMRs were identified by intersecting all four pairwise comparisons between two different conditions with two biological replicates and were used for further analysis (Figure 6B).

DNA Methylation analyses

Average DNA methylation level of each DMR was calculated by using CpG methylation levels estimated by methylCRF. Principal component analysis was performed on the average DNA methylation levels of all union set of 2,419 DMRs with high reproducibility. The chromatin states of each DMR was determined by intersecting either core 15 chromHMM states or expanded 18 states of fetal keratinocytes (ID# 057 and 058) from Roadmap Epigenomics Consortium data.

Prediction of cumulative population doublings

To predict cPDs, genomic coordinates of 6 CpG sites were determined by using sequences surrounding the CpG sites. For cg03891191, two genomic locations were found in hg19 genome assembly and average of methylCRF values of two locations were used. Cumulative population doublings were calculated by using the following equation (Koch and Wagner, 2013): pcPD = 45.89 + (23.63 × cg02332525) + (31.61 × cg17453778) + (−53.70 × cg03891191) + (14.86 × cg01459453) + (−23.94 × cg01999333) + (−10.34 × cg16431978).

QUANTIFICATION AND STATISTICAL ANALYSIS

IncuCyte Zoom was used for live imaging and automatic cell counts in 96-well or 384-well plates. Automatic cell count was facilitated by counting nuclear-localized RFP. For multi-well plate assay, each condition had 3–4 replicates, and all error bars correspond to SD. GraphPad Prism 7.0 was used for non-parametric tests. Unpaired data was compared using Kruskal-Wallis test. Differences were considered significant if p < 0.05. All collected data were included for the quantification and the statistical analysis.

In the transcriptome study, regularized-logarithm transformation (rlog) was performed for count data to stabilize the variance across the mean for PCA analysis and heatmap. Differentially expressed genes were identified by comparing two different media conditions, and only genes with log2FC > 1 and adj.p < 0.01 were considered significant. In the DNA methylation study, DMRs with a q-value (false discovery rate) <1 × 10−5 were selected for each pairwise comparison (Figure S10A). Highly reproducible set of DMRs were identified by intersecting all four pairwise comparisons between two different conditions with two biological replicates and were used for further analysis (Figure 6B).

DATA AND SOFTWARE AVAILABILITY

The accession number for the data reported in this paper is Gene Expression Omnibus (GEO): GSE103759, and contains the subseries GEO: GSE103756 (RNA-seq), GSE103757 (MeDIP-seq), and GSE103758 (MRE-seq). All data generated in this study have been visualized in WashU Human Epigenome Browser (Zhou et al., 2011) and is publicly available in the following url: http://epigenomegateway.wustl.edu/browser/?genome=hg19&datahub=http://wangftp.wustl.edu/~hlee/EpiX/EpiX_hg19.json.

KEY RESOURCES TABLE

REAGENT or RESOURCE	SOURCE	IDENTIFIER
Antibodies
Rabbit polyclonal anti-MUC5AC	Santa Cruz Biotechnology	Cat# sc-20118; RRID: AB_2146854
Mouse monoclonal anti-Acetylated Tubulin	Sigma-Aldrich	Cat# T7451; RRID: AB_609894
Rabbit polyclonal anti-Ki67	Abcam	Cat# ab15580; RRID: AB_443209
Mouse monoclonal anti-KRT14	Thermo Fisher Scientific	Cat# MA5-11599; RRID: AB_10982092
Rabbit monoclonal anti-YAP	Cell Signaling Technology	Cat# 14074; RRID: AB_2650491
Alexa Fluor 594 Phalloidin	Thermo Fisher Scientific	Cat# A12381; RRID: AB_2315633
Rabbit monoclonal anti-P63	Cell Signaling Technology	Cat# 13109; RRID: AB_2637091
ZO-1 Monoclonal Antibody (ZO1-1A12), Alexa Fluor 488	Thermo Fisher Scientific	Cat# 339188; RRID: AB_2532187
Occludin Monoclonal Antibody (OC-3F10), FITC	Thermo Fisher Scientific	Cat# 33-1511; RRID: AB_2533102
Bacterial and Virus Strains
NucLight Red Lentivirus (EF-1 Alpha Promoter, Puromycin selection)	Essen Biosciences	Cat# 4476
Biological Samples
Human foreskin tissue, fresh	Zenbio	Cat# T-FS
Human nasal brushing samples from CF patients	UNC CF Clinics	N/A
Chemicals, Peptides, and Recombinant Proteins
Tocriscreen Stem Cell Toolbox	Tocris	Cat# 5060
A83-01	Sigma-Aldrich	Cat# SML0788
Y-27632 dihydrochloride	Enzo Life Sciences	Cat# ALX-270-333
Isoprenaline hydrochloride	Sigma-Aldrich	Cat# I5627-5G
(−)-Blebbistatin	Selleck Chemicals	Cat# S7099
IPA-3	Tocris	Cat# 3622
GSK-429286	Sigma-Aldrich	Cat# SML0231
Critical Commercial Assays
Telomere length measurements by quantitative PCR service	Risques Lab at University of Washington	http://depts.washington.edu/risques/
G-Band Karyotyping service	Cell Line Genetics	N/A
Whole genome sequencing service	Novogene	N/A
Ussing assays service	ChanTest, a subsidiary of Charles River Laboratories	N/A
Histology service (Paraffin section, H&E and alcian blue staining)	VitroVivo Biotech LLC	N/A
TGF-β ELISA service	University of Maryland Cytokine Core Lab	N/A
RT2 Profiler PCR Array Human Cellular Senescence	QIAGEN	Cat# PAHS-050ZA-6
In vivo tumorigenicity service	Noble Life Sciences	N/A
Deposited Data
Raw and analyzed RNA-seq, MeDIP-seq and MRE-seq data	This paper	GSE103759 containing the Subseries GSE103756 (RNA-seq), GSE103757 (MeDIP-seq) and GSE103758 (MRE-seq)
Analyzed whole genome sequencing data	This paper	Figure 3; Table S4
Experimental Models: Cell Lines
Human bronchial epithelial cells, P1, from healthy and CF donors	UNC CF Center Tissue Procurement and Cell Culture Core	N/A
NHBE-Bronchial Epi Cells for B-ALI	Lonza	Cat# CC-2540S
CF-DHBE - Diseased Bronchial Epi. Cells (CF)	Lonza	Cat# 196979
Human Epidermal Keratinocytes, adult (HEKa)	Thermo Fisher Scientific	Cat# C0055C
Human Epidermal Keratinocytes, neonatal (HEKn)	Thermo Fisher Scientific	Cat# C0015C
HMEC-Human Mammary Epithelial Cells	Lonza	Cat# CC-2551
LNCap Clone FGC Cell Line human	Sigma-Aldrich	Cat# D-073-1ML
Human normal prostate epithelial cells	Georgetown University	N/A
Experimental Models: Organisms/Strains
Nude mice (Noble Life Sciences)	Jackson Lab	N/A
NSG mice (Noble Life Sciences)	Jackson Lab	N/A
Oligonucleotides
qRT-PCR primers for TERT (forward primer, 5′-TGACACCTCACCTCACCCAC-3′, reverse primer, 5′-CACTGTCTTCCGCAAGTTCAC-3′ and Taqman probe (5′-ACCCTGGTCCGAGGTGTCCCTGAG-3′)	This paper	N/A
Software and Algorithms
STAR, v.2.5.3a	Dobin et al., 2013	https://github.com/alexdobin/STAR
featureCounts (Subread v.1.5.3)	Liao et al., 2014	http://bioinf.wehi.edu.au/featureCounts/
DESeq2 bioconductor package v.1.14.1	Love et al., 2014	https://bioconductor.org/packages/release/bioc/html/DESeq2.html
Metascape	Tripathi et al., 2015	http://metascape.org
BWA-MEM v.0.7.10	Li, 2013	http://bio-bwa.sourceforge.net/
methylQA	Li etal., 2015	http://methylqa.sourceforge.net/
methylCRF	Stevens et al., 2013	http://methylcrf.wustl.edu/
methylMnM	Zhang et al., 2013	http://bioconductor.org/packages/release/bioc/html/methylMnM.html
Predict cumulative population doublings using DNA methylation of 6 CpG sites	Koch and Wagner, 2013	Figure 6 and Table S7
Other
IncuCyte ZOOM System	Essen Biosciences	N/A
Data visualization in WashU Human Epigenome Browser	This paper	http://epigenomegateway.wustl.edu/browser/?genome=hg19&datahub=http://wangftp.wustl.edu/~hlee/EpiX/EpiX_hg19.json