Prospective isolation of NKX2-1–expressing human lung progenitors derived from pluripotent stem cells

Little is known about the early stages of human lung development, preventing an understanding of whether successful healing from adult lung injury involves recapitulation of embryonic mechanisms and limiting approaches for generating lung progenitors from pluripotent cells in vitro. Inbred mouse models have begun to define the mechanisms regulating lung specification and patterning, but how this breadth of work applies to human lung development is unknown. Current claims suggest that all cells of the postnatal mammalian lung epithelium derive from embryonic NKX2-1⁺ progenitors; however, scant literature exists formally testing this hypothesis either in mice or humans. Support for this paradigm derives mainly from the observation that Nkx2-1 is the first gene locus known to be activated in cells of the endodermal lung primordium (1, 2). NKX2-1–null mutant mice have hypoplastic lungs that fail to mature, and human children with NKX2-1 mutations develop respiratory insufficiency, hypothyroidism, and neurological impairment (3), but these observations do not necessarily indicate that all lung epithelial cells derive via NKX2-1⁺ progenitor intermediates. Since lung lineage specification is thought to occur in relatively few endodermal cells during a narrow developmental time period in vivo, it has been difficult to gain access to these cells in human embryos or to follow their cell fate decisions in real time. Hence, we sought to interrogate the earliest moments of human lung lineage specification by engineering an in vitro system that would allow the isolation and differentiation of pure populations of NKX2-1⁺ putative human lung progenitors. Given the known capacity of mouse pluripotent stem cells (PSCs) to form all cell types, including lung lineages, after transfer into mouse blastocyst embryos and the known broad differentiation repertoire of human PSCs in vitro, we based this system on the in vitro differentiation of PSCs.

Initial published attempts at deriving lung epithelium from PSCs relied on the presence of drug-resistance genes or used incompletely defined media (4–6), resulting in inefficient induction of selected lung markers. Subsequently, a number of groups, including our own, had more success by broadly attempting to recapitulate the key milestones of embryonic lung development in vitro through the exogenous addition of sequential combinations of growth factors (7–11). By differentiating iPSCs into definitive endoderm, patterning this endoderm via inhibition of TGF-β and BMP signaling (7), and then adding various combinations of Wnts, FGFs, BMPs, and retinoic acid (Ra) these groups demonstrated the in vitro derivation of cultures expressing a broad array of lung epithelial markers. However, the characterization of the cells derived at different stages of these protocols suggested that heterogeneous cell types were present (12, 13). The most recently published directed-differentiation protocols describe variable efficiencies of induction of NKX2-1⁺ cells from embryonic stem cells (ESCs) or induced PSCs (iPSCs) ranging from approximately 36% to 86% (9, 10, 14). Such heterogeneity limits the utility of these cultures for downstream applications and has caused uncertainty as to whether subsequent lung lineages derive directly from these early endodermal NKX2-1⁺ precursors. To derive more mature lung cell types from iPSCs, some groups have employed prolonged in vitro cultures, murine kidney capsule or subcutaneous transplantations, or coculture with lung mesenchyme (LgM) (9–11, 14–16). These results suggest that at some point during differentiation of iPSCs, progenitor intermediates with competence for forming mature lung cells likely emerge. However, it has not previously been possible to isolate these cells for characterization or to properly test their differentiation repertoire.

To perform a detailed profiling of candidate human lung progenitors we targeted the NKX2-1 human locus with a fluorescent reporter, enabling isolation of the earliest identifiable putative lung-lineage-committed cells derived from PSCs. In addition, these same reporter iPSC lines facilitate the derivation and purification of alternate developing human progenitors that express NKX2-1, such as endodermal thyroid–like and ectodermal forebrain–like lineages. After directed differentiation of these PSCs in defined media designed to promote lung rather than thyroid or forebrain development, we demonstrate that the NKX2-1⁺ endodermal population is highly enriched in undifferentiated (primordial) progenitors that are competent at expressing a broad repertoire of lung epithelial marker genes, supporting the paradigm that the human lung epithelium derives from embryonic NKX2-1⁺ progenitors. We provide population-based as well as single-cell global transcriptomic profiles that define developmental stage–specific gene signatures of iPSC-derived lung progenitors. These signatures suggest that an evolutionarily conserved lung developmental program exists in both mice and humans. Furthermore, these signatures reveal that NKX2-1⁺ human lung progenitors can be prospectively isolated from patient-specific iPSCs based on the cell surface phenotype CD47^hiCD26^lo. Ultimately, the detailed characterization and purification of human lung progenitors presented here should provide access to an inexhaustible supply of these cells for disease modeling, future cell-based therapies, and basic developmental studies.

An NKX2-1^GFP reporter enables purification of human lung, thyroid, and forebrain lineages. In order to generate a tool for the identification and purification of candidate human lung progenitors, we used gene editing technologies to target an enhanced green fluorescence reporter (GFP) gene to the endogenous human NKX2-1 locus in multiple human PSC lines. Prior reports of targeting GFP to the NKX2-1 locus in human PSCs for the derivation of forebrain lineages resulted in NKX2-1 haploinsufficiency (17). Hence, pursuing a strategy designed to retain intact expression of targeted loci without haploinsufficiency (Figure 1A and Supplemental Figure 1A; supplemental material available online with this article; https://doi.org/10.1172/JCI89950DS1), we targeted an exon3-2A-GFP cassette to the second intron of NKX2-1 using either transcription activator–like effector nucleases (TALENs) or CRISPR-Cas9 tools deployed in ESCs (H9) or in our previously published iPSC lines: cystic fibrosis patient–specific C17 iPSCs (18) and normal BU3 iPSCs (Figure 1 and Supplemental Figure 1) (19). The resulting NKX2-1^GFP reporter PSC clones (hereafter H9NKX2-1^GFP, C17NKX2-1^GFP, and BU3NKX2-1^GFP) demonstrated successful mono- and bi-allelic integration of the donor template by PCR (Supplemental Figure 1B). For further profiling we selected 1 homozygous targeted clone for each of the TALENs-targeted lines (H9 and C17) as well as 1 homozygous CRISPR-Cas9–targeted BU3 clone (Supplemental Figure 1, B and C).

Figure 1

Purification of human NKX2-1⁺ lineages derived from ESCs/iPSCs using NKX2-1^GFP reporters. (A) Gene editing strategy based on TALENs or CRISPR technology to target a GFP reporter to the human NKX2-1 locus to engineer NKX2-1^GFP iPSC/ESC lines. See also Supplemental Figure 1. (B) Schematic overview of in vitro directed differentiation of ESCs/iPSCs into NKX2-1⁺ lineages: endodermal lung or thyroid epithelia versus ectodermal forebrain. (C) Representative timeline of GFP expression measured by flow cytometry during lung-directed differentiation (C17). (D) Immunostaining of day 15 lung-directed differentiation for NKX2-1, GFP (anti-GFP), and EPCAM. Scale bar (left panels): 100 μm. Right panel is a zoom in (white dashed-line box). Scale bar: 25 μm (C17). (E) Phase contrast and fluorescence microscopy of C17 iPSC-derived cells generated in the lung (day 15), thyroid (day 27), and forebrain (day 10) protocols. See also Supplemental Figure 2. (F) Flow cytometric analysis of lung (day 15), thyroid (day 17), and forebrain (day 14) protocols (C17). (G) Fold change, compared with day 0, of mRNA expression of sorted NKX2-1^GFP+ and NKX2-1^GFP– cells from those time points by RT-qPCR; quantified as 2^(–ΔΔCT), n = 3 (C17).

To differentiate each targeted PSC line, we tested protocols previously developed by us and others for the in vitro directed differentiation of human PSCs into the 3 lineages known to express NKX2-1: neural/forebrain (20, 21), lung (7–9), or thyroid (19) (Figure 1B). Protocols for lung and thyroid both required generating anterior foregut–like endoderm followed by the addition of Chir99021 (CHIR), BMP4, KGF, FGF10, and Ra for lung versus BMP4 and FGF2 for thyroid. Consistent with prior publications (9), the percentage of NKX2-1⁺ cells on day 15 of the lung-directed differentiation was typically 41% ± 21% (mean ± SD) when using the ESC line RUES2 (data not shown). However, when other cells lines including H9 and C17 were differentiated the percentage of NKX2-1–expressing cells was initially less than 1% (data not shown). Thus, we developed a methodology for the optimization of lung differentiation for each cell line by altering the duration of endoderm induction, the cell density of replated endoderm, and the duration of TGF-β/BMP inhibition in order to augment the percentage of NKX2-1⁺ cells emerging by day 15 (Supplemental Figure 1H and Supplemental Figure 2, A and B). Employing the new lung differentiation protocol optimized for each clone, GFP expression was first detected between days 8 and 10 of differentiation and the percentage of GFP⁺ cells peaked between days 12 and 16, with an efficiency of 25.6% ± 9.5% for C17NKX2-1^GFP (Figure 1, C and F), 29.5% ± 4.4% for BU3NKX2-1^GFP, and 19.8% ± 13.2% for H9NKX2-1^GFP (mean ± SD, data not shown). Day 15 of differentiation was selected for further characterization of GFP⁺ cells in the lung protocol. Immunostaining for cytoplasmic GFP and nuclear NKX2-1 protein indicated faithful and specific expression of the GFP reporter in NKX2-1⁺ cells (Figure 1D). We confirmed that NKX2-1 protein levels were not significantly perturbed by our targeting strategy by comparing homozygously targeted BU3 NKX2-1^GFP iPSCs with nontargeted parental control iPSCs (Supplemental Figure 1G). Sorting GFP⁺ cells from each of the 3 differentiation protocols enriched for cells expressing NKX2-1 mRNA (Figure 1, F and G). GFP⁺ cells from thyroid and forebrain protocols were selectively enriched in expression of early lineage markers: PAX8, HHEX, and FOXE1 for thyroid and OTX2, OTX1, SOX1, PAX6, and SIX3 for forebrain (Figure 1G, Supplemental Figure 2, C–E and data not shown). The NKX2-1^GFP+ cells isolated on day 15 in the lung protocol were NKX2-1⁺ PAX8^– (Supplemental Figure 2E) and importantly lacked detectable expression of markers of lung epithelial differentiation (e.g., SFTPC or SCGB1A1; data not shown), raising the possibility that they may have undergone lung lineage specification but were still undifferentiated or primordial. Lung-specific progenitor markers or transcriptomic signatures are not known at this primordial developmental stage. Hence, to test the hypothesis that day 15 NKX2-1^GFP+ cells represented lung epithelial progenitors, we next interrogated their competence for subsequently expressing markers of more differentiated lung epithelium.

Purified iPSC-derived NKX2-1^GFP+ cells exhibit lung progenitor potential. Employing our human lung–directed differentiation protocol and the C17NKX2-1^GFP iPSC line, we tested whether differentiated lung epithelial cells derive directly from NKX2-1⁺ progenitors. We sorted NKX2-1^GFP+ versus NKX2-1^GFP– cells on day 15 and plated each population (vs. unsorted controls) in 3D Matrigel in serum-free media supplemented with factors we and others have previously shown to support lung epithelial differentiation (8, 9): CHIR, KGF, and FGF10 for 7 days followed by the addition of dexamethasone, cAMP, and IBMX until day 36 (Figure 2A). We observed the outgrowth of proliferating cell aggregates over the next 2 to 3 weeks (hereafter referred to as organoids) (Figure 2A). Unsorted day 15 cells plated as clumps gave rise to lobular organoids with GFP⁺ and GFP^– areas on day 36 (Figure 2A), and GFP⁺ cells could be followed in real time in these unsorted cultures by time-lapse photography (Figure 2B and Supplemental Video 1). In contrast, when unsorted cells were plated as single-cell suspensions, simpler spherical organoids formed (data not shown). Immunostaining of the unsorted organoids demonstrated areas of monolayered NKX2-1⁺EPCAM⁺ epithelium surrounding an inner lumen but also areas and organoids that were NKX2-1^– (Figure 2A). Sorted day 15 GFP⁺ cells gave rise to GFP⁺ aggregates in 3D culture, with 50.1% ± 17.6% (mean ± SD; n = 6 runs) of the progeny remaining GFP⁺ by flow cytometry on day 36 (GFP⁺ outgrowth) (Figure 2C). Sorted day 15 GFP^– cells remained GFP^– on day 36 (99.8% ± 0.2%) (GFP^– outgrowth) (Figure 2C). The GFP⁺ outgrowth formed predominantly NKX2-1⁺EPCAM⁺ spheroids, whereas the GFP^– outgrowth formed EPCAM⁺ and EPCAM^– organoids that were uniformly NKX2-1^– (Figure 2A). The sorted NKX2-1^GFP+ progenitors on day 15 comprised the entirety of cells competent at subsequently expressing the lung-specific marker SFTPC by day 36 (Figure 2C and data not shown), suggesting this population contained lung progenitors. The GFP⁺ outgrowth was also highly enriched for cells competent at expressing lung markers SFTPB and MUC1, although these markers are known to have less lung specificity than SFTPC (Figure 2C). Immunostaining confirmed discrete populations of cells expressing SFTPB and MUC1 proteins in the GFP⁺ outgrowth (Figure 2D). In contrast, TP63⁺ cells were present in both GFP⁺ and GFP^– outgrowths, but were more prevalent in GFP^– outgrowths, suggesting TP63 is not a lung-specific marker in this system. Consistent with this interpretation, TP63⁺ cells in the GFP⁺ outgrowth coexpressed NKX2-1, whereas TP63⁺ cells in the GFP^– outgrowth did not express NKX2-1 (Figure 2D and Supplemental Figure 3, A and B). Reverse transcription quantitative PCR (RT-qPCR) confirmed significantly higher levels of TP63 and the esophageal marker, PITX1, in the GFP^– outgrowth, raising the possibility that the NKX2-1^–TP63⁺ cells might represent alternative foregut derivatives, such as developing esophageal epithelium (Supplemental Figure 3A) (22).

Figure 2

Sorted iPSC-derived NKX2-1^GFP+ cells exhibit lung progenitor potential and ability to form epithelial spheroids in 3D culture. (A) Schematic overview of organoid generation with representative phase contrast and GFP fluorescence microscopy images (day 25–28) as well as immunostaining (day 36) for NKX2-1 and EPCAM proteins. Cell nuclei are counterstained with DAPI. Each panel shows outgrowth in 3D culture of structures arising from iPSC-derived cells that were either unsorted or sorted on day 15 as GFP⁺ versus GFP^– populations (C17). Scale bars: 100 μm (left and center columns) and 20 μm (right column). (B) Time-lapse microscopy (merged GFP fluorescence and phase contrast) of unsorted organoids over 25 hours. Arrows indicate epithelial organoids undergoing induction of the GFP reporter in real time (C17). See also Supplemental Video file. (C) FACS quantification on day 36 of the percentage of cells expressing GFP in the outgrowth wells shown in A (C17). Data indicate individual biological replicates (squares and triangles) with mean ± SD, n = 6 biological replicates. ****P ≤ 0.0001 by Student’s t test. Fold change [RT-qPCR; 2^(–ΔΔCT)] in mRNA expression on day 36 compared with day 0 for each GFP⁺ versus GFP^– outgrowth compared with fetal lung control tissue (n = 3 biological replicates) (C17). (D) Immunostaining of GFP⁺ outgrowth organoids on day 36 for MUC1, SFTPB, and TP63 (C17). Nuclei are counterstained with DAPI. Scale bars: 10 μm.

Fetal LgM augments distal lung differentiation in iPSC-derived lung organoids. Having established that NKX2-1⁺ primordial progenitors could be induced to upregulate markers of lung epithelial lineages without mesenchymal coculture support, we next sought to determine whether these progenitors might also respond to developmental cues provided by primary embryonic LgM. Lung epithelial-mesenchymal interactions are essential for lung epithelial growth, branching, and differentiation (23–25). For example, separating and recombining rat embryonic lung epithelium and mesenchyme (recombinants) has previously revealed both the importance of LgM and the stage-specific plasticity of the developing lung epithelium in response to mesenchymal signals (23). Hence, we asked whether iPSC-derived NKX2-1⁺ cells are competent at responding to developing mouse LgM and if distal lung epithelial gene expression might be induced by distal LgM compared with bronchial mesenchyme or standard directed-differentiation conditions without mesenchyme (Figure 3A). Lung organoids generated from unsorted day 15 human iPSCs were cultured in 3D conditions until day 21. Either NKX2-1^GFP+ or NKX2-1^GFP– areas were microdissected and recombined with E12 mouse lung distal or bronchial mesenchyme (Figure 3A). After 5 to 7 further days in culture as recombinants, we observed continued growth of the human GFP⁺ aggregates in response to mouse distal LgM, continued robust expression of nuclear human NKX2-1 protein, expression of proliferation marker Ki67, and no detectable expression of thyroid markers (Figure 3, B–D and Supplemental Figure 3C). Importantly, we noted induction of cytoplasmic human pro-SFTPC and LPCAT1 protein expression in the majority of the human iPSC–derived cells by immunostaining, with validation of human SFTPC mRNA expression by both in situ hybridization as well as RT-qPCR (Figure 3, B–D and Supplemental Figure 3, D and E). Levels of SFTPC induction were higher in organoids recombined with distal LgM than in those continued solely through directed differentiation without recombination (Supplemental Figure 3D). Since E12 mouse LgM presumably lacks the signals needed for full maturation of the alveolar epithelium, which normally begins in the mouse at E18.5, human cells in these recombinants did not exhibit robust induction of transcripts associated with mature lamellar body biogenesis, such as LAMP3 (data not shown), and did not appreciably display lamellar body–shaped inclusions by microscopy, as expected. Bronchial mesenchyme did not induce SFTPC expression or proximal lung epithelial markers (SOX2, TP63, or SCGB3A1) in NKX2-1^GFP+ cells (Supplemental Figure 3D and data not shown). In addition, control recombinants generated using GFP^– organoids were not competent at inducing either human NKX2-1 or SFTPC expression (Figure 3C and Supplemental Figure 3D). These results suggest that human iPSC–derived NKX2-1⁺ lung progenitors respond to developing distal lung mesenchymal cues, findings in keeping with our prior observations that recombining rat distal LgM with isolated early distal primary lung epithelium induces distal alveolar marker gene expression (25).

Figure 3

Mouse-human recombinant cultures demonstrate that fetal mouse lung mesenchyme augments distal lung differentiation in iPSC-derived human lung organoids. (A) Schematic of microdissecting and combining iPSC-derived lung organoids with E12 mouse lung mesenchyme (LgM). (B) Light microscopy of the same recombinant on days 1, 3, and 5 of in vitro culture. White dashed lines indicate boundaries of mouse LgM and human iPSC-derived epithelium. (C) SFTPC mRNA expression (purple) assessed by in situ hybridization using an anti–human SFTPC probe to stain recombinants generated with GFP⁺ versus GFP^– human iPSC–derived organoids recombined with distal mouse LgM. Scale bar: 500 μm. (D) Immunostaining of NKX2-1^GFP+/distal LgM recombinants for NKX2-1, Ki67, LPCAT1, and pro-SFTPC proteins (with zoom). Brown = immunoperoxidase product after DAB exposure. Scale bars: 200 μm. (E) Representative phase microscopy of day 32 organoids grown from day 15 GFP⁺–sorted progenitors plated at limiting dilution. Shown beneath each image are the cell numbers plated per well of a 96-well plate on day 15. Scale bar: 500 μm. (F) Fold change of human SFTPC mRNA expression in day 32 organoids, generated from sorted day 15 NKX2-1^GFP+ cells at concentrations ranging from 15,000 to 240 cells per well, compared with day 0 by RT-qPCR using the 2^(–ΔΔCT) method. Lines with error bars indicate mean ± SD, n = 8 biological replicates except for sample 240, where only 4 samples had appreciable RNA. B–D were performed with C17, E and F with BU3.

The expression of SFTPC in response to LgM in the majority of epithelial cells analyzed raised the question of whether SFTPC⁺ cells were being derived by the selective outgrowth of rare distal lung–competent day 15 precursors versus the possibility that distal lung–competent progenitors might be common within the NKX2-1⁺ day 15 population. To distinguish these possibilities, we purified day 15 NKX2-1^GFP+ cells and replated increasingly dilute numbers of cells (15,000 cells down to 240 cells per well of a 96-well plate) for further directed differentiation from day 15 to day 32 in 3D cultures without mesenchymal coculture support. This limiting-dilution assay should result in declining SFTPC competence with dilution if only rare distal progenitors are present within the day 15 NKX2-1⁺ population (Figure 3, E and F). On day 32 we observed that lower cell numbers, plated at limiting dilution, resulted in stable to increased SFTPC mRNA expression, consistent with the existence of common rather than rare distal lung–competent progenitors within the NKX2-1⁺ day 15 population and suggesting inhibition of distal differentiation in increasingly dense replating conditions in epithelial-only sphere outgrowths.

Global gene expression kinetics of early human lung development modeled by directed differentiation of human PSCs. We next sought to define the fundamental programs of early human NKX2-1⁺ lung progenitors and their global gene expression kinetics during the course of in vitro directed differentiation. We prepared time series microarray expression profiles representing the following 5 key stages of iPSC lung-directed differentiation (Figure 4A): undifferentiated iPSCs (day 0), definitive endoderm (day 3), anterior foregut–like endoderm (day 6), sorted NKX2-1^GFP+ and NKX2-1^GFP– primordial progenitors (day 15), and sorted NKX2-1^GFP+ and NKX2-1^GFP– differentiated cells (day 28). For positive and negative controls we included primary distal fetal lung epithelial cells (21 weeks of human gestation) and forebrain-like iPSC-derived neural NKX2-1^GFP+ cells (as shown in Figure 1B), respectively (Supplemental Table 1).

Figure 4

Time series global transcriptomic profiling of human iPSCs undergoing lung-directed differentiation. (A) Schematic overview of the lung-directed differentiation with time points and populations analyzed by microarrays. (B) Principal component analysis (PCA) of global transcriptomes of the biological triplicates from time points shown in A. (C) Unsupervised hierarchical clustering by dendrogram of the samples shown in A, based on the top ~1,000 transcripts differentially expressed as determined by ANOVA across all 27 samples. (D) Heatmap of the top 20 transcription factors differentially expressed in each direction (10 up and 10 down) between neural NKX2-1^GFP+ and day 15 lung NKX2-1^GFP+ (15+) populations (ranked by fold change; filtered by FDR < 0.01). (E) Heatmap of the top 10 genes differentially expressed between day 15 lung NKX2-1^GFP+ (15+) and day 15 lung NKX2-1^GFP– (15–; ranked by fold change; filtered by FDR < 0.01). (F) Heatmap of the expression of known markers of neuroectoderm, endoderm, and lung epithelium (separated into progenitor stage and distal/alveolar vs. proximal/airway epithelium). Scale = row-normalized log2 expression. All differentiation samples were from C17 iPSCs.

Principal component analysis (PCA) indicated that the global transcriptome of NKX2-1^GFP+ lung cells was easily distinguished from NKX2-1^GFP+ neural cells (Figure 4B). Unsupervised hierarchical clustering of all 27 samples based on the top ~1,000 transcripts differentially expressed as determined by ANOVA (1,032 genes at P < 5 × 10^–13) revealed that day 15 and 28 cells prepared in the lung-directed differentiation protocol clustered closer to distal fetal lung epithelial controls than to endoderm or neural NKX2-1⁺ cells (Figure 4C). The transcriptional profile of neural NKX2-1^GFP+ cells compared with lung NKX2-1^GFP+ cells (day 15) included 4,329 differentially expressed transcripts (FDR-adjusted P < 0.01) and a distinct set of transcription factors including SIX3, DBX1, OTX1, OTX2, FOXD1, PAX6, and LHX5 (ranked by fold change [FC], filtered by FC > 5, FDR < 0.01, and gene ontology [GO] classification GO:0003700: “transcription factor activity, sequence-specific DNA binding”; Figure 4D), further emphasizing the marked differences between these early NKX2-1⁺ forebrain and lung progenitors. Directed-differentiation protocols have previously used FOXA2 expression as a marker to define an NKX2-1⁺ population as endodermal and TUJ1 to indicate neuroectodermal fate, but we found that FOXA2 is expressed in both neuronal and lung populations (Figure 4F), and that TUJ1 is not highly expressed in the neuronal NKX2-1⁺ population, suggesting neither marker is useful in distinguishing NKX2-1⁺ neural from lung lineages. In contrast, the top 10 transcription factors differentially expressed in day 15 lung NKX2-1^GFP+ versus neuronal NKX2-1^GFP+ cells include 6 genes expressed in the developing lung in vivo (GRHL2, ELF3, GATA6, HNF1B, HOXA1, and FOXA1 (Figure 4D) (26–28), suggesting a constellation of transcripts better able to distinguish these 2 iPSC-derived populations.

To interrogate the differences between day 15 NKX2-1^GFP+ and NKX2-1^GFP– cells at the primordial progenitor stage, we ranked the top 10 differentially expressed genes by FC (Figure 4E). We found that NKX2-1 as well as neighboring long noncoding RNAs (lncRNAs) SFTA3 and NKX2-1AS were highly upregulated in the GFP⁺ population. Also in this top 10 list were genes (BMP3, CRH, and SPOCK3) previously described in lung development (Figure 4E, with validation by RT-qPCR in S4B) (28–34). The finding that SFTA3 (aka NKX2-1–associated noncoding intergenic RNA [NANCI]) is the top differentially expressed transcript in the genome distinguishing day 15 NKX2-1^GFP+ cells is in keeping with recent reports that in developing mouse lungs this transcript is coexpressed with Nkx2-1 and shares the same regional and temporal expression pattern (35). To the best of our knowledge, 4 genes in this list (PALMD, FAM189A2, GRM8, and WDR49) have not been previously identified in the lung epithelium.

From our microarray datasets we selected known markers of definitive endoderm, forebrain, and 24 genes of known importance in the developing lung epithelium and profiled their expression patterns over the course of directed differentiation in comparison with human primary fetal epithelial lung control cells (HFL) (Figure 4F). Endodermal markers, such as GATA4, GATA6, SOX17, NODAL, and FOXA2 were upregulated early during endodermal differentiation, with retained expression of GATA6 and FOXA2 in day 15 and day 28 NKX2-1^GFP+ cells. In contrast, the transcripts NKX2-1, SFTA3, SOX9, and FOXP family members were low or absent prior to day 6, and their clear emergence in the day 15 GFP⁺ population is consistent with their published expression during the early lung progenitor period in mouse lung development (35, 36). These findings, together with the lack of mature lung marker gene expression in day 15 GFP⁺ cells (low SCGB1A1, SCGB3A2, TP63, SFTPB, and SFTPC), further suggest the day 15 GFP⁺ population represents a relatively undifferentiated or primordial lung progenitor population, as has been observed in early NKX2-1⁺ progenitors in developing mouse cells in vivo (28). In contrast to day 15, by day 28 the GFP⁺ population had begun to express markers known to be enriched in maturing alveolar epithelial cells (ETV5, CLDN18, LPCAT1, MUC1, SFTPB, and low SFTPC) or in airway epithelia, such as basal (TP63), secretory (SCGB3A2, MUC5B, MUC5AC, and AGR2), and neuroendocrine (ASCL1) cells. PDPN, which has occasionally been referred to as a PSC-derived type 1 pneumocyte marker in prior publications (8, 9), was actually expressed in day 0 as well as day 15 GFP⁺ cells, consistent with its expression patterns in developing mice where it is robustly expressed in both the foregut endoderm and the developing pseudoglandular lung epithelium prior to the emergence of type 1 cells (37).

Next we sought to identify unbiased gene signatures of primordial (day 15) and maturing (day 28) NKX2-1^GFP+ cells. We generated lists of the top 100 differentially expressed genes (ranked by FC, filtered by FDR < 0.01) of each sample across multiple comparisons and identified a common gene set for each sample (Supplemental Figure 4 and Supplemental Table 2). The day 28 GFP^– population was enriched for diverse but predominantly liver (APOA2, FGB, AFP, CDH17, and TF) and intestinal (CDX2, CDH17, and GIF) markers (Supplemental Table 2). In addition to NKX2-1, SFTA3, CPM, NFIB, and CRH, which are all expressed in primordial lung progenitors, the maturing lung cells (day 28 GFP⁺) expressed higher levels of SCGB3A2, SFTPB, TP63, ICAM1, IL8, and ITGB6 (Supplemental Figure 4 and Supplemental Table 2). SCGB3A2 was the most differentially expressed transcript of 23,786 probesets ranked by FC (day 28 GFP⁺ vs. GFP^– groups; FC = 76.6; FDR-adjusted P = 1.6 × 10^–9). SFTPC was upregulated by day 28 (GFP⁺), but not yet at levels equivalent to HFL, findings in keeping with our recombinant experiments, which suggested that current differentiation protocols without the use of primary mesenchyme have not yet been optimized for efficient and full distal alveolar maturation (Figure 3, C and D, Supplemental Figure 3D, Figure 4F, and Supplemental Table 2).

Transcription factors play critical roles in organogenesis including in lung development (27). To identify candidate genes that control human lung specification and development, we screened for enrichment of transcription factors or regulators of transcription prior to lung specification (day 6) and at different stages of lung maturity (day 15 NKX2-1^GFP+, day 28 NKX2-1^GFP+, and HFL; Figure 5). To identify genes of interest, we ranked the significantly differentially expressed genes by FC (FC > 4; FDR < 0.01) and filtered genes based on GO classification for transcription factor activity (GO:0003700, “transcription factor activity”). The majority of the most highly differentially expressed genes in day 6 anterior foregut–like endoderm were known transcription factors of the foregut endoderm previously described in Xenopus and mouse model systems: HHEX, GATA3, GATA4, FOXC1, EOMES, OTX1, OTX2, ISL1, and PITX2 (28, 38, 39) (Figure 5, Supplemental Figure 4C, and Supplemental Table 3). In comparison with day 6, the day 15 NKX2-1^GFP+ population was enriched for many transcription factors known to be present in the developing mouse lung (JUN, MECOM, SOX2, HES1, HOXA1, NKX2-1, FOXA1, ELF3, ELF5, NFIB, and FOXP2) (26–28, 38–44) (Figure 5, Supplemental Figure 4D, and Supplemental Table 3). In addition to HOXA1, a number of other HOX genes were upregulated in day 15 NKX2-1^GFP+ samples (HOXA4 and HOXC4). SHH, essential for normal lung development in mice (27), was the most highly expressed gene in the day 15 samples in our analysis (FC = 60, FDR = 4.8 × 10^–11). Day 28 NKX2-1^GFP+ cells expressed higher levels of transcription factors associated with basal cells (TP63) and neuroendocrine cells (ASCL1). Taken together, these time series data provide unbiased stage-dependent signatures of the putative transcriptomic programs of human lung progenitors and their differentiated progeny as they emerge during developmental directed differentiation. Moreover, these signatures reveal that many evolutionarily conserved transcription factors, previously observed in developing Xenopus and mouse lung endoderm in vivo, also are differentially expressed in the iPSC human lung development model system.

Figure 5

Transcriptomic signatures of iPSC-derived anterior foregut endoderm and lung progenitors, focused on genes associated with transcription factor activity. Heatmap of the top transcription factors or genes with transcription factor activity GO terms that are differentially expressed (filtered by fold change > 4) between successive stages of lung differentiation starting from day 6 and compared with human fetal lung. Vertical text box and heatmap with black outline identify the signature of each stage being analyzed. *Indicates the gene was previously described in foregut and/or lung development (see Supplemental Table 3 for further publication details). See also Supplemental Figure 4. Day 6, day 15, and day 28 time points were derived from C17 iPSCs.

Single-cell RNA sequencing and surface marker profiling of day 15 iPSC-derived lung progenitors. While transcriptomic profiles of purified groups of cells allows a deep understanding of the genetic program of the NKX2-1⁺ progenitor population, it does not allow interrogation of the heterogeneity of these programs at the individual cell level. Hence, we next sought to profile the transcriptomic programs of individual iPSC-derived cells at the day 15 stage of lung differentiation employing the C17 NKX2-1^GFP targeted line as well as the untargeted iPSC line, BU3 (19). We performed RNA sequencing (RNA-Seq) of 84 BU3 iPSC-derived cells without any cell sorting and 69 C17 iPSC-derived cells, sorted based on NKX2-1^GFP+ expression (Figure 6A).

Figure 6

Single-cell RNA sequencing of sorted and unsorted iPSC-derived cells reveals NKX2-1⁺ lung and NKX2-1^– non-lung lineages and suggests markers for their identification. (A) Schematic of the single-cell capture and global RNA sequencing of sorted C17 NKX2-1^GFP+ and unsorted BU3 iPSCs on day 15 of lung-directed differentiation. (B) PCA of the top 150 most variant genes of all sequenced cells reveals 4 cell clusters, color coded to match the clusters shown in panel C. (C) Heatmap of gene expression (global Z score) with unsupervised hierarchical clustering of all 153 cells (x axis) and the top 150 most variant genes (y axis). Dendrograms as well as colored boxes indicate 4 cell clusters (CC1–4; matching colors shown in panel B). y-axis dendrograms and thick black lines indicate 3 gene clusters (GC1–3). Key genes indicated on right. CC1 = orange, CC2 = blue, CC3 = red, CC4 = green. (D) Top 10 genes correlated with NKX2-1 expression. (E) Unsupervised cell clustering using Monocle’s pseudotime spanning tree analysis reveals 7 cell states. Individual cells are labeled in subsequent panels by NKX2-1 level, genetic background (iPSC clone), or cell cycle (mitosis), respectively. (F) Pseudotime plots of expression levels of NKX2-1, CD47, SOX9, NFIB, FGB, and APOA2 with cells colored based on the 7 states determined in E. (G) Immunostaining of RUES2-derived day 15 cells for NKX2-1 (red) and SOX9 (green) nuclear proteins. Nuclei counterstained with DAPI. Scale bars: 25 μm.

PCA (Figure 6B) as well as unsupervised hierarchical clustering analysis (Figure 6C) both suggested that 4 broad clusters of cells (hereafter CC1-4) were present on day 15, and these cell clusters could be distinguished based on 3 broad gene clusters (GC1-3). Notably, 3 of the 4 clusters (CC1, CC2, and CC4) expressed high levels of NKX2-1 (Figure 6C), whereas cluster CC3 exhibited an absence of transcripts encoded by the NKX2-1 locus or its neighboring locus, SFTA3. NKX2-1–expressing clusters were most robustly distinguished by whether they exhibited mitotic (CC1 and CC2) or noncycling (CC4) gene signatures. For example, CC1 and CC2 were highly enriched for the expression of genes associated with mitosis or cytokinesis (e.g., KIF11, KIF14, KIF22, KIF23, CDC20, and AURKB). Thus, CC1 was labeled mitotic and CC2 was labeled premitotic based on slightly lower expression levels of these markers in the latter cluster, whereas CC3 and CC4 did not appear to be in active cycle at the time of cell capture for analysis. Importantly, NKX2-1⁺ cells clustered together (CC4) regardless of whether they were sorted GFP⁺ C17 iPSCs or unsorted BU3 iPSCs. Furthermore, only 1 GFP⁺ sorted cell could be found misclustering among the 26 cells that comprised the NKX2-1–negative cluster (CC3), and as expected 25 out of 26 cells found in this NKX2-1–negative cluster derived from the unsorted BU3 iPSCs (Figure 6C).

We undertook 3 approaches to interrogate the gene expression differences that distinguished each cell cluster. First, we used unsupervised hierarchical clustering of the top 150 most variant genes (Figure 6C; y-axis dendrograms; hereafter gene clusters GC1–3). We found that GC1 was highly enriched for cell-cycle regulation genes (including AURKB, BIRC5, BUB1, CCNB1, CCNB2, CENPE, CENPF, KIF11, KIF14, KIF22, KIF23, KIF22c, MELK, and TOP2A; Figure 6C), further supporting the interpretation that changes in genes of cytokinesis and cell cycle dominate the first level of clustering of day 15 cells. However, 2 additional distinct gene clusters were also apparent, most notably GC2, including NKX2-1, SFTA3, NFIB, CD47, WNT5A, CPM, and LAMA3; and GC3 which was associated with the NKX2-1–negative/SFTA3-negative cells of CC3 (Figure 6C). This analysis suggested GC2 genes as potential markers associated with NKX2-1⁺ lung cells, and GC3 genes as potential markers associated with non-lung (NKX2-1–negative) cells. Consistent with this observation, the top 10 genes most highly correlated with NKX2-1 expression across individual cells were GC2 genes, including CD47, SFTA3, CPM, and LAMA3 (Figure 6D). In contrast, GC3 was enriched in liver-lineage genes (APOA2 and FGB) as well as nonspecific mesenchymal genes (COL19A1 and S100A10). We have previously reported that in iPSC-derived hepatic cells, FGB represents the most upregulated transcript in the genome during hepatic-directed differentiation (45, 46). Furthermore, in postnatal human tissues both APOA2 and FGB are transcripts specifically enriched in liver cells (gtexportal.org). Significantly, all NKX2-1–negative cells of putative hepatic lineage (20 out of 20 APOA2⁺ cells; Figure 6C) were solely composed of unsorted BU3 iPSCs, indicating that sorting on the NKX2-1^GFP+ marker successfully depleted any contaminating hepatic cells in this protocol. The majority of cells (23 of 37) in the mitotic and premitotic groups clustered with NKX2-1⁺CD47⁺ cells when hierarchical clustering was run again after cell cycle genes were removed (Supplemental Figure 5A).

Second, we applied the Monocle computational algorithm (47) to our single-cell dataset in an effort to identify, in an unbiased manner and irrespective of cell cycle, cell subtypes, or intermediate states. Because a particular challenge in single-cell RNA-Seq experiments is the high cell-to-cell variation observed in most genes, including key developmental regulators, during differentiation (48–51), Monocle was developed to improve the resolution of individual transcriptomes and allow the ordering of cells by potential progress through a biological process without relying on known lineage markers (47). Hence, we used Monocle to order all day 15 cells in pseudotime: an abstract, semiquantitative measure of progress through a biological process (Supplemental Figure 5B). In this analysis, pseudotime represents a computational, high-dimensional ordering of the transcriptional spectrum of differentiating cells, accounting for the likelihood that day 15 cultures contain diverse cell types at various stages of differentiation. In light of the dominant effect of proliferation on the first-level analysis (Figure 6C), we excluded cell-cycle genes from the Monocle analysis. Unsupervised cell clustering revealed 7 Monocle states (labeled states 1–7; Figure 6E and Supplemental Figure 5B), composed of 24, 8, 10, 57, 21, 13, and 12 cells, respectively. Next, we labeled cells in each state based on expression levels of NKX2-1 (high vs. no/low expression; Figure 6E) and determined that cells in states 5 and 6 clustered separately because they expressed no or low levels of NKX2-1, whereas cells in states 1, 2, 4, and 7 were almost entirely NKX2-1 high. In contrast, most states were neither defined by cell cycle effects (mitotic state) nor genetic background (cell origin) of each iPSC line, with the exception of states 5 and 6, which were entirely composed of unsorted BU3 cells, as expected (Figure 6E). Consistent with our clustered heatmap results (Figure 6C), we found that states 1, 2, 4, and 7 cells expressed higher levels of NKX2-1, CD47, and SFTA3, consistent with a lung signature, whereas state 5 cells were enriched for APOA2 and FGB, consistent with a fetal liver signature (Figure 6, Supplemental Figure 5D, and data not shown). The lineage identity of state 6 cells, enriched for transcription factors including MSX1, EGLN3, and OTX2, was uncertain. The presence of discrete NKX2-1⁺ states suggested some degree of either temporal or lineage heterogeneity within the overall NKX2-1 population. For example, SOX9 expression varied across cells in CD47⁺ NKX2-1⁺ states. State 1 and 2 cells were significantly enriched for SOX9 expression, whereas state 4 cells were not. State 3 cells were highly enriched for SOX9 but expressed lower levels of NKX2-1. Using Monocle to examine genes that follow similar expression trends when cells are ordered in pseudotime, we identified genes with increasing expression towards a lung phenotype (including SOX9, NKX2-1, CD47, NFIB, LAMA3, and SFTA3) and conversely genes with increasing expression towards a liver phenotype (DLK1, AFP, MSX1, FTL1, FN1, FGB, and APOA2) (Figure 6F and Supplemental Figure 5, C and D). The expression of SOX9 at variable but easily detected levels in the majority of NKX2-1⁺ putative lung cells was confirmed at the protein level by immunostaining (Figure 6G), and only a minor subset of NKX2-1⁺ cells expressed high levels of the proximal airway patterning marker SOX2 without SOX9 (either by immunostaining or by supervised hierarchical clustering of single-cell transcriptomes; Supplemental Figure 6A). This predominance of the distal progenitor marker, SOX9, in day 15 NKX2-1⁺ cells is consistent with the efficient distal alveolar differentiation competence of the NKX2-1⁺ progenitor population observed in our recombinant cultures in Figure 3. Importantly, mature distal or proximal markers (e.g., SFTPC and SCGB1A1) were not detected in any cell on day 15 (Supplemental Figure 6A).

Finally, we performed repeat ANOVA with hierarchical cell and gene clustering of the 97 cells that were not in active mitosis (focusing solely on CC3 and CC4; Figure 6C and Supplemental Figure 6B). This analysis identified 4 cell subgroups (SG1–4; Supplemental Figure 6B). The largest subgroup of cells, SG2, expressed key genes of the early developing lung (SFTA3, NFIB, and WNT5A) and, in keeping with our population-based transcriptomic profiles, they lacked detectable expression of markers of lung maturation (SFTPC, SFTPB, SCGB3A2, ASCL1, FOXJ1, or SCGB1A1) (Supplemental Figure 6B). In the 3 remaining minor subgroups (SG1, 3, and 4; comprising predominantly NKX2-1–negative BU3 cells), we found differentially expressed genes suggestive of non-lung endoderm or undetermined identity. Gene sets most significantly correlated with SG1 and SG3 were consistent with hepatic lineages (Supplemental Figure 6B). Significance testing of SG2 versus SG3 demonstrated that CD47 was the most highly differentially expressed gene in SG2 (ranked by either P value or correlation with NKX2-1 expression; Supplemental Figure 6C), followed by IGFBP5, SFTA3, EIF1AY, LAMA3, CPM, SOX9, and LMO7 (Supplemental Figure 6C). LMO7 is a known target of FGF10 that is upregulated in early developing mouse lung epithelium (52), and NKX2-1, SFTA3, CPM, and SOX9 are all known to be enriched in developing mouse and human lung epithelia.

To determine whether key markers identified in our single-cell RNA-Seq or microarray analyses are expressed in the developing human lung epithelium in vivo, we analyzed available microarray data of human fetal lungs ranging from 53 to 154 days of gestation (53). Consistent with our PSC in vitro model system, we observed increasing in vivo expression with time of known lung differentiation markers (SFTPC, SFTPB, and LAMP3), absence at any time point of non-lung markers (APOA2 and CDX2), and early, unchanging expression of NKX2-1, CD47, NFIB, HOXA3, and JUN (Supplemental Figure 6D). In addition, we confirmed early developmental CD47 protein expression in NKX2-1⁺ epithelial cells in vivo by immunostaining week 10 human fetal lung (Supplemental Figure 6E).

Taken together, our results provided an improved understanding of the heterogeneity of iPSC-derived cells emerging with lung-directed differentiation, supported the utility of NKX2-1^GFP+ sorting to deplete non-lung endodermal lineages that contribute to this heterogeneity, and suggested transcripts associated with NKX2-1⁺ cells in this in vitro model system.

Prospective isolation of iPSC-derived NKX2-1⁺ primordial lung progenitors by CD47^hi cell sorting. Because our single-cell RNA-Seq profiles revealed that CD47 was the transcript in the genome most highly correlated with NKX2-1 (Figure 6D), we sought to determine whether NKX2-1⁺ primordial progenitor cells might be prospectively isolated based on cell surface protein expression of CD47 without the need for a GFP knockin reporter. Using both immunofluorescence microscopy as well as FACS of day 13–15 unsorted PSCs (C17, BU3, and RUES2 lines), we observed that the brightest CD47⁺ cells selectively coexpressed NKX2-1 nuclear protein as well as the NKX2-1^GFP reporter (Figure 7, A and B and Supplemental Figure 7).

Figure 7

Cell surface profiling and prospective isolation of iPSC-derived NKX2-1⁺ primordial lung progenitors by CD47^hiC26^lo cell sorting. (A) Day 15 iPSCs after lung-directed differentiation, immunostained for CD47 and NKX2-1 proteins (C17). Nuclei are counterstained with DAPI. Scale bars: 50 μm. (B) Flow cytometry dot plots of 4 cell surface markers identified in a screen of 243 surface markers on day 15 of lung-directed differentiation (C17); CD26 is depleted while CD47, ALCAM, MUC1, and CPM exhibit higher expression in the GFP⁺ population. (C) Schematic of experimental data for C, D, and E. Flow cytometry dot plots of live day 13 cells (BU3) indicates staining with isotype control antibodies (left panel) or antibodies against CD47 and CD26. Sort gates identify presorted cells (gray box) versus a CD47^hiCD26^lo population (red box) or a CD47^lo population (black box) profiled in D and E. (D) Fold change of NKX2-1 mRNA (left graph, n = 4) in each indicated day 13 population, and SFTPC mRNA (right graph, n = 3) expression in the outgrowth of each indicated population on day 36 compared with day 0 iPSCs by RT-qPCR; 2^(–ΔΔCT). Data indicate the mean ± SD. *P ≤ 0.05, **P ≤ 0.01, ****P ≤ 0.0001 by Student’s t test. (E) GFP expression quantified in each of the gates shown in C: day 13 presort = 62% NKX2-1^GFP+; CD47^hiCD26^lo = 97% GFP⁺; CD47^lo = 8% GFP⁺. Lower panel is FACS of day 36 outgrowth of each indicated day 13 sorted population: GFP⁺, CD47^hiCD26^lo versus CD47^lo. (F) Phase contrast and fluorescence microscopy (GFP) of day 42 organoids derived from day 13 sorted GFP⁺, CD47^hiCD26^lo, and CD47^lo populations from E. Scale bars: 100 μm. (G) Confocal microscopy of outgrowth organoids (C17), sorted on day 13 based on CD47^hiCD26^lo and analyzed on day 44 by coimmunostaining for NKX2-1 (green) and pro-SFTPC (purple). The individual panels for this merge image are contained in Supplemental Figure 7E. Scale bars: 25 μm.

In independent experiments we screened day 15 iPSC-derived NKX2-1^GFP+ progenitors by FACS using a panel of 243 antibodies, validating that CD47 was in the top 4 cell surface markers most associated with GFP⁺ expression (Figure 7B). Notably, this screen also confirmed that: (a) NKX2-1^GFP+ cells are EPCAM positive (Supplemental Figure 7A); (b) CPM, recently published as a marker of iPSC-derived NKX2-1⁺ cells (16) and highly associated with NKX2-1 in our single-cell RNA-Seq (Figure 6D), indeed costains most NKX2-1^GFP+ cells (Figure 7B), although it is also associated with NKX2-1^– hepatic cells (54) that emerge at low levels in this protocol (Figure 6F and Supplemental Figures 5B and 6D); and (c) ALCAM (CD166) and MUC1 (CD227) are 2 additional candidate markers that selectively identify NKX2-1^GFP+ cells at this stage (Figure 7B). Importantly, our antibody screen also suggested CD26 as a negative-selection marker since the brightest CD26⁺ cells were lower in expression of the NKX2-1^GFP reporter (Figure 7B). When sorting day 13–15 cells based solely on CD47^hiCD26^lo gating, we found significant enrichment for NKX2-1⁺ cells: 89% ± 4.1% of cells (mean ± SD; n = 11 runs) expressed NKX2-1 nuclear protein as well as the NKX2-1^GFP reporter compared with CD47^lo cells, which were depleted of NKX2-1 expression (Figure 7, C–F). These findings were validated for both human ESCs (RUES2) as well as multiple iPSC lines (C17, BU3, 100-3, RC202, and RC204; refs. 18, 19, 55 and data not shown), despite varying efficiencies of NKX2-1⁺ induction in any given differentiation run (Supplemental Figure 7, C and D). For example, regardless of whether a line differentiated to lung with low or high efficiency (e.g., 13% vs. 56% NKX2-1⁺ in 3 separate runs for C17; or 50% for RUES2), in each case CD47^hiCD26^lo gating provided significant enrichment in NKX2-1⁺ cells (7-fold vs. 2-fold enrichment for C17 and 2-fold for RUES2), resulting in populations approximately 90% pure for NKX2-1 expression in each run. Furthermore, human ESCs or human iPSCs sorted solely on CD47^hiCD26^lo gating produced predominantly NKX2-1⁺ spheroids in 3D culture that expressed lung differentiation markers including pro-SFTPC protein and SFTPC mRNA at levels similar to those of sorted NKX2-1^GFP+ cells (Figure 7, D and G).

Our results indicate that iPSC-derived lung epithelial cells originate from identifiable NKX2-1⁺ progenitors. Through the use of an NKX2-1–targeted GFP reporter these progenitors can be sorted and then further differentiated without mesenchymal coculture support in 3D Matrigel culture. Importantly, human NKX2-1⁺ progenitors derived with our methods undergo efficient distal SFTPC⁺ differentiation and proliferation after recombinant culture with primary distal embryonic mouse LgM. We have labeled these NKX2-1⁺ cells, which emerge between days 8 and 15 of iPSC differentiation, primordial progenitors because they express a transcriptome that includes the earliest transcripts known to emerge during the endodermal and primary lung bud stages of mammalian development (NKX2-1, SFTA3, SOX9, and SOX2), but they are otherwise lacking in transcripts associated with differentiated/maturing lung epithelia, most of which emerge during the later pseudoglandular stage of lung development.

Our findings support a paradigm in which the human lung epithelium derives directly from NKX2-1⁺ endodermal progenitors rather than from alternate cells, because sorting human NKX2-1⁺ cells at the primordial stage highly enriches for cells competent at further differentiating into lung epithelia, while depleting this population significantly depletes cells competent at forming lung. Given the difficulty in accessing and tracking live human fetal cells in vivo during the earliest stages of lung development (approximately 1 month of human gestation), our in vitro model enables the purification, tracking, and visualization of cells undergoing the earliest moments of human lung cell fate decisions, a time period in human lung development that remains elusive to scientific study.

We provide evidence that the genetic control of early human lung development is similar to that in mouse. Indeed, our finding that human iPSC-derived lung epithelial progenitors respond to inductive differentiation cues provided by developing mouse LgM suggests that an evolutionarily conserved biology is common to early mouse and human lung development. Furthermore, our stage-dependent transcription factor signatures for developing lung, thyroid, and forebrain revealed by our in vitro iPSC model provides important validation in a human system of many of the gene ensembles previously identified in mice. For example, we found that many genes and transcription factors of murine lung development are expressed in iPSC-derived human lung progenitors, including NKX2-1, SHH, FOXA2, FOXA1, GATA6, SOX2, SOX9, GRHL2, IRX1, IRX2, NFIA, NFIB, FOXP2, HNFIB, ELF3, and ELF5. Our developmental stage–dependent signatures also suggest novel genes requiring further study, and the ability to employ the human iPSC model system should now provide a tractable developmental human system to examine the roles of these and other genes in human lung lineage specification. Partnered future work focused on this early period of lung lineage specification in vivo in mice is likely to provide a further understanding of the phenotype and biology of primordial lung progenitors across mammalian species, and should enable isogenic head-to-head comparisons of iPSC-derived progenitors with their in vivo counterparts. It is important to point out that the differentiation kinetics of our human iPSC in vitro model is generally faster than that observed in vivo in developing human endoderm and lungs. Thus, we cannot exclude the possibility that cells in the in vitro developmental model take a slightly different or alternate path towards lung cell fates, and we propose that further comparisons of our cells to in vivo developing cells, both mouse and human, will be important to address this possibility.

Given the primordial nature of early NKX2-1⁺ lung progenitors, whether derived from iPSCs/ESCs or emerging in vivo in embryos (2, 56), specific lung progenitor markers have not been previously identified with certainty to enable their prospective isolation. While FOXP2 in mice has been proposed as a lung primordial marker (56) and CPM has been proposed as a cell surface marker for sorting NKX2-1⁺ cells derived from human iPSCs (16), the previous lack of any available tool for specifically tracking or purifying live NKX2-1⁺ cells has left uncertainty regarding the specificity of those markers for prospective lung progenitor isolation. Our profiling of NKX2-1^GFP+ primordial progenitors by microarrays, single-cell RNA-Seq, and FACS-based screens reveals a cell surface phenotype, CD47^hiCD26^lo, which can be used to prospectively isolate NKX2-1⁺ progenitors when derived from iPSCs in culture. Our findings validate the utility of CPM, recently published by Gotoh et al., to also serve as a cell sorting marker (16), although CPM is also expressed in NKX2-1^– hepatic cells which emerge in this lung-directed-differentiation protocol. While FOXP2 is enriched in the NKX2-1⁺ progenitor population, our results indicate that it lacks lung specificity and is also expressed elsewhere in the protocol, for example, in day 15 NKX2-1^– cells (Supplemental Table 1).

Like other published markers for the purification of iPSC-derived endoderm (e.g., CKIT/CXCR4), CD47 is broadly expressed in many tissues. However, in the iPSC model system it has particular utility as a marker that allows sorting of NKX2-1⁺ lung progenitors with ~90% purity based on its unexpectedly high levels of cell surface expression compared with other cells. CD47 is a broadly expressed cell surface glycoprotein with diverse roles in cellular processes including apoptosis, proliferation, and migration. The extracellular domain of CD47 acts as a thrombospondin-1 (TSP-1) receptor but also interacts with integrins and SIRPα (57). CD47 is expressed in lung epithelial cells in vivo and in vitro where it has a role ascribed to regulating leukocyte migration into the lung (58). Future work will be needed to determine if CD47 has a unique functional role during early lung development.

Just as our results provide an increased understanding of the level of heterogeneity present in iPSCs undergoing initial lung lineage specification in culture, further work is needed to interrogate the increasing heterogeneity that appears to emerge with each subsequent lung differentiation step. Similar to our NKX2-1^GFP tool, engineering multicolored reporters that become activated in more differentiated lung epithelial lineages should gradually facilitate this understanding and enable purification of subsets of increasingly mature cells in order to understand and overcome this obvious heterogeneity. These approaches should facilitate modulations of later lung differentiation stages in order to efficiently pattern cells into proximal versus distal lineages and their downstream progeny.

In summary, we have purified human lung progenitors derived from iPSCs, and these cells are reminiscent of early stages of lung developmental differentiation. Our profiling of these cells as well as their precursors and progeny during the time course of directed differentiation has resulted in an understanding of their global transcriptomic programs at the single-cell level and provides a validated set of cell surface markers and transcription factors selectively enriched in these cells. It should now be possible to test whether cells of similar primordial lung progenitor phenotype remain in the lung postnatally or can be rederived in patients during responses to injury. Given the broad differentiation repertoire of the primordial progenitors, we anticipate that access to pure populations of these cells should facilitate basic developmental studies as well as clinical applications focused on disease modeling, drug development, and potentially future regenerative therapies.

Human PSC maintenance and gene editing. Previously published PSC lines iPSC17 (18), BU3 (19), and WA09 (H9 ESC) were maintained in feeder-free conditions on Matrigel (Corning) in mTeSR1 (Stem Cell Technologies) and passaged with Gentle Cell Dissociation Reagent (Stem Cell Technologies). In order to generate NKX2-1^GFP reporter PSC lines, TALENs or CRISPR-based technologies were employed as detailed in the supplement to introduce a DNA double-stranded break into the second intron of NKX2-1 (Supplemental Figure 1A). A donor matrix containing a splice acceptor, NKX2-1 exon 3, 2A-eGFP, and loxP-flanked PGK-puroΔTK selection cassette was integrated by homologous recombination and targeted PSC clones resistant to puromycin underwent Cre-mediated excision of the loxP-flanked puroΔTK selection cassette (Supplemental Figure 1B), followed by confirmation of cassette excision, karyotyping, and characterization as detailed in the supplement.

Human iPSC directed differentiation. Neuroectodermal NKX2-1^GFP+ cells were generated using STEMDdiff Neural Induction Medium (STEMCELL Technologies) according to the manufacturer’s protocol. On day 6 of differentiation, 2 μM purmorphamine (Stemgent) was added to the base media. NKX2-1⁺ cells were sorted on days 12–14. Thyroid NKX2-1^GFP+ cells were derived using our recently published protocol (19) and lung NKX2-1^GFP+ cells were generated by adapting published protocols (8, 9). For both lung and thyroid differentiations we first induced definitive endoderm using a STEMdiff Definitive Endoderm Kit (STEMCELL Technologies) according to the accompanying protocol. After approximately 72 to 84 hours, we harvested and analyzed cells by flow cytometry for efficiency of definitive endoderm induction by the coexpression of the surface markers c-KIT and CXCR4 (see supplemental experimental procedures). After definitive endoderm induction cells were plated in small clumps at approximately 50–150,000 cells/cm² on Matrigel-coated plates in complete serum-free differentiation media (CSFDM) supplemented with 2 μM dorsomorphin (Stemgent) and 10 μM SB431542 (Tocris) for 72 hours. Y-27632 (10 μM, Tocris) was added for the first 24 hours. To specify thyroid epithelium, differentiation medium was changed on day 6 to CSFDM supplemented with 250 ng/ml rhFGF2 (R&D Systems), 100 ng/ml rhBMP4 (R&D Systems), and 100 ng/ml heparin sodium salt (Sigma-Aldrich) according to our published methods (19). To specify lung epithelium, differentiation medium was changed on day 6 to CFKBRa: CSFDM supplemented with 3 μM CHIR99021 (Tocris), 10 ng/ml rhFGF10, 10 ng/ml rhKGF, 10 ng/ml rhBMP4 (all from R&D Systems), and 50–100 nM Ra (Sigma-Aldrich) (9). Day 15 cells were dissociated with 0.05% trypsin (Thermo Fisher Scientific) followed by resuspending as small clumps in CSFDM supplemented with 3 μM CHIR99021, 10 ng/ml rhFGF10, and 10 ng/ml rhKGF (CFK medium) and plated on freshly-coated Matrigel (Corning, 354277) plates. Y-27632 (10 μM) was added to CFK medium for the first 24 hours. On day 22 the medium was changed to CFK+DCI: CFK medium plus 50 nM dexamethasone (Sigma-Aldrich), 0.1 mM 8-bromoadenosine 3′,5′-cyclic monophosphate (8-Br-cAMP) sodium salt (Sigma-Aldrich), and 0.1 mM 3-isobutyl-1-methylxanthine (IBMX) (Sigma-Aldrich).

Sorting iPSC-derived lung progenitors and organoid generation. On day 15, live cells identified by propidium iodide (PI) exclusion were sorted by flow cytometry based on GFP expression for downstream applications including RNA analysis and generating organoids. To generate unsorted organoids, day 15 cells were dissociated with trypsin and pelleted cells were resuspended in Matrigel (Corning, 356230). CFK medium was then added to each well, supplemented with 10 μM Y-27632 medium for the first 24 hours. To generate sorted NKX2-1^GFP+ or NKX2-1^GFP– organoids, we resuspended the relevant sorted populations in Matrigel at a density of 50,000 cells per 50 to 100 μl Matrigel and allowed to gel as above. Further methods are detailed in the supplement.

Single-cell RNA-Seq analysis of day 15 iPSC-derived lung progenitors. Day 15 NKX2-1^GFP+ and BU3 unsorted cells were generated using the lung protocol, dissociated, and sorted as described above. Fluidigm C1 integrated fluidics circuits (IFCs) were used to capture individual live cells, lyse, convert polyA⁺ RNA into full-length cDNA, amplify and generate cDNA according to the manufacturer’s protocol (“Using C1 to Generate Single-Cell cDNA Libraries for mRNA Sequencing”, Fluidigm, PN 100-7168). Paired-end 50-bp reads were aligned with an average of 2.8 million reads per cell per end. Aligned, mapped reads were statistically analyzed using methods detailed in the supplement. The clustering, PCA, and significance testing were performed using SCICAST (details and a walkthrough can be found at https://github.com/iandriver/SCICAST) with additional hierarchical clustering linkage, Pearson’s correlation coefficients, ANOVA, and FDR-adjusted P value calculation methods detailed in the supplement. Unbiased cell clustering was performed and cells ordered in pseudotime using Monocle 2 (47) (see supplemental experimental procedures). All raw data gene expression files can be downloaded from the NCBI Gene Expression Omnibus (GEO GSE96106).

Isolation of primary human fetal lung epithelium. Week 10 and 21 human lung tissues or their purified cell derivatives were isolated in the Guttentag lab using methods detailed in the supplement and previously published (59, 60).

Microarray analysis. Biological triplicates of all samples except human fetal lung were prepared. Biological duplicates from 1 embryo (uncultured naive lung epithelium) and a singlicate from a different embryo (differentiated AT2 cells) were prepared for the human fetal lung sample controls. Global gene expression in all 27 samples was analyzed by Affymetrix GeneChip Human Gene 2.0 ST arrays using methods and computational analyses detailed in the supplement. Differential gene expression with respect to experimental group across all samples was assessed by performing a 1-way ANOVA with correction for multiple hypothesis testing using the Benjamini-Hochberg FDR. All raw data gene expression files can be downloaded from the NCBI Gene Expression Omnibus (GEO GSE83310).

Recombination with mouse embryonic LgM. Recombinations were performed essentially as previously described (25). Briefly, small GFP⁺ or GFP^– fragments of day 21 human iPSC organoids were recombined with 10 to 12 pieces of mouse E12 LgM manually dissected free of any epithelia as described (25). The LgM rudiments were teased into close apposition to the human fragments with microsurgery knives (Fine Science Tools, Inc.). After overnight culture to promote tissue adherence, the recombinants were transferred to the surface of an 8-μm pore size Whatman nuclepore filter and cultured for 5 to 7 days in BGJb medium containing 20% FBS, 0.2 mg/ml vitamin C (Sigma-Aldrich), and 5 μg/ml recombinant mouse amino-terminal SHH (R&D Systems) to promote mesenchyme viability (61). The recombinants were maintained for 7 days, with medium changes every other day. Dexamethasone (50 nM) was added to the medium for the final 48 hours to promote lung epithelial differentiation.

RT-qPCR. RNA extracts were converted to cDNA and analyzed during 40 cycles of real-time PCR using Taqman probes (Applied Biosystems) according to methods detailed in the supplement. Relative gene expression, normalized to 18S control, was calculated as FC in 18S-normalized gene expression, over baseline, using the 2^(–ΔΔCT) method. Unless otherwise specified in the text, baseline, defined as FC = 1, was set to undifferentiated stem cell levels, or if undetected, a cycle number of 40 was assigned to allow FC calculations.

Statistics. Statistical methods relevant to each figure are outlined in the accompanying figure legend. Unless otherwise indicated unpaired, 2-tailed Student’s t tests were applied to 2 groups of n = 3 or more samples, where each replicate (n) represents either entirely separate differentiations from the PSC stage or replicates differentiated simultaneously and sorted into separate wells. A P value less than 0.05 was considered to indicate a significant difference between groups.

Study approval. The IRB of Boston University approved the generation and differentiation of human iPSCs with documented informed consent obtained from participants. Human lung tissues were obtained in the Guttentag laboratory under protocols originally reviewed by the IRB at the Children’s Hospital of Philadelphia and subsequently reviewed by Vanderbilt University. Experiments involving mouse lung recombinant cultures were approved by the IACUC of Cincinnati Children’s Medical Center.

For details of the methods used for immunostaining, mRNA in situ hybridization, and genomic DNA Southern blotting please see the supplemental procedures.