Conserved and Divergent Features of Mesenchymal Progenitor Cell Types within the Cortical Nephrogenic Niche of the Human and Mouse Kidney

Abstract

Cellular interactions among nephron, interstitial, and collecting duct progenitors drive mammalian kidney development. In mice, Six2⁺ nephron progenitor cells (NPCs) and Foxd1⁺ interstitial progenitor cells (IPCs) form largely distinct lineage compartments at the onset of metanephric kidney development. Here, we used the method for analyzing RNA following intracellular sorting (MARIS) approach, single-cell transcriptional profiling, in situ hybridization, and immunolabeling to characterize the presumptive NPC and IPC compartments of the developing human kidney. As in mice, each progenitor population adopts a stereotypical arrangement in the human nephron-forming niche: NPCs capped outgrowing ureteric branch tips, whereas IPCs were sandwiched between the NPCs and the renal capsule. Unlike mouse NPCs, human NPCs displayed a transcriptional profile that overlapped substantially with the IPC transcriptional profile, and key IPC determinants, including FOXD1, were readily detected within SIX2⁺ NPCs. Comparative gene expression profiling in human and mouse Six2/SIX2⁺ NPCs showed broad agreement between the species but also identified species-biased expression of some genes. Notably, some human NPC-enriched genes, including DAPL1 and COL9A2, are linked to human renal disease. We further explored the cellular diversity of mesenchymal cell types in the human nephrogenic niche through single-cell transcriptional profiling. Data analysis stratified NPCs into two main subpopulations and identified a third group of differentiating cells. These findings were confirmed by section in situ hybridization with novel human NPC markers predicted through the single-cell studies. This study provides a benchmark for the mesenchymal progenitors in the human nephrogenic niche and highlights species-variability in kidney developmental programs.

The self-renewal and commitment of NPCs is finely balanced by complex reciprocal signaling networks through Fgf, Bmp, Gdnf, Wnt, Notch, Fat4, and Hippo signaling pathways.^4–15 Their control of NPC fate, self-renewal, proliferation, and survival is directed by a number of transcription factors operating within NPCs, including Six2, Sall1, Osr1, Pax2, and Hox11 paralogs.^16–20 On induction, a subset of NPCs around each ureteric epithelial branch tip undergoes a mesenchymal-to-epithelial transformation forming an epithelial nephron precursor, the renal vesicle. The mouse kidney forms approximately 16,000 nephrons,²¹ all initiated over a 12–13-day period of development.²¹ In contrast, the final nephron count varies widely in the human kidney, with 1,000,000 a reasonable estimate, and all nephrogenesis initiated and completed in a period of 30–32 weeks.²² Interestingly, low nephron counts are linked to kidney disease.²² Cessation of mouse nephrogenesis is marked by the exhaustion of the nephron progenitor pool^21,23 and this is likely the case in the developing human kidney.

In conjunction with the commitment of NPCs, there is a progressive commitment of adjacent IPCs to different interstitial compartments in the mouse kidney.² This process is less well understood, although temporal fate mapping studies indicate a progressive recruitment of IPCs along the expanding radial axis of the kidney into distinct vascular-associated (pericytes and mesangial cells) and tubule-associated interstitial fibroblast populations in cortical and medullary regions.^2,24,25 At the molecular level, IPCs or their derivatives activate a number of transcriptional regulators that distinguish IPCs from adjacent NPCs: these include Foxd1, Meis1, and Pbx1.^10,26,27 The role of these is largely unclear, although genetic analysis shows Foxd1 and other IPC proteins play a key role in regulating interactions with NPCs and collecting duct precursor cells (CDPCs).^{9,10,14,28,29}

Our understanding of mammalian nephron and interstitial progenitor types is almost entirely on the basis of rat and mouse models. In the human kidney, initial studies have pointed to molecular differences between mouse and human NPCs in the expression of SIX family members, key factors in the specification and maintenance of NPCs.^16,30,31 Recent advances in generating kidney-like structures from human pluripotent stem cells highlight the need to understand human kidney progenitor types and their differentiated cellular derivatives to characterize and optimize in vitro strategies.^32–36 Here, we employed a variety of approaches to examine NPC and IPC compartments in the developing human fetal kidney. These data yield new insights into human kidney development and provide a valuable resource to guide in vitro efforts to engineer normal kidney structures.

Results

Differences and Similarities in Anchor Gene Expression Patterns in the Nephrogenic Zone

Mouse studies have identified Cited1 and Six2 as transcription factor–encoding genes expressed specifically by NPCs^3,37 and each is an anchor gene for the NPC compartment.³⁸ NPCs are surrounded by IPCs that in the mouse control NPC self-renewal and differentiation^9,14,29 and branching growth of the CDPC population.²⁸ Two well characterized transcriptional regulators identifying the mouse IPC compartment are Foxd1 and Meis1. Each is present in IPCs but not NPCs; however, Foxd1 is IPC specific within this lineage, whereas Meis1 extends into IPC derivatives outside of the nephrogenic zone.^2,26,39,40

We examined expression of human orthologs of these well characterized mouse NPC and IPC markers in the developing human kidney at weeks 14–15. As in the mouse, CITED1 and SIX2 were strongly expressed within mesenchymal cells capping the ureteric epithelial branch tips, the likely human NPC population (Figure 1, A and B). However, whereas Cited1 transcripts were restricted to NPCs in the mouse, CITED1 expression extended into differentiating pretubular aggregates in the human kidney. Further, Six2 RNA extends into early NPC derivatives, pretubular aggregates, and renal vesicles in the mouse,⁴¹ but in the human SIX2 expression was detected much later, within proximal regions of the S-shaped body (Supplemental Figure 1, D and E).

Figure 1.: In situ hybridization labeling for nephron compartment marker genes. (A-F) show expression for genes as indicated on fields. Left-hand and right column fields display in situ hybridization labeling of cryo-sectioned human week 14–15 kidneys. Sections show peripheral nephrogenic niches and interlobular nephrogenic niches (left and right, respectively). Red, blue, and black dashed lines indicate nascent nephrons, cap mesenchyme, and ureteric bud epithelium, respectively. PTA, pretubular aggregate; RV, renal vesicle; SSB, S-shaped body. Scale bar 50 μm.

Examining FOXD1 and MEIS1, we observed a FOXD1+/MEIS1+ population of peripheral mesenchymal cells similarly positioned to mouse IPCs, that are likely human IPC counterparts (Figure 1, C and D). Surprisingly, expression of both genes also extended into adjacent NPCs and early NPC derivatives, although expression of both genes was weaker in the NPC population (Figure 1, C and D). FOXD1 was also detected in podocytes consistent with a separate role for Foxd1 in podocyte programs from mouse kidney studies.¹⁰Sall1 and Wt1 encode zinc finger–containing transcription factors critical for kidney development expressed in both NPCs and IPCs in the mouse kidney with highest levels in the NPC population.^17,42,43 Human counterparts of both genes showed a mouse-like expression in the likely human NPC and IPC populations (Figure 1, E and F). In all material examined, no differences in gene expression were observed between peripheral and interlobular regions of the human kidney.

To determine whether overlapping gene expression profiles resulted in cotranslation of SIX2, CITED1, MEIS1, and FOXD1 mRNAs in NPCs, we performed immunolabeling studies on week 8 and 16 human kidneys comparing these data with E15.5 and P2 mouse kidneys. These developmental stages were chosen for reasons discussed previously⁴⁴ as they represent two stages of active nephrogenesis during and after ureteric branching.^21,23 In the mouse nephrogenic niche, Six⁺/Cited1⁺ cells cluster around Krt8⁺ ureteric epithelial branch tips (Figure 2A). High Six2 levels were observed in NPCs and Six2 was present at lower levels in anatomically distinct pretubular aggregates (Figure 2B), whereas Cited1 was restricted to NPCs, as predicted from in situ hybridization data (Figure 2C) and previous studies.⁴¹ In the human nephrogenic niche, SIX2⁺/CITED1⁺ cells were more broadly distributed around epithelial branch tips (Figure 2D), with a less marked difference in SIX2 levels in pretubular aggregates (Figure 2E), with detectable SIX2 and CITED1 extending into renal vesicles (Figure 2F and data not shown).

Figure 2.: Nephron and interstitial progenitor markers mix and persist into epithelializing nephrons. (A and D) Immunofluorescent stains for CITED1 and SIX2 in mouse and human kidneys. Insert in (D) shows CITED1 protein in the human RV (scale bar, 10 µm). (B, C, E, and F) Quantitative analyses of signal intensity distribution for CITED1 and SIX2. (G–L) Immunofluorescent analysis for FOXD1 and SIX2 and intensity correlation plots for these. White, blue, and red dashed lines indicate ureteric bud epithelium, cap mesenchyme, and nascent nephrons, respectively. Scale as indicated on fields. CM, cap mesenchyme; IM, interstitial mesenchyme; N, nephron; PTA, pretubular aggregate; RV, renal vesicle; UB, ureteric bud.

Analysis of Foxd1 showed Foxd1⁺ IPCs surrounding Six2⁺ NPCs in the developing mouse kidney; no Foxd1 was detected in the NPC population (Figure 2, G–I). At P2 Foxd1 was detected at very low levels around the nephrogenic niche (Supplemental Figure 1, F and G). In the human kidney, a strong FOXD1+ putative IPC population surrounded SIX2+ NPCs; however, FOXD1 was present in SIX2+ NPCs (Figure 2, J–L), albeit at lower levels (13% lower than in IPCs). MEIS1 was also detected within human NPCs and mouse NPCs showed low levels of Meis1 at both E15.5 and P2 (Supplemental Figure 1, A–C). In summary, human and mouse kidneys differ in the broader extent of coactivation of IPC-associated regulatory factors within the NPC population and the persistence of NPC-associated regulatory factors into differentiating nephron components.

RNA Sequence Analysis of Purified Human and Mouse Nephron and Interstitial Progenitors Suggests Divergences and Similarities in Regulatory Pathways

Previous studies have attempted to obtain transcriptional profiles of human NPCs utilizing an ITGA8-directed antibody enrichment protocol to compare transcriptional profiles between mouse and human NPC compartments. These approaches identified the transcription factor SIX1 as a specific component of the human NPC population during periods of active nephrogenesis. These NPC-enriched populations showed significant contamination from other mesenchymal cell types. To obtain a more specific human SIX2+ NPC profile, we applied the MARIS (method for analyzing RNA following intracellular sorting) approach.⁴⁵ In this, a cortical mesenchymal kidney isolate was fixed and permeabilized, then immunostained with anti-Six2/SIX2 antibodies. Mouse (Six2+) and human (SIX2+) NPC-enriched populations were purified by FACS, mRNA isolated, and RNA sequencing performed to obtain NPC expression profiles.

Initially, we performed control experiments using a transgenic Six2GFP reporter mouse strain which labels Six2+ nephron progenitors with nuclear GFP.³ E16.5 Six2GFP-expressing kidneys were gently dissociated to release cells from the cortical nephrogenic niche and Six2GFP+ and Six2GFP− cells were either directly isolated by FACS or subjected to MARIS, then transcriptionally profiled (Figure 3A). Mouse cells processed for MARIS with a Six2 antibody are hereafter referred to as mMARIS-Six2+ or mMARIS-Six2−. The two isolation techniques were similar in cell content: approximately half of the cortical cell preparation (56%) were Six2-GFP+ and a similar fraction (66%) were mMARIS-Six2+. Six2GFP+ and mMARIS-Six2+ RNA-seq datasets showed a strong correlation (R²=0.97), indicating that MARIS generates a comparable transcriptional profile to FACS isolation of viable, GFP-labeled Six2+ cells (Figure 3B).

Figure 3.: Transcriptional profiling of mouse and putative human NPCs assisted by intracellular staining of Six2/SIX2 followed by FACS (MARIS). (A) Separation of Six2+ cell population from dissociated mouse (m) embryonic kidney cortex cells by either FACS of Six2GFP reporter line (middle) or Six2 MARIS (right). (B) Gene-level correlation of normalized mRNA-Seq reads between NPC profiles generated by Six2 reporter line (Six2GFP+) and Six2 MARIS (mSIX2+). (C) Overlap (left) between NPC-specific genes identified by differential gene expression analysis between Six2GFP+ versus Six2GFP− (TPM_Six2GFP+>5, TPM_Six2GFP+/TPM_Six2GFP−>3, P<0.05), or between mSix2+ versus mSix2− (TPM_mSix2+>5, TPM_mSix2+/TPM_Six2−>3, P<0.05). Results (middle) of GO term enrichment analysis of the indicated gene sets, with representative ones (right) from each set of genes. (D) Separation of SIX2+ cell population from dissociated human (hu) fetal kidney cortex. (E) Gene-level correlation of normalized mRNA-seq reads between human and mouse NPC profiles obtained by MARIS; human (orange) or mouse (cyan) enriched genes were indicated. (F) Top three GO terms enriched from the human (top) and mouse (bottom) enriched genes. Ab, antibody; enrich., enriched; norm., normalized vs., versus.

As expected, both the Six2GFP+ and mMARIS-Six2+ cells expressed nephron progenitor markers such as Phf19, Cited1, Osr1, and Six2 (Figure 3C, Supplemental Table 1). The key difference between these two approaches was in the detection of early nephron induction markers such as Wnt4 and Pax8: these were weakly expressed in Six2-GFP+ cells but not mMARIS-Six2+ NPCs. Thus, mMARIS-Six2+ sorted cells likely displayed a more progenitor-like profile, presumably a reflection of the bias in setting a window for Six2 detection in the MARIS that selects for the higher Six2 levels relative to live FACs of Six2-GFP+ NPCs. In addition, perdurance of GFP in Six2-GFP+ NPCs may capture a small population of the earliest induced NPCs. As expected, gene ontology (GO)–term analyses showed mMARIS-Six2 and Six2-GFP samples were enriched for genes associated with kidney development. Collectively, these data demonstrate that the MARIS strategy can generate a robust NPC transcriptional signature.

To profile SIX2+ human NPCs, we performed a brief cortical dissociation of human fetal kidneys (week 16) to release mesenchymal cell types, then performed MARIS to isolate SIX2+ cells (Figure 3D). Approximately 70% of the human cortical cell population was positive for SIX2, a comparable number to mouse cortical isolations. An initial comparison of the hMARIS-SIX2+ RNA-seq data to the previously generated ITGA8+ NPC-enriched cell profile³⁰ showed a good correlation (R²=0.81; Supplemental Figure 2A). Further analysis showed the hMARIS-SIX2+ sample displayed a higher expression of nephron progenitor markers (Supplemental Figure 2C) and lower expression of genes expressed by differentiating cells, except for PAX8, or in epithelializing nephron structures (Supplemental Figure 2D). The hMARIS-SIX2+ RNA profile also showed minimal contamination by blood and vascular endothelial cell types (Supplemental Figure 2E), and reduced expression of genes indicative of cells within the ureteric epithelium that underlies the NPC niche (Supplemental Figure 2F). In conclusion, the hMARIS-SIX2+ RNA-seq profile matched expectations for a highly enriched NPC population.

Overall, hMARIS-SIX2+ and mMARIS-Six2+ expression profiles showed a significant correlation in their gene expression profiles (R²=0.61) (Figure 3E). Genes common to both mouse and human nephron progenitors included known progenitor markers such as SIX2, CITED1, PHF19, OSR1, SALL1, and EYA1. Several genes expressed in nephron progenitors in both species still displayed variations in absolute expression levels; for example, CITED1 was expressed at lower levels compared with Cited1 (transcripts per million [TPM] values of 54 human versus 259 mouse) and, as expected, SIX1 showed an extreme difference³⁰; Six1 was undetectable in mouse NPCs (TPM values of 27 human versus 0 mouse).

To specifically identify genes with differential expression profiles between human and mouse NPCs, we looked for genes expressed at levels >5 TPM with an expression level differential of three-fold or greater. By these criteria, 1230 genes were enriched in human NPCs and 1087 genes in mouse NPCs (see Supplemental Table 2). GO-term analyses on the mouse-enriched genes indicated a strong upregulation of genes involved in oxidative phosphorylation and mitochondrial function, and the regulation of cell proliferation (Figure 3F). Other genes strongly enriched in human progenitors included components of retinoic acid signaling CRABP2 (TPM values of 349 human versus 54 mouse), cell-adhesion complexes CDH24 (TPM values 284 human versus 9 mouse), and genes linked to human disease and congenital disorders: DAPL1, linked to renal neoplasia (Klomp et al.⁴⁶; DAPL1 TPM values 155 human versus 0 mouse), and COL9A2, linked to Stickler syndrome, which is associated with renal agenesis (Baker et al.⁴⁷; COL9A2: TPM values 139 human versus 2 mouse).

Conversely, several genes were more strongly represented in mouse NPCs, including: Crym, a previously identified anchor gene for mouse NPCs that is not detected in human NPCs (Rumballe et al.²³; TPM values of 205 mouse versus 0 TPM in human), and Capn6 (TPM values of 248 mouse versus two human). In situ data and data from a transgenic mouse model confirm mouse NPC activity of Crym (GUDMAP ID 22105 and 14077) and Capn6 expression in cap mesenchyme and nascent nephrons has been documented in mouse kidney studies.⁴⁸ Variation was also observed in related genes with potential overlapping activities. As an example, Rspo1, a modulator and agonist of WNT-signaling,⁴⁹ was expressed at six-fold higher levels compared with RSPO1 (TPM values of 84 mouse versus 13 human); in contrast, RSPO3 and Rspo3 displayed comparable levels (TPM values of 20 human versus 19 mouse). A redundancy between Rspo1 and Rspo3 actions could underlie the absence of a phenotype in Rspo1 mutants.^50,51

To validate expression predictions from the MARIS data, we selected 17 genes enriched in either human or mouse NPCs for further characterization by SISH (Figures 4 and 5, Supplemental Table 3). Of these, we could detect expression for 16 in either the mouse or human kidney, or both. The exception was Fgf20/FGF20, which is predicted to be mouse NPC specific, but was below the limits of SISH detection (TPM values of 10 mouse versus 0 human). Matching the differential expression predictions, HIP1R, UNC5B, LYPD1, DAPL1, ECEL1, COL9A2, WASF3, and TNFRSF19 were expressed at high levels in human NPCs, and Crym, Serpinf1, Slc12a2, Foxd2, Rspo1, and Capn6 in mouse NPCs (Figures 4 and 5). Interestingly, RSPO1, CRYM, CAPN6, and SLC12A2, four human homologs of mouse-specific nephron progenitor–enriched genes, were actually expressed in developing human nephrons at the S-shaped body stage, a profile not observed in the developing mouse kidney (Figure 4, I, J, L, and M). Conversely, mouse homologs of human-specific nephron progenitor–enriched genes Unc5b, Col9a2, and Pcdh15 were also expressed in other cellular compartments of the mouse kidney (Figure 5, D, G, and K). In summary, SISH analysis verified predicted species-specific differences in expression profiles, indicating that other predicted differences are likely to be valid, and identified additional species-specific differences in the expression of this cohort of genes in other kidney structures.

Figure 4.: In situ hybridization labeling for human and mouse enriched nephron progenitor genes. (A–O) In situ hybridization labeling of cryo-sectioned human week 16 kidneys. In situ labeling as indicated on fields. Inserts show enlarged regions from main fields. Scale bars as indicated on fields.

Figure 5.: In situ hybridization labeling for human and mouse enriched nephron progenitor genes. (A–O) Complementary in situ hybridization labeling of cryo-sectioned mouse E15.5 and P2 mouse kidneys complementary to Figure 4. In situ labeling as indicated on fields. Scale bars as indicated on fields. Scale bars in magnified inserts are 20 µm.

Data in Figures 1 and 2 indicate that human orthologs of mouse genes distinguishing mouse IPCs from mouse NPCs were also expressed within the human NPC population. To evaluate if this is a wider trend and to compare human and mouse IPC gene expression profiles we adopted a two-pronged approach. First, we extracted a mouse IPC-enriched gene expression profile by isolating mouse IPCs using the Foxd1-GCE strain, where Cre-ERT2 is expressed from the Foxd1 locus^2,24 in combination with the Rosa26tdTomato reporter line.⁵² Pregnant mice were injected with tamoxifen at E13.5 and Tomato+ cells were isolated at E14.5 for RNA sequencing. To be able to compare this expression profile to mouse NPCs we isolated mouse NPCs using the Six2-GFP reporter strain as described above. Second, to identify a human IPC RNA profile we performed MARIS colabeling with MEIS1 (FOXD1 antibodies were not compatible with this procedure) and SIX2 antibodies on preparations of cortical, mesenchyme cell–enriched human kidney isolates at weeks 13–15 (Supplemental Figure 3A), generating RNA-seq profiles for MEIS1+/SIX2− (IPC-enriched), MEIS1+/SIX2+ (NPCs), and cortex cells (Supplemental Figure 3B).

We first contrasted genes enriched in either human or mouse IPCs (hIPC and mIPC) to their respective cortex RNA profiles (Supplemental Figure 3C). hIPCs and mIPCs enriched genes showed a low correlation (R=0.36) as expected because of the broader distribution of MEIS1 compared to Foxd1. We focused the analysis to identify genes differentially expressed between human IPCs and NPCs (Figure 6A). As anticipated, the human IPC and NPC fractions both expressed FOXD1 and MEIS1 (TPM 38 versus 41 and 120 versus 110, respectively), whereas SIX2 and CITED1 where confined to NPCs (TPM 2 versus 191 and 1 versus 49).

Figure 6.: Transcriptional profiling of human and mouse IPCs. (A) Gene-level correlation of normalized mRNA-seq reads between human (hu) IPC and NPC. (B) Gene-level correlation of normalized mRNA-seq reads between mouse (ms) IPC and NPC. (C) Immunostaining of interstitial markers in mouse and human kidneys as specified on fields. (D) Breakdown of mouse (top) or human (bottom) genes expressed in IPC or NPC by their relative expression in the two cell types. Genes enriched in one of the cell types satisfy TPM>5 and fold change >3. Other expressed genes are categorized as “non-DE.” (Middle) Pie chart shows breakdown of mouse IPC-enriched genes by their relative expression between human IPC and NPC. DE, differentially expressed enrich., enriched; norm., normalized; UB, ureteric bud epithelium; vs., versus.

Applying similar thresholding criteria as in earlier MARIS data, 503 genes showed enriched expression in hIPCs versus hNPC cell fractions (Figure 6A, Supplemental Table 4), including genes associated with extracellular matrix or matrix interactions (ITGA9, ITGA1, COL3A1), transcription (GATA3), and cell signaling (PDGFRB). Conversely, 534 genes were specifically enriched in NPC versus IPC fractions. These included well characterized NPC marker genes such as SIX2, CITED1, and EYA1; genes identified earlier in SIX2 MARIS comparisons such as PCDH15, LYPD1, and ECEL1; and novel gene predictions including ELAVL4, FAT3, and CRABP2. Enrichment of PDGFRA, PDGFRB, and PBX1 within hIPCs was confirmed through immunolabeling studies (Figure 6C). In agreement with expression profiles, PDGFRA and PDGRFB were only detected in IPCs (expression extends also into likely IPC interstitial/stromal derivatives) whereas PBX1 was present in both hIPCs and hNPCs but at markedly elevated levels within hIPCs (Figure 6C).

To determine whether human orthologs of mouse IPC markers were expressed more broadly in NPCs, we first identified a full set of genes whose expression was enriched in mouse mIPCs compared with mNPCs (Figure 6B, Supplemental Table 5); this gave a set of 647 genes, including Foxd1 and Meis1, and other recognizable interstitial markers. We next determined whether the human orthologs of these 706 genes were expressed in hIPC and hNPC enriched fractions, and the relative expression between each cell population (Figure 6D, Supplemental Table 6). Twenty-seven percent of human orthologs were enriched in hIPCs displaying a similar expression to the mouse. Twenty-two percent of genes were coexpressed in both hNPCs and hIPCs including MEIS1, FODX1, SMOC2, and ROR2. Smoc2 is expressed broadly in the cortical nephrogenic interstitium of the mouse but is largely absent from NPCs.⁵³SMOC2 was expressed at similar levels in hNPC and hIPC cell fractions (TPM 10 versus 14, respectively). In contrast, Ror2 expression is linked to mouse NPCs⁵⁴ and, although ROR2 was found at higher levels in hNPCs, ROR2 transcripts were also present in IPCs (TPM 29 versus 12, respectively). Three percent of the mIPC markers were not expressed in hIPCs but were expressed in hNPCs, having potentially shifted expression from the interstitial to nephrogenic lineage (e.g., CRABP2), whereas 49% of genes were not expressed above the cutoff threshold (TPM 5) in the human IPC fraction. Collectively, the data indicate a significant disparity between the transcription profiles of human and mouse IPCs.

Cellular Diversity of Human Nephron Progenitors

Single-cell RNA sequencing (scRNAseq) can potentially reveal cellular heterogeneity that is difficult to define with other procedures. We applied scRNAseq using the 10× Genomics platform⁵⁵ to profile 2750 predominantly mesenchymal cell types from the cortical nephrogenic niche of the week 16 fetal human kidney. Twelve cell populations emerged from unsupervised clustering analyses using Seurat⁵⁶ (Figure 7A). Clusters were identified by known marker genes for each population (Supplemental Table 7). Five cell population clusters belonged to the interstitial lineage, whereas the nephrogenic lineage was represented by three cell populations; the remaining populations included vascular endothelial cells, distinct proliferating cell compartments, and cells of the immune system (Figure 7, A and B, Supplemental Table 7).

Figure 7.: Single-cell transcriptional profiling of human nephrogenic niche cells. (A) tSNE plot displaying principal component analysis of approximately 2800 human kidney cortex cells from week 16 kidney. Cell identities as indicated on figure by gene expression. Dashed line demarks the nephron progenitors. (B) Cluster hierarchies inferred from differential gene expression and GO-term analyses of top 50 differentially expressed genes per cluster. (C) tSNE plot displaying principal component analysis of cells from cluster 4 in (A). (D) Cluster hierarchies inferred from differential gene expression and GO-term analysis of top 50 differentially expressed genes per cluster as seen in (C). (E) tSNE plots displaying gene expression levels in cells. (F–J) Gene expression plots for novel and established NPC markers. Genes as indicated on plots. tSNE, t-Distributed Stochastic Neighbor Embedding.

Here, we focused on the NPC compartment of these data. NPCs were identified as cell population 4 on the basis of CITED1, SIX1, LYPD1, and DAPL1 expression (Figure 7A). To scrutinize cellular diversity within cell population 4, we re-examined this cluster (Figure 7C, Supplemental Table 8). Four NPC subclusters emerged, which segregated into four distinguishable cell populations; we termed these 4A, 4B, 4C, and 4D (Figure 7C). NPCs (4A) expressed TMEM100, CITED1, and MEOX1. A second cell population (4B) differentially expressed ID1, MEG3, and DAPL1, and induced NPCs (4C) expressed CCND1, PAX8, and LHX1, and contained a small subpopulation of differentiating MAFB and PODXL expressing cells likely initiating podocyte differentiation (Figure 7J). Proliferating cells (4D) expressed CENPF, MKI67, and TOP2A (Figure 7, C–E and I) and likely represented a mixture of several cell types aggregated by their shared strong cell-cycle profile, because TOP2A/MKI67/CENP- expressing cells included subsets of cells expressing PAX8, SIX1, and CITED1 (Figure 7, E–I).

In the progression of mouse nephrogenesis, Cited1+/Six2+ self-renewing nephron progenitors transition into Cited1−/Six2+ cells, a cell state primed for differentiation, then to an induced, committed Cited1−/Six2−/Pax8+ nephron-forming cell state.^5,41 The expression domains for CITED1, SIX2, and SIX1 recapitulated that expected for self-renewing, primed, and committed NPCs (Figure 8A) and their expression profiles (Figure 7F) suggest they closely correspond to the cell populations identified as 4A to 4C, which define overlapping expression domains. In situ hybridization showed that COL9A2, PCDH15, UNC5B, and ECEL1, as identified earlier (Figures 3 and 4), were mainly expressed in 4A cells, whereas WASF3, DAPL1, PHF19, and TNFRSF19 were expressed also in 4B, and LYPD1 most strongly in 4C (Figures 8A and 7, F–H). Each domain showed significant overlap. TMEM100 and ROBO2 were predicted to be coexpressed in population 4A and follow a CITED1-like pattern (Supplemental Table 8). TMEM100 and ROBO2 displayed a restricted localization within NPCs, but unlike CITED1 their expression did not persist into early nephron-forming stages (Figure 8C), although much later in nephron development ROBO2 was upregulated in podocytes. ROBO2 was also present in the interstitial lineage.

Figure 8.: Validation of NPC cell populations and exploration of novel NPC marker genes. (A) Genes identified in NPC subclusters stratify into distinct gene expression patterns; section in situ hybridization for genes as specified on fields, clusters as specified. (B) Gene expression plot for TMEM100, ROBO2, and CITED1. (C) Immunofluorescent staining for TMEM100, CITED1, ROBO2, and KRT8/19 in human fetal kidney. Scale bar as indicated. RV, renal vescile UB, ureteric bud epithelium

Discussion

Here, we examined the conserved and divergent features of the human and mouse nephrogenic niche using single-cell sequencing, MARIS sequencing, RNA sequencing, in situ hybridization, and immunohistochemistry. In summary, we find that cells in the human nephrogenic niche show significant divergence from their mouse counterparts and that the boundaries between NPC and IPC lineages follow different rules to those defined in the mouse. We focus our discussion to the differences and similarities of human and mouse NPCs and IPCs and the effect these may have on the nephrogenic niche.

Comparative Analysis of Human and Mouse NPCs

A number of genes have been identified genetically as having important roles within the NPC compartment of the developing mouse kidney. These include Six2, Eya1, Osr1, Gas1, Itga8, Fgf20, and Pax2, which are expressed within NPCs and not IPCs; and Wt1 and Sall1, which are expressed in both NPCs and IPCs, but at elevated levels in the former.^{4,16–18,43,57–60} Human equivalents of these genes showed broad conservation consistent with conserved roles from mouse to man. Interestingly, we observe conservation in gene expression profiles for other highly NPC restricted genes, such as Cited1 and Phf19, that have no observable function in the mouse kidney (J.A. McMahon, unpublished data),⁶¹ indicating conservation in regulatory programs that do not appear to underlie a functional role in all mammalian species. Whether there is a distinct role for either gene in human NPCs remains to be determined.

Although several functionally important genes showed conserved expression, our comparative analyses of human and mouse NPCs highlighted a large number of genes (1230 and 1087, respectively) enriched in NPCs of each species. We confirmed the predictions hold for all 16 genes detectable by SISH, indicating that there are likely to be many more genes with bona fide expression differences in these datasets. These findings beg the question, what distinct biologic processes might be at play within mouse and human NPCs?

Two critical processes are the regulation of progenitor self-renewal and differentiation, the balance of which ultimately determines the final number of nephrons formed. In the mouse, both require Wnt9b signaling in NPCs; the Wnt9b ligand is secreted by the subjacent ureteric epithelium.^6,62 Several mouse Wnt9b target genes have been suggested, including Cited1, Btbd11, Etv5, Gdnf, and Itga8⁶; many of these are functionally important in the mouse NPC.^13,59,63 Strikingly, human NPCs displayed multiple examples of putative WNT9B targets that could not be detected, such as CDH4, SLC45A3, SORBS2, CLDN9, PLA2G7, and SLC12A2, suggesting differences in regulatory mechanisms between mouse and man.

GO analyses of species-enriched genes in mouse and human NPCs suggested increased oxidative metabolism in mouse. Consistent with these data, basic cellular processes, such as metabolism, formulate a large portion of interspecies gene expression disparity.⁶⁴ Although a functional role for metabolism in NPCs has not been directly addressed, deletion of p53 in mouse NPCs results in aberrant cell metabolism and reduced NPC numbers.⁶⁵ A temporally condensed nephrogenic program, as seen in the mouse compared with human (Lindström et al.⁴⁴), may increase metabolic demands on NPCs.

Of note, during the validation of gene expression differences, we found Rspo1, Crym, Capn6, and Slc12a2 to be enriched in mouse but not in human NPCs. RSPO1, CRYM, CAPN6, and SLC12A2 were instead specifically expressed in narrow segments within the human S-shaped body. The transposition of gene expression from one compartment to another suggests either the loss of a requirement or the necessity of a different cell type. Alternatively, these genes are simply a readout of broader underlying changes to pathways and may represent nonessential genes that fluctuate in expression without functional consequences. The exact nature of all of these differences requires further investigation to determine if they relate to biologic function.

Gene Expression in Human Nephron and Interstitial Progenitors Adheres to Different Rules from that in the Mouse

Our data suggested that human orthologs of mouse IPC marker genes (Foxd1/FOXD1 and Meis1/MEIS1) were not restricted to the IPC lineage, as they are in the mouse, but were also expressed in human NPCs. To examine this in greater detail we performed RNA profiling on mouse and human IPCs and NPCs and demonstrated that only 27% of genes that we categorized as mouse IPC marker genes were enriched in human IPCs. The remaining mouse IPC marker genes (73%, 472 genes) would therefore be categorized as displaying nonconserved expression patterns.

It has recently been shown that mouse NPCs maintain their identity and prevent lineage switching to interstitial cell types via a Pax2-dependent mechanism.²⁰ When Pax2 is genetically removed from mouse NPCs, they upregulate IPC-enriched genes such as Col1a1, Col3a1, Col1a2, Anax2, and Dcn.²⁰ Given that Pax2 is implicated in interstitial cell fate repression, we scrutinized the expression of PAX2 to determine if it is downregulated in human NPCs. However, human PAX2 was robustly expressed in NPCs (NPC TPM values of 364 human versus 159 mouse). COL1A1 and ANXA2 were not expressed in either NPCs or IPCs in the human but we did find expression of COL3A1, COL2A2, and DCN in IPCs and not NPCs, as expected from the mouse. It is therefore unlikely that the increased expression of IPC genes in human NPCs is a result of a PAX2-dependent mechanism.

The change in FOXD1 expression is noteworthy due to the known function of Foxd1 in regulating, directly or indirectly, the patterning of the kidney capsule, collecting duct, and nephron.^10,26FOXD1 mRNA and protein levels were very similar comparing human NPCs and IPCs, suggesting FOXD1 within NPCs could fundamentally alter signaling dynamics within the nephrogenic niche. In mouse IPC’s Foxd1 is required to regulate both Dcn and Fat4, which encode a leucine-rich proteoglycan and a membrane-bound signaling factor, respectively.^9,66Foxd1 and Fat4 mutants both display an expansion of NPCs due to a failure of NPC commitment, a phenotype that closely resembles a gross ablation of the IPC compartment.^9,14,29 Whereas DCN expression remains enriched in IPCs compared with NPCs as in the mouse (TPM IPCs: 11 versus 45; TPM NPCs: 4 versus 2), interestingly, human FAT4 is expressed at much lower levels than its mouse counterpart in IPCs (TPM values of 4 versus 21), suggesting a potential alteration in the FAT4 signaling axis. Because reduced FAT4 signaling is predicted to enhance progenitor expansion, such a mechanism could contribute to a larger, longer-lived human NPC population.

Single-Cell Analyses Reveal Population Complexities in the Human Cap Mesenchyme

scRNAseq is likely to play an important role in defining cell diversity and providing evidence for diversity-generating processes in human kidney development where genetic approaches, a mainstay of mouse studies, are not possible. A current mouse-centered model of NPC differentiation suggests NPCs differentiate from a self-renewing Cited1+ state and progress through an intermediate state primed for differentiation in response to a combination of Wnt, Bmp, and Hippo signaling.^5,9,62 The human NPC population displays comparable diversity as judged by single-cell transcriptional profiling. We validated predicted cell clusters by examining the expression of 14 genes by in situ hybridization (Figure 8) and showed that gene expression patterns can be categorized within the expression domains as defined by CITED1, SIX2, and SIX1. Our data agree with scRNAseq data from the mouse which indicates that the NPC population displays low diversity.⁶⁷ The focus is now on determining the pathways and genes that control the differentiation cascade during induction and the differentiation trajectories that generate specific cell states within developing nephron precursors.

Concise Methods

Complete methods are included as Supplemental Material. Further information as pertinent to in situ hybridization, immunolabeling, and microscopy are as described previously in this series of papers (Lindström et al.⁴⁴). Here, we elaborate in detail on methods relating to MARIS, single-cell sequencing, and image quantification as pertinent to the work described here. We also included details for the specimens used.

Human Kidney Specimens

Consented, anonymized, human fetal tissue was obtained from elective terminations following review of the study by Keck School of Medicine of the University of Southern California’s Institutional Review Board. Kidney samples ranging in age from 13 to 18 weeks of gestation were supplied from the Children’s Hospital of Los Angeles and the University of California, San Francisco. Gestational age was determined per guidelines specified by the American College of Obstetricians and Gynecologists using ultrasound, heel-to-toe, and crown-to-rump measurements following published Carnegie Stages.^68–70 Stages as stated in the manuscript indicate age of embryo or fetus from point of conception/fertilization. Samples from the Children’s Hospital of Los Angeles were received immediately after elective terminations and transported on ice at 4°C in 10% FBS, 25 mM Hepes, high-glucose DMEM (SIGMA). Samples from the University of California, San Francisco, were transported similarly but by overnight courier. Given the anonymized nature of the specimens, no further information was available regarding the specimens or the normalcy of the pregnancy.

Animals

All animal work was reviewed and institutionally approved by Institutional Animal Care and Use Committees at the University of Southern California and performed according to institutional guidelines. Timed matings were set up to recover embryos and neonates at the appropriate age. The Foxd1-GCE strain (B6;129S4-Foxd1tm2(GFP/cre/ERT2)Amc/J) was generated as previously described (Humphreys et al.²⁴). The Rosa26tdTomato reporter line (B6.Cg-Gt(ROSA)26Sortm14(CAG-tdTomato)Hze/J) was obtained from JAX.⁵² Heterozygous Foxd1-GCE animals were crossed with female Rosa26tdTomato homozygous females. Pregnant females were injected with 3 mg Tamoxifen per 40 g at E13.3 and kidneys collected at E14.5. The Six2TGC line was generated as previously described.³ Male heterozygous Six2TGCs animals were crossed with female Swiss Webster mice and embryos collected at E16.5.

Image and Sample Quantification

2D Immunofluorescent Analyses

Frozen and sectioned samples were stained as described previously (Lindström et al.⁴⁴) to detect SIX2, SIX2, LEF1, FOXD1, and CITED1. KRT8 and β-laminin were used as structural markers to determine the location of the ureteric epithelium and nephrons. Images were captured with a 63× objective on a Leica SP8. Data were captured as 8-bit images. IMARIS 8.2 (Bitplane) was used for quantification of nuclear antibody signals. The Spot function was used to manually add circular spots to mark all nuclei on the image frame using DAPI-highlighted nuclei as a reference. Because of the convoluted shape of nuclei in 2D sections and the circular shape of the Spot function, we used multiple smaller spots to represent single nuclei to ensure accurate quantitation and coverage across the nuclei. Spots were grouped into three cell populations: (1) cap mesenchyme (SIX2⁺ cells), (2) cortical interstitium (FOXD1⁺ cells), and (3) all other cells. To compare the mean intensity of spots we first normalized the mean intensity values, taking into consideration the background signal and the maximum signal for each channel. To do this we measured the intensity for each channel throughout all spots and identified the lower fifth percentile intensity mean (background), as well as the maximum value. Each spot’s intensity mean was thereafter normalized as follows:This transforms the intensity mean of each spot onto a 0–100 scale with the fifth percentile equaling 0 and the maximum being 100, respectively. To plot the normalized intensity of spots against their position within the cap mesenchyme population we marked the most cortical point of the cap mesenchyme and utilized this as point 0. A line was extended in a medullary direction parallel to the ureteric epithelium around which the cap mesenchyme was located.

RNA Sequencing Data

All RNA sequencing data are provided at the Gene Expression Omnibus; GEO accession numbers: GSE102378 (mouse RNA-seq data), GSE102230 (human RNA-seq data), and GSE102596 (human single-cell RNA-seq data). Full details for the number of samples can be found for each submission. In brief, MARIS or conventional RNA sequencing was performed on five kidneys for SIX2 MARIS, two kidneys for MEIS1/SIX2 MARIS, three kidneys for mouse Six2 MARIS, three kidneys for mouse Six2-GFP sequencing, and two kidneys for mouse Foxd1 sequencing. One kidney was dissociated for scRNAseq. The GEO submission comprises 39 RNA-seq libraries for RNA-seq and one multicell library for scRNAseq.

MARIS Staining and FACS

The MARIS staining and FACS procedure was performed as described in Hrvatin et al.⁴⁵ with the following modifications. Human and mouse cortical nephrogenic zone cells were digested from E16.5 embryonic mouse kidneys or 13–18-week fetal human kidneys using 10 mg/ml pancreatin (P1625; Sigma) and 2.5 mg/ml collagenase A (11 088 793 001; Roche) enzyme mixture and filtered through a 40-µm filter (352340; BD Falcon) as described in.⁷¹ Cell fixation, washing, permeabilization, and centrifugation were performed as described in Hrvatin et al.⁴⁵ using the following solutions with all subsequent steps performed at 4°C. Fix buffer: 4% PFA (Electron Microscopy Sciences), 0.1% saponin (47036; Sigma-Aldrich) in molecular grade PBS (Ambion) supplemented with 1∶100 RNasin Plus RNase Inhibitor (N2615; Promega). Wash buffer: PBS containing 0.2% BSA (Gemini Bio-Products), 0.1% saponin, 1∶100 RNasin Plus RNase Inihibitor. SIX2 Primary antibody (MBS610128; Mybiosource) and MEIS1/2/3 Primary antibody (39795; Active Motif) staining of cells at 1:5000 dilution was carried out while rocking overnight at 4°C in staining buffer containing PBS with 1% BSA, 0.1% saponin, and 1∶25 RNasin Plus RNase Inhibitor. Cells were washed and stained with donkey anti-rabbit Alexa 488 (A-21206; Thermofischer) and goat anti-mouse IgG1 555 (A21127; Thermofischer) secondary antibody for 45 minutes. Subsequent washing and FACS sorting were performed at a concentration of 5–10 M cells/ml with sort buffer containing PBS, 0.5% BSA, and 1∶25 RNasin Plus RNase Inhibitor. Cells were sorted on the FACSAria I and II (BD Biosciences) using FACSDiva software. Sorting gates were set with reference to negative controls with no primary antibody stain. The sorting efficiency was maintained at >90%. Cells were collected in tubes that were coated with a small amount of sorting buffer.

For FACS sorting of Six2TGC and TdTomato positive cells, the mouse kidneys were dissected and dissociated in the same enzymatic solution as described above. The cells were not fixed but instead immediately FAC sorted for GFP or tdTomato.

RNA Isolation of MARIS

The RNA isolation was performed as described in Hrvatin et al.⁴⁵ with the following modifications. FACS-collected SIX2+, SIX2−, MEIS1+SIX2+, MEIS1+SIX2−, and MEIS1−SIX2− cells were pelleted by centrifugation at 3000 × g for 10 minutes at 4°C. Total RNA was isolated using the RecoverAll Total Nucleic Acid Isolation kit (AM1975; Ambion), starting at the protease digestion step with protease incubation time of 1 hour at 50°C, and inactivated at 80°C for 15 minutes. Cell lysates were frozen at −80°C overnight and extracted for RNA according to the manufacturers recommended protocol.

For the nonfixed mouse cells, the RNA was isolated as previously described (Lindström et al.⁴⁴).

RNA-Seq Analysis

mRNA-seq libraries were synthesized with Kapa stranded mRNA-Seq kit, and were sequenced on an Illumina NextSeq500 platform at USC Epigenome Center. All mRNA-seq reads were aligned to the mouse or human reference genome (mm10 or hg38) using TopHat2.⁷² Quantification of RNA-seq reads to generate RPKM was performed by Partek Genomics Suite software, version 6.6 (St. Louis, MO). TPM was calculated by dividing RPKM value by ratio of sequencing reads from the corresponding library that were mapped to exon regions of the genome. We identified differentially expressed genes as those satisfying the following three criteria: (1) P value <0.05 from statistical tests performed by DESeq2⁷³; (2) >3-fold difference of average normalized read counts between the groups compared; and (3) average TPM>5 in at least one of the groups. GO analysis was performed using PANTHER classification system⁷⁴ (http://pantherdb.org/). We ranked the relevance of GO terms by fold enrichment of number of observed genes over number of expected genes. The GO terms with binomial P value >0.01 were omitted due to statistical insignificance. We analyzed the variability between all of the MARIS RNA-seq data for NPCs from both the huSIX2+MARIS and the huSIX2+/MEIS1+ MARIS data and found that correlation between samples was high in six of seven samples (R² range 0.93–0.98). Replicate 1 from the huSIX2+MARIS displayed lower correlation to the other samples (R²=0.72–0.78). This variability may have arisen at various points: (1) each replicate RNA sample was extracted from a different human fetal kidney with no known, but presumed, genetic variability, in addition to samples originating from a range of close developmental stages; (2) RNA from replicates 2–5 for the huSIX2+MARIS exhibited lower quality, as measured by RNA integrity, due to difficulties in library construction/sequencing consequent to low RNA content, thereby indicating that replicate 1 may be higher quality; and (3) the total amount of mapped reads from replicates 2–5 is approximately 25% less than replicate 1, which might have contributed to decreased sample complexity.

Single-Cell Sequencing

Cell Preparation

Cells were dissociated as described for the MARIS protocol from a week 16 kidney and live cells sorted by FACS using DAPI (Thermo Fisher Scientific) and DRAQ5 (Thermo Fisher Scientific) to select against dead cells and for live cells, respectively. Seventy-eight percent of cells were live and intact, indicating robust isolation methods. Seven thousand live cells were input into a 10× Genomics Chromium device expecting the capture of 4000 cells. Illumina ready sequenceable libraries were then generated using the 10× Chromium single cell 3′ RNA-seq protocol. Subsequently, sequencing was carried out on the Illumina NextSEquation 500/550 platform with the goal of obtaining at least 50,000 reads per cell. Three thousand seven hundred thirty-one valid barcodes (“cells”) were recovered after filtering and UMI counting.

Sequence Mapping

Mapping was performed using the CellRanger software version 1.3.1 through the CellRanger count command. We used STAR version 2.5.1b to map the second end of the FASTQ reads to the human genome version GRCh37.p13 and uniquely mapped reads were counted using the Ensembl GTF annotation as reference. A total of 70.1% of the reads had unique mapping, which corresponded to a total of 3731 valid barcodes.

Quality Control

To filter out potential doublets and low-quality cells we calculated three quality measures for each individual cell:

1. The Good–Turing estimate of observed expression⁷⁵ given by S=1−n₁/N, where n₁ is the number of genes with one mapped read and N is the total number of reads in the cell. Saturation ranged from 40% to 100%. We chose to keep only cells with S>0.6.
2. The percentage of mitochondrial gene expression. We filtered out any cell with >5% of the total expression mapped to genes annotated to come from mitochondrial DNA.
3. The deviation from a read-UMI fitted curve: We expect the number of observed genes to increase linearly with the number of reads for cells that have not attained full saturation. We fitted a line between the number of nonzero genes and number of reads and filtered out cells whose residuals were >5 SDs from the line.

A total of 2750 cells were kept after filtering through these three criteria, indicating that 73% of sequenced cells were of high quality.

Analysis of Week 16 scRNAseq Dataset

We used the Seurat R package⁵⁶ version 1.4 for further analysis of the remaining cells. The MeanVarPlot function with default parameters was used to find a subset of genes whose variability is above the expected technical noise. We found 582 such genes, which were further used for principal component analysis (PCA function). To find significant PCs, we used the JackStraw test⁷⁶ and kept the first 24 PCs, which had P<10⁻⁴. These PCs were used for clustering using the FindClusters function with k=30 nearest neighbors. We found 12 clusters whose identities were further validated by the AssessNodes function, which builds a random forest classifier for each split node in the cluster hierarchy. The highest out-of-bag error we found was 9%, which indicates that all clusters have a clear identity. Differential expression was performed with the likelihood ratio proposed by McDavid et al.⁷⁷ and implemented in the FindAllMarkers function, in which genes inside a cluster are compared with the expression in all cells outside of the cluster. We set the minimum average difference between inside and outside clusters to 0.15 and no minimum average expression threshold.

Analysis of Cluster 4

We repeated the aforementioned procedures with only the subset of 318 cells that were assigned to cluster 4. We found 685 variable genes and the four first principal components to be statistically significant (P<10⁻⁴). Clustering with four principal components yielded four subclusters whose maximum out-of-bag error was 6.5%.

Disclosures

None.