Research ArticleHEALTH AND MEDICINE

Molecular characterization of the human kidney interstitium in health and disease

See allHide authors and affiliations

Science Advances  10 Feb 2021:
Vol. 7, no. 7, eabd3359
DOI: 10.1126/sciadv.abd3359

Abstract

The gene expression signature of the human kidney interstitium is incompletely understood. The cortical interstitium (excluding tubules, glomeruli, and vessels) in reference nephrectomies (N = 9) and diabetic kidney biopsy specimens (N = 6) was laser microdissected (LMD) and sequenced. Samples underwent RNA sequencing. Gene signatures were deconvolved using single nuclear RNA sequencing (snRNAseq) data derived from overlapping specimens. Interstitial LMD transcriptomics uncovered previously unidentified markers including KISS1, validated with in situ hybridization. LMD transcriptomics and snRNAseq revealed strong correlation of gene expression within corresponding kidney regions. Relevant enriched interstitial pathways included G-protein coupled receptor. binding and collagen biosynthesis. The diabetic interstitium was enriched for extracellular matrix organization and small-molecule catabolism. Cell type markers with unchanged expression (NOTCH3, EGFR, and HEG1) and those down-regulated in diabetic nephropathy (MYH11, LUM, and CCDC3) were identified. LMD transcriptomics complements snRNAseq; together, they facilitate mapping of interstitial marker genes to aid interpretation of pathophysiology in precision medicine studies.

INTRODUCTION

The structure and components of the renal interstitium underlie its functional and pathophysiologic role in the kidney. These components include a diverse assembly of cells such as fibroblasts, vascular smooth muscle cells (VSMCs), and immune cells (1, 2). Noncellular constituents include the extracellular matrix, proteins, and exosomes (3). Beyond providing a scaffold to support the glomeruli (GLOM) and tubules, components of the interstitium regulate processes such as the endocrine functions of erythropoietin and renin (46), extracellular and blood pressure homeostasis (7, 8), and adaptation to inflammation (911).

In contrast to the nephron and renal vasculature, the transcriptomic signature of the renal interstitium remains largely uncharacterized in health and disease. We sought to uncover the gene expression signature of the interstitium to better characterize the contribution of its resident cells to kidney health and pathology. We postulated that the aggregate interstitial transcriptomic signature would complement and enhance our understanding of individual renal cell type signatures acquired from single nuclear RNA sequencing (snDrop RNAseq) (12). In an overlapping set of nephrectomies and renal biopsy specimens, spatial context was provided for these single nuclear signatures by performing laser microdissection (LMD) of the interstitium without the nephron segments, but including capillaries, small vessels, and lymphatics. Last, we compared transcriptomic signatures of these reference nephrectomies to those of diabetic biopsy specimens to understand the effect of chronic injury on interstitial cell expression markers. The spatial information afforded by LMD allowed the analysis of interstitial cell expression in the context of its heterogeneous microenvironment, a valuable asset to assess expression changes in disease.

RESULTS

Expression of known markers

Samples of the interstitium, GLOM, proximal tubules (PTs), thick ascending limbs (TALs), and collecting ducts (CDs) were acquired by LMD from each of nine reference nephrectomy samples with no notable pathological changes. Known markers of each nephron subsegment have been described and serve as landmarks to identify cell types within the kidney (13). Total RNAseq revealed significantly increased expression of these expected markers in the corresponding subsegmental dissections (Fig. 1). A relative absence of marker expression was observed in neighboring subsegments, indicative of a precise sample collection with minimal contamination from adjacent structures (fig. S1).

Fig. 1 Known segmental nephron markers in the LMD interstitium compartment.

(A) Dissection of the renal interstitium revealed increased expression of collagen markers (COL1A1, COL6A2, and COL6A3) compared with the other subsegments. (B) The interstitial dissection lacked expression of glomerular markers such as nephrin (NPHS1), podocin (NPHS2), and podocalyxin (PODXL). (C) Proximal tubular expression of megalin (LRP2), cubilin (CUBN), and the sodium-hydrogen exchanger regulatory factor 3 (PDZK1) was higher as compared with the interstitium. (D) Uromodulin (UMOD), NKCC2 (SLC12A1), and calcium-sensing receptor (CASR) expression was higher in the thick ascending loop than the interstitium. (E) The interstitium had minimal expression of CD markers such as aquaporin 2 (AQP2), the amiloride-sensitive sodium channel subunit beta (SCNN1), and solute carrier family 4 member 9 (SLC4A9). Immunofluorescence staining with DAPI (nuclei), phalloidin (F-actin, brush border), Tamm-Horsfall protein (THP) (TAL), and peanut agglutinin (collecting duct) were used along with morphology to dissect each specific compartment. Dot plots are provided to visualize individual data points, but all P values correspond to the whole-transcriptome analysis provided in table S1 to maintain consistency throughout the study. ***P <0.05 after FDR correction.

Differential expression between subsegments and the interstitium

We sought to delineate the most differentially expressed genes (DEGs) in the interstitium and each nephron subsegment from the reference samples. To facilitate comparison, an expanded list of known markers is presented in Fig. 2A, illustrating their specificity for the subsegments with which they are associated. In contrast, Fig. 2B depicts the most significant genes that were differentially expressed between each nephron subsegment. The top 15 up-regulated genes for each subsegment are sorted by significance after removing noncoding genes and micro RNA species (Fig. 2C). There was an overlap between the known interstitial markers and the top DEG interstitial markers. For example, three of six known markers resided within the top 15 DEGs, and all six had a false discovery rate (FDR)–corrected level of significance of P < 0.05.

Fig. 2 Known markers and top DEGs in nephron compartments.

(A) Heatmap illustrating relative gene expression of expected subsegmental markers in the subsegments of interest from reference tissue. (B and C) Heatmap depicting the expression of the top 15 protein-coding differentially expressed markers for each subsegment sorted by significance. The identity of each gene is provided in (C). (D to G) Scatter plots indicating differential expression of selected interstitial markers (COL6A2, collagen 6A2; AEBP1, adipocyte enhancer-binding protein 1; SYNM, synemin; TPSB2, tryptase beta 2) identified in the unbiased analysis with representative images of immunohistochemical staining; each marker is obtained from the Human Protein Atlas. Dot plots are provided to visualize individual data points, but all P values correspond to the whole-transcriptome analysis provided in table S1 to maintain consistency throughout the study. ***P <0.05 after FDR correction.

At a significance threshold of P < 0.05, we found 4626 DEGs in the interstitial compartment as compared with other nephron segments. This included 136 genes at an FDR threshold of 0.05, and 607 genes at an FDR threshold of 0.1 (table S1). The top protein-coding DEG of the interstitium was COL6A2 (collagen 6A2; P = 2.4 × 10−7), known to be expressed in the interstitium and, to a lesser extent, the glomerulus (Fig. 2, A and D). Among the top DEGs, additional underrecognized markers of the renal interstitium were identified, beyond the predicted markers found in Fig. 2A. These additional up-regulated genes align with cell types known to reside in the interstitium. For example, expression of a stromal/structural protein adipocyte enhancer-binding protein 1 (AEBP1; P < 0.0001), the pericyte marker synemin (SYNM; P < 0.001), and an immune/mast cell marker tryptase beta 2 (TPSB2; P < 0.001) were all enriched in the renal interstitial compartment. These unbiased markers were recognized as interstitial markers by differential gene expression following LMD, and protein expression was localized to the interstitium with immunohistochemistry from the human protein atlas (Fig. 2, D to G). Thus, the differential expression analysis guided identification of additional specific interstitial markers and was supported by immunohistochemical staining in the kidney.

Interstitial pathway analysis

To uncover pathways enriched in the renal interstitium, we juxtaposed its expression signature to the mean expression of the remaining subsegments (GLOM, PT, TAL, and CD). A resulting transcriptogram illustrates the differential expression, significance, and pathway enrichment for the protein-coding genes (Fig. 3A). A comprehensive list of DEGs (table S1) and pathways (table S2) is provided. The enriched pathways included guanine nucleotide–binding protein (G protein)–coupled receptor (GPCR) ligand signaling (P = 0.0019) and collagen biosynthesis (P = 0.0029) (Fig. 3, B and C). Among the key genes driving enrichment of GPCR signaling, the taste receptor, TAS2R38 (P = 2.4 × 10−5), and KISS1 (kisspeptin; P = 8.0 × 10−5) were both found to have increased expression in the interstitial compartment as compared with the GLOM and tubules (table S1). The protein of TAS2R38 localizes to the interstitium upon immunohistochemistry (fig. S2A), although it appears in the glomerular endothelium as well. KISS1 also localized to the interstitium, and this is described in greater detail below. Collagen VI and collagen I, as well as the prolyl hydroxylases, like P3H1 and CRTAP, important in collagen assembly, contributed to the enrichment of collagen biosynthesis pathways in the interstitium. P3H1 and CRTAP were expressed in our interstitial dataset and contributed to the pathway enrichment, but neither gene reached an FDR-corrected level of significance. Expression was confirmed in the interstitium with immunohistochemistry and colocalized with fibroblasts by in situ hybridization, but the expression of these two genes was not specific and found in most tubules on immunohistochemistry and additional nonfibroblast cells by fluorescence in situ hybridization (FISH) (fig. S2). In contrast, collagen 4 was not found in the interstitium but was expressed in GLOM.

Fig. 3 Network and pathway analysis of gene expression enriched in the renal interstitium.

(A) A transcriptogram illustrating pathway enrichment of the interstitial compartment as compared with the mean expression of all tubular and glomerular compartments. Protein-coding genes are ordered along the x axis according to their network interactions. The top plot and middle plot illustrate the average expression and P value of the network region. The bottom plot provides pathway enrichment based on expression and significance. (B) S-curve plot of genes contributing to the GPCR ligand-binding pathway. (C) S-curve plot of genes contributing to the collagen biosynthesis and modifying enzyme pathway.

Deconvolution of the interstitial expression signature

To better understand the components of the renal interstitial expression signature, we compared the compartmental expression pattern of the interstitium to publicly available single nuclear RNA sequencing (snDrop seq) data (12). Single nuclear data were obtained for 15 reference samples, including 6 of the same reference nephrectomies that underwent LMD. Unsupervised clustering led to the identification of 30 distinct cell clusters (Fig. 4A), which could be reidentified as identical to the clusters in a prior analysis (12) using known marker gene expression (fig. S3A). Gene expression profiles of cell types known to reside in a specific compartment correlated most strongly with the transcriptomic signature of the corresponding microdissected compartment (Fig. 4B). For example, the podocyte, mesangial cell, and glomerular capillary expression signatures most strongly correlated with the dissected GLOM. In contrast to the dissected compartments, an RNAseq bulk expression analysis of undissected nephrectomy cross sections was less specific, revealing correlation with most of the snRNAseq clusters.

Fig. 4 Deconvolution and alignment of LMD expression with single nuclear RNA sequencing.

(A) Uniform Manifold Approximation and Projection (UMAP) of 30 snRNAseq clusters. (B) Cluster identities are provided for all 30 clusters, along with correlation plots between the snRNAseq clusters and compartmental transcriptomic signatures. Selected genes with increased expression in the interstitial compartment are displayed in the corresponding dot plot. In a separate analysis, bulk expression was used to calculate the expression correlation in the bulk heatmap column. (C) UMAP illustrating subclustering of the fibroblast population into three subpopulations. (D) Dot plot revealing the expression of marker genes associated with each of the three identified fibroblast subclusters. (E to H) Feature plots displaying selected genes associated in the fibroblast subclusters: ACTA2 as a marker of myofibroblasts, PDGFRA corresponding to fibroblast subcluster 3, and SPP1 and MYL12B corresponding to fibroblast subcluster 2.

The renal interstitial expression signature held the strongest correlation with the VSMC/pericyte, fibroblast, and nonglomerular endothelial cell clusters. The top DEGs of the interstitial compartment were similarly found to be expressed most strongly in these clusters. For example, ADCY5 and SYNM were more highly expressed in the VSMC and pericyte cluster 27, while spondin 1 (SPON1) and collagen VI (COL6A3) were most frequently expressed in fibroblast cluster 28. We also explored whether subclustering of cells in fibroblast cluster 28 at an increased resolution might yield additional clusters corresponding to important fibroblast subtypes, such as myofibroblasts (Fig. 4, C to H). This high-resolution subclustering led to three distinctive clusters: (i) a platelet-derived growth factor receptor A (PDGFRA)–positive cluster-positive cluster, (ii) an osteopontin/myosin light chain–positive cluster, and (iii) a cluster with a variety of retained epithelial cell markers. This third cluster may include cells undergoing epithelial-to-mesenchymal transition. Actin alpha 2 (ACTA2)–positive fibroblasts/myofibroblasts were distributed among all three fibroblast subclusters, and a distinct myofibroblast cluster could not be observed in this dataset even at higher resolutions. A subset of cells coexpressed both ACTA2 and PDGFRA in the snRNAseq dataset. Receptor ligand interactions between fibroblasts, proximal tubular cells, and other interstitial cell types were assessed, and significant differences were not identified in cross-talk with other cell types for each of the three fibroblast subclusters (fig. S3, B and C).

Feature plots in fig. S3 (D to L) reveal the distribution of cells expressing other top genes from LMD regional transcriptomics mapped onto the snDrop RNAseq Uniform Manifold Approximation and Projection (UMAP). These genes corresponded to endothelial cells, immune cells, VSMCs, and fibroblasts. While there exists considerable overlap between technologies, unique DEGs were identified by both platforms.

To further understand the origin of DEGs unique to the LMD renal interstitial signature, we carried out pathway-based deconvolution using the snDrop RNAseq dataset. Using a dimensionality reduction approach (see Materials and Methods), the 30 snDrop RNAseq clusters were regrouped (or reclustered) into major cell type categories (e.g., all nonglomerular endothelial cell clusters were combined into a single endothelial cell type), and the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways were attributed to each cell type based on enrichment of DEGs from the snDrop RNAseq dataset. The member genes of these pathways were then used as filters to screen the unique LMD DEGs as described in Materials and Methods. In this analysis, three conditions had to be met: (i) a unique DEG in the LMD interstitial compartment; (ii) a member of a pathway enriched in the fibroblast, VSMC, endothelial, or immune cell types; and (iii) localized to the interstitium on immunohistochemistry (Table 1, part A, and fig. S2, A to I). Using this approach, 13 genes were identified.

Table 1 Pathway-based deconvolution of interstitial genes and expression changes in disease.

REF, reference samples; DN, diabetic nephropathy; FIB, fibroblast; END, endothelial cell; VSMC, vascular smooth muscle cell or pericyte; IMM, immune cell.

View this table:

Among the unique genes in the aggregate interstitial expression signature, we identified several extracellular matrix proteins (COL6A2, COL1A1, etc.) in the analysis. Expression of these markers was not present in the original snDrop RNAseq dataset. In addition, known markers of immune cells were identified in the aggregate signature, which were likely too lowly represented in the snDrop RNAseq dataset to form unique clusters (GZMB, NCF2, etc.). Last, the deconvolution approach also identified interstitial expression of the GPCRs TAS2R38 and KISS1, discussed earlier. Pathway enrichment for these genes included relevant pathways for endothelial cells, VSMCs, and fibroblasts.

Localization of unique differentially expressed interstitial genes

To improve the confidence with which novel interstitial genes were localized, we conducted a series of FISH- and antibody-based immunofluorescence studies. Localization of COL6A1, CELA3A, and KISS1 (Fig. 5) was tested with known markers for endothelial (CD31), VSMC/pericyte (TAGLN), fibroblast (PDGFRA), and myofibroblast (ACTA2) cells. Immunohistochemistry images acquired from the human protein atlas localized all three markers to the renal interstitium. By in situ hybridization, COL6A1 and CELA3A specifically colocalized with both PDGFRA and ACTA2. Similar to the observations in the snRNAseq dataset (Fig. 4), a subset of fibroblasts/myofibroblasts coexpressed PDGFRA and ACTA2 by in situ hybridization (Fig. 5). Potential explanations for this coexpression will be delineated in the discussion. KISS1 specifically colocalized with the VSMC marker, TAGLN. As expected, none of the three genes (COL6A1, CELA3A, and KISS1) were localized in cells with tubular epithelial cell morphology. The fibroblast, VSMC, and endothelial markers did not colocalize with epithelial cells. FISH for TAS2R38 was undertaken, but only a small number of cells expressed TAS2R38 in the kidney. As a result, conclusions on colocalization for this gene could not be drawn beyond its immunohistochemistry localization to the interstitium.

Fig. 5 Localization of identified interstitial cell genes with in situ hybridization and immunofluorescence.

Colocalization of COL6A1 (collagen 6A1), KISS1 (kisspeptin 1), and CELA3A (elastase 3A) with known markers of endothelial cells, VSMCs, and fibroblasts. (A, F, and K) Immunohistochemical staining of COL6A1, KISS1, and CELA3A, as illustrated in the Human Protein Atlas. (B to E, G to J, and K to O) Single-molecule FISH colocalization. The 5-μm measurement bar of “E” corresponds to all FISH images unless otherwise noted. COL6A1 (B; purple) colocalized with fibroblast markers PDGFRA (C; red) and ACTA2 (D; green) and in a merged image of all channels (E; *colocalization of all three markers]; KISS1 (G; purple) colocalizing with VSMC marker TAGLN (H; red), and minimal expression of the endothelial cell marker CD31 (I; green) and with merged image in (J) (#colocalization of purple and red); CELA3A (L; purple) colocalized with PDGFRA (M; red) and ACTA2 (N; green); merged image (O; white arrows indicate colocalization of markers, and the yellow arrowhead indicates an RNA clump artifact, which should not be mistaken for colocalization). (P to S) DAPI staining of tubular nuclei is seen in a tubular cross section (P and R) with selected nuclei at higher magnification (Q and S) illustrating the absence of targeted RNA molecules associated with the tubular cell. (T to W) Immunofluorescence staining illustrating the presence of CELA3A protein within interstitial cells and lack of colocalization with endothelial cells, further supporting CELA3A expression being associated with fibroblasts. (T and U) CELA3A staining in the interstitial compartment. BV, blood vessel; DT, distal tubule; G, glomerulis; PT, proximal tubule. Phalloidin was used to visualize F-actin, and DAPI was used to stain nuclei. (V and W) CELA3A (red) not colocalized with CD31 stain (in blue). Arrowheads indicate CELA3A-positive cells.

Antibody-based immunofluorescence for CELA3A localized protein expression to the interstitium, neither colocalizing with epithelial or endothelial cells. This pattern of expression is consistent with in situ hybridization data revealing expression of this gene in the fibroblast.

Differential expression in diabetic kidneys

A potential limitation of single-cell and single nuclear RNA sequencing technologies is the reliance upon marker genes for cell identification. Understanding how expression of marker genes changes in disease is paramount to cell type identification. We hypothesized that transcriptomics obtained from LMD would allow a tissue-guided, spatially registered identification of markers in the interstitium that either change or remain stable in disease.

The histopathologic phenotype of diabetic nephropathy often includes significant interstitial fibrosis (14). We assessed the differential expression of genes in the interstitium of renal biopsies of patients with diabetic nephropathy, compared with the reference nephrectomies. The age and glomerular obsolescence of these groups were similar; however, the area affected by interstitial fibrosis and tubular atrophy (IFTA) in the diabetic samples was greater (45.7% versus 20.0%) by design. Other histopathologic and clinical characteristics are included in table S3. Our goal was not to uncover the pathogenesis of early diabetic nephropathy from samples matched by IFTA; instead, we strived to understand the relevance of transcriptomic markers in a common human type of chronic kidney disease (CKD) and to understand the effect of interstitial fibrosis on the interstitial signature. Additional insights gathered regarding diabetic nephropathy were welcomed, but not the central goal.

A number of genes were differentially expressed between the reference and diabetic interstitium (table S4). Marker genes of interstitial cell types were ranked on the basis of significance and categorized as up-regulated, down-regulated, or unchanged in diabetic specimens with selected examples given in Fig. 6. Examples of markers that are up-regulated in diabetic nephropathy specimens include plexin A2 (PLXNA2) expressed in endothelial cells and the Cbl proto-oncogene (CBL) expressed in macrophages (table S5). Of note, CBL is known to suppress macrophage activation in insulin resistance (15). The expression of certain markers remained unchanged in the setting of diabetic nephropathy, including NOTCH3 and JAG1 in VSMCs, HEG1 in endothelial cells, and ANGPT2 in fibroblasts (table S5). While not a marker in snDrop RNAseq, expression of TPSB2 to identify mast cells was found in both health and disease. Genes with down-regulated expression in diabetic nephropathy included MYH11 in VSMCs, LUM in fibroblasts, and CCDC3 in endothelial cells (Table 1, part B) (16).

Fig. 6 Network and pathway analysis of DEGs in healthy and diseased renal interstitium.

(A) S-curve plot of DEGs between diabetic nephropathy and reference samples contributing to extracellular matrix organization. (B) S-curve plot of genes contributing to small-molecule catabolic processes. (C to H) Scatter plots of example genes with differential expression in reference and diabetic samples. Examples include up-regulation of plasminogen (PLG) and fatty acid binding protein 1 (FABP1) and down-regulation of collagen 11A2 (COL11A2), dimethylarginine dimethylaminohydrolase 2 (DDAH2), pyruvate kinase L/R (PKLR), and lumican (LUM). (*P < 0.05 uncorrected, **P < 0.1 FDR-corrected, and ***P < 0.05 FDR-corrected). (I) Transcriptogram illustrating pathway enrichment between the diabetic nephropathy and reference interstitial compartments.

Pathway analysis of the interstitium in reference and diabetic kidneys

The gene expression signatures of the reference and diabetic interstitial samples were compared to uncover enriched pathways, as defined by KEGG, Reactome, and Gene Ontology (GO) terms (table S6). Most genes in diabetic samples had reduced expression as compared with the reference samples (Fig. 6, C to F). Among the most distinctly affected pathways were small-molecule catabolism (P = 0.0008) and extracellular matrix organization (P = 0.002). The most significant DEG in the small-molecule catabolism pathway was DDAH2 (P = 9.5 × 10−8), a gene that encodes dimethylarginine dimethylaminohydrolase 2. DDAH2 metabolizes asymmetric dimethylarginine, a metabolic by-product that contributes to uremic vasculopathy (17). Contributing to extracellular matrix organization enrichment in diabetic nephropathy, COL11A2 was down-regulated (1.6 × 10−9), but ITGB4 was up-regulated (P = 0.001). Plasminogen (PLG) was also up-regulated (P = 1 × 10−4). Angiostatin, a proteolytic product of PLG, inhibits angiogenesis and has been shown to be up-regulated in the serum of diabetics (18).

DISCUSSION

In this study, we characterized the expression signature of the renal interstitium obtained from LMD. This aggregate signature was deconvolved using orthogonal snDrop RNAseq signatures of fibroblasts, VSMCs, nonglomerular endothelial cells, and immune cells. Markers were spatially localized to the renal interstitium during LMD, both validating existing snDrop RNAseq markers and uncovering a small set of genes, which were not detected in the snDrop RNAseq platform. A subset of markers was colocalized to specific cell types with FISH. Furthermore, markers were subcategorized as those that maintained expression in disease or had differential expression in the interstitium of patients with diabetic nephropathy. Given the significance of interstitial fibrosis in the progression of kidney disease, elucidating the composition of the interstitium and its role in health and disease is of paramount importance. The interstitium is not simply a passive compartment housing fibroblasts, pericytes, and other cell types, but a space where hormonal processes, extracellular homeostasis, and adaptation to inflammation take place. Understanding its role in progression of fibrosis, resistance to injury, and regeneration of neighboring cells on a patient-by-patient basis will enable integration with precision medicine and an individualized therapeutic approach to treat disease. Therefore, with its depth, the data contained in this study will become an important reference useful to further expand our understanding of the role of the renal interstitium in maintaining homeostasis and modulating response to disease.

The use of tissue markers is a well-established practice to identify cells or regions of interest. Numerous markers specific to cell types of the nephron have been identified (13, 19). For example, NPHS1 (nephrin) and NPHS2 (podocin) are known to be predominantly expressed in GLOM, while AQP2 (aquaporin 2) is characteristic of the CD. In the interstitium, extracellular matrix proteins, such as collagens and prolyl hydroxylases, have been used as markers. In humans, many of these markers have been identified on the basis of assessment of tissue without known pathology. As we have shown, marker expression can change in disease; therefore, this study lends confidence to the use of certain markers in both healthy and pathological tissues. For example, the expression of TPSB2, a mast cell marker (20), did not change significantly in the setting of diabetic nephropathy.

Single-cell and single nuclear sequencing of the kidney has enhanced our understanding of the transcriptomic signature of each renal cell type (21, 22). One critique of single-cell and single nuclear RNA sequencing is that the cells lack spatial information. However, we found that the markers used to cluster nuclei of renal cells mapped reliably. Regional transcriptomics increased confidence in the spatial localization of these markers to the interstitium. In a complementary way, the use of snDrop RNAseq data increased the confidence of mapping marker genes uniquely expressed in the LMD dataset to a cell type of origin. A small number of additional markers were identified through regional transcriptomics or, alternatively, found to be down-regulated in disease.

RNA in situ hybridization allowed precise localization of the novel markers we identified through LMD to the interstitium. For example, KISS1 rather specifically colocalized with a VSMC marker, TAGLN. COL6A1 and CELA3A both colocalized with fibroblast markers. We did observe coexpression of PDGFRA and ACTA2 on in situ hybridization for a subset of fibroblasts. We speculate that this finding is related to the sensitivity of the RNAscope technology, wherein gene expression was observed for both genes, but cells may not have expressed both proteins. A second alternate hypothesis is that these cells were undergoing a transition between a fibroblast and a myofibroblast state, as has been reported by the Benjamin Humphreys laboratory and others (2325)). This hypothesis is supported by the snRNAseq data, which also identified coexpression of PDGFRA and ACTA2 in a subset of cells within the fibroblast cluster. Regardless, the multitude of orthogonal datasets, including pathway-based deconvolution, immunohistochemistry, antibody-based immunofluorescence, and in situ hybridization, all support the fibroblast/myofibroblast as the origin cell(s) for these genes.

A limitation of this study is its relatively small sample size, which precludes incorporation of clinical data into these transcriptomic models. Considerable heterogeneity exists among humans, and dozens of specimens are required to glean insights into the relationship between clinical features and the underlying biology of disease. In addition, nephrectomy samples used in this study, while valuable, are never completely free of pathology as some are acquired from deceased donors or the unaffected regions of tumor nephrectomy tissue. However, these particular samples were reviewed by pathologists and found to have minimal glomerular obsolescence, interstitial fibrosis, or tubular atrophy. The reference nephrectomies did not have histologic evidence of diabetic or hypertensive nephropathy. Unfortunately, a full clinical dataset of all comorbidities and medication usage was unavailable for the reference nephrectomies, as it was in the diabetic samples. Last, a limitation inherent in all LMD is that it represents an enriched compartment of cells and tissue, not the specific cells themselves. Neighboring cells may be inadvertently captured along with regions of interest. However, the correlation analyses of our LMD compartments with snDrop RNAseq suggest relatively specific collections. Given the limitations noted, we are cautious to avoid drawing sweeping conclusions about the pathogenesis of diabetic nephropathy from this study. Instead, the value of this study lies in the localization of transcriptomic signals as well as discerning the stability of cell type–specific markers in reference and diseased tissue.

LMD facilitates interrogation of whole tissue sections adherent to a slide, allowing for excision of regions of interest based on the morphology, staining, or proximity to other structures. Consequently, the cell type, location, and volume of dissected area within a specific region of the tissue provide important insight into the functionality of the organ. Beyond contrasting reference and disease samples in this study, the spatial information afforded by LMD provides the opportunity to compare regions within the same kidney. In this study, we found significant enrichment of immune signals in the interstitium as compared with the GLOM. This is explained by the fact that our reference samples did not have glomerulonephritis but may have had interstitial immune cell infiltration from hypoxia during harvest. Therefore, regional LMD-based transcriptomics can determine the site of immune signaling activity, which will be relevant in investigating the pathobiology of kidney disease. Future directions will include comparisons of the cortical and medullary interstitium, as well as areas with and without significant inflammation or fibrosis.

Transcriptomic analysis of the kidney has progressed considerably in the past decade. Advances in RNAseq technology have enabled acquisition of whole-genome expression signatures in archived kidney tissue samples (26). In this study, we showed the feasibility and advantages of obtaining the transcriptomic signature of the interstitium in reference nephrectomy and diabetic biopsy specimens. These insights are crucial to deepen our understanding of the interstitial space in maintaining homeostasis and modulating response to disease. The LMD-based transcriptomic approach identified cell type marker genes for which expression changes during kidney disease. Therefore, by using tissue landmarks to guide the dissection, LMD transcriptomics can be very useful in providing a gold standard for gene expression changes in the interstitium. In addition, using snDrop RNAseq data to assist in deconvolution of the aggregate interstitium signature, we show that LMD regional transcriptomics and snDrop RNAseq are complementary methods that allow for cross-examination of the transcriptome.

MATERIALS AND METHODS

Tissue source

Nine reference kidney samples (from tumor-free regions of nephrectomies or deceased donors) and six diabetic nephropathy biopsy specimens were obtained from the Kidney Precision Medicine Project (KPMP) or Biopsy Biobank Cohort of Indiana. This study was approved by the Institutional Review Board (IRB no. 1906572234) of Indiana University.

LMD and RNA isolation

Twelve-micrometer frozen sections were cut from tissue blocks preserved in optimal cutting temperature (OCT) medium. The sections were adhered to PPS [poly(p-phenylene sulfide)] membrane slides (Leica, Buffalo Grove, IL). Slides were processed using a “Rapid Stain” protocol (27, 28) involving (i) fixation with ice-cold acetone for 30 s twice (Sigma-Aldrich, cat. no. 270725-1L), (ii) two 30-s washes with ribonuclease (RNAse)–free phosphate-buffered saline (PBS) (VWR, cat. no. K812-500ML), (iii) antibody application in 10% bovine serum albumin (BSA) (VWR, cat. no. 0332-100G) for 5 min, (iv) two 30-s rinses with 10% BSA, and (v) quick slide air drying. Up to three of the following antibodies facilitated compartment isolation: (i) Tamm-Horsfall protein (THP), directly conjugated to Alexa Fluor 555, 1:50 concentration, to identify TALs (29); (ii) peanut agglutinin lectin (Vector Laboratories, cat. no. FL1071), 1:50 concentration, to identify CDs (30); (iii) phalloidin (Oregon Green 488, Thermo Fisher Scientific, cat. no. O7466), 1:25 concentration, to visualize F-actin and brush border on PTs, as well as GLOM; (iv) 4′,6-diamidino-2-phenylindole (DAPI) (Thermo Fisher Scientific, cat. no. 62248), to visualize nuclei, and alternatively, instead of THP, we use a megalin/LRP2 antibody (Abcam, cat. no. ab76969), 1:50 concentration, followed by Alexa Fluor 555 secondary antibody for 5 min, to identify PTs (31).

Slides were dissected on a Leica LMD6500 system with pulsed ultraviolet (UV) laser and ±0.5-μm accuracy for a 20× objective. Each subsegment was carefully dissected following the perimeter of a tubule or glomerulus. The interstitium was isolated following the outline of the area located between the tubules and along the tubular basement membrane, all contained within cortical space. For the purpose of this work, we extend the definition of the interstitium to include both small vascular and lymphatic vessels in addition to fibroblasts, dendritic cells, immune cells, and particles of extracellular matrix. Dissected tissue was collected in an autoclaved RNAse-free microcentrifuge tube containing RNA extraction buffer from the Arcturus PicoPure RNA Isolation Kit (Applied Biosystems, cat. no. KIT0204). A minimum area of 500,000 μm2 was collected for each kidney subsegment: interstitium, GLOM, PT (S1 and S2), TAL, and CD. A detailed LMD protocol can be found on protocols.io (www.protocols.io/view/laser-microdissection-8rkhv4w). RNA was isolated using the PicoPure RNA Isolation Kit according to the manufacturer’s instructions.

Whole-transcriptome sequencing of ultralow input total RNA

An Agilent Bioanalyzer assessed RNA quality. A DV200 (percent of transcripts of greater than 200 nucleotides in length) of 25% or more was the minimum threshold. At least 0.5 ng of total RNA was used for each sample’s library preparation. Ribosomal RNA was removed using the standard RiboGone–Mammalian Kit protocol (Clontech, cat. no. 634847) to enable expression measurements of fragmented RNA from laser microdissected specimens (26, 32, 33). Double-stranded complementary DNA was synthesized and amplified using the SMARTer Universal Low Input RNA Kit protocol (Clontech, cat. no. 634938) (26). Each resulting barcoded library was quantified with quality assessed by Bioanalyzer. Multiple libraries were pooled in equal molarity. Eight microliters of 100 pM pooled libraries was loaded into lanes of an Illumina HiSeq 4000. Approximately 100 million reads were generated per library. Sequencing quality control metrics are provided in table S3.

Mapping and expression

Mapping of reads to the hg38 reference genome was performed using RNA-seq aligner STAR (v2.5.2b) (34). After preprocessing and deduplication (35), read counts were quantified with featureCounts (subread v.1.5.0) (36) applying a Q10 threshold for mapping quality. Intronic reads were not counted toward expression. Total read counts mapping to each gene were generated with edgeR (37), normalized, and converted to expression ratios. Reference and diabetic datasets were quantile normalized jointly. The value of each normalized read count was increased by a value of 2, and a ratio was calculated on the basis of the following formula: X = (Expregion1)/mean [(Expregion2 + Expregion3 + Expregion4 + Expregion5)]. “Exp” was the quantile-normalized read count of each gene, region 1 refers to the region of interest (e.g., interstitium), and regions 2 to 5 refer to the other LMD subsegments (Glom, PT, TAL, or CD). All values were log base 2 transformed. Bulk expression was calculated in an analogous fashion with log base 2 transformation in a separate but parallel analysis. LMD compartment expression was not divided by bulk expression because each LMD compartment was a component of the bulk.

To quantify the proportion of exonic coverage for each gene, intron-to-exon (I/E) coverage ratios were calculated according to the following equation: I/E = [(number of intronic counts) × (length of counts in nucleotides)/(sum of the length of all intronic regions in nucleotides)]/[(number of exonic counts) × (length of counts in nucleotides)/(sum of the length of all exonic regions in nucleotides)]. Although intronic reads are a feature of ribosomal depletion RNA sequencing, they are not counted toward composite transcript expression. As a metric of sequencing fidelity, the ratio of I/E read coverage was assessed, and more than 13,000 genes had more than 70% exonic or junctional coverage in each reference sample, corresponding to an I/E ratio less than 3/7 or 0.428. The I/E ratio of housekeeping genes, known interstitial markers, and key focus genes discussed in this manuscript are exhibited in fig. S4.

Transcriptogram and pathway analysis

Genome-wide expression of data was visualized by transcriptogram methodology as previously described (38). In brief, protein-coding genes are ordered according to their protein-protein interaction network, defined from a STRING database with a confidence score >0.800. A Monte Carlo algorithm that minimizes a cost function draws groups of connected genes together. These groups of connected genes tend to share biological functions, clustering them in closer proximity along the x axis. Significance is calculated by a two-tailed Welch’s t test after 500 random permutations are conducted to estimate the FDR. Enrichment profiles are generated for each GO term, KEGG pathway, or Reactome pathway using clusterProfiler or ReactomePA (39). The top panel of the figure includes a gray horizontal plot serving as the comparator baseline with the indigo plot as the direction of differential expression. The middle panel represents significance of differential expression. The bottom panel illustrates pathway enrichment by the dimension and height of the respective peaks.

Single nuclear analysis

Six KPMP reference nephrectomies underwent both LMD regional RNA sequencing and nuclear isolation for snDrop RNAseq (12). Publically available data were acquired [Gene Expression Omnibus (GEO), GSE121862] for these samples and nine additional samples from the Washington University Kidney Translational Research Center. Nuclei data for all 15 samples were filtered for quality control before their deposition in the GEO as previously described (12). Nuclei were reclustered using PAGODA2 (https://github.com/hms-dbmi/pagoda2) and displayed as a UMAP using Seurat 3.

Following the methodology of Lake et al. (12), the 2000 genes with highest dispersion were included in the principal components analysis. An infomap community detection algorithm assisted in reclustering after considering the 30 closest neighbors measured through a cosine distance in an approximate k-nearest neighbor graph (top 150 principal components). Thirty clusters were obtained, which were identical to those previously reported. Clusters were annotated based on cell type marker genes. Visualization of expression in dot plots and feature plots was performed in Seurat 3.

Expression correlation

For correlation analyses between LMD compartment expression data and snDrop RNAseq in Fig. 4B, averaged gene expression values were obtained from the GEO count matrix for cluster marker genes. Pearson correlation coefficients between the expression of those genes in LMD and snDrop RNAseq were visualized as a heatmap using the ggplot2 package in R. Bulk expression was compared with snDrop RNAseq in a separate, but parallel, analysis. For fig. S1E, multivariate normal distribution was modeled for the gene expression of the top eight marker genes of the respective segment. The probability density function, including evaluating the likelihood of each microdissected sample being classified as a particular nephron compartment, is presented in a heatmap.

Pathway-based deconvolution

Differential expression analysis of single nuclei clusters (N = 30) was performed using a Wilcoxon rank sum test on genes detected in over 10% of nuclei in a cluster. The gene markers for each cluster were selected using Seurat 3 with P value and positive log fold change thresholds of 0.01 and 0.25, respectively. Marker lists were reduced and merged into four major cell type groups (endothelial, vascular, fibroblast, and immune) for the clusters most highly correlated with LMD interstitial expression by Pearson correlation. The combined marker set for each cell type group was tested for KEGG pathway enrichment using clusterProfiler (at P < 0.05). Significant pathways were queried to select genes exclusive to each cell type group and cross-referenced with the LMD data. To improve the confidence of linking a marker gene to a cell type, each gene was (i) uniquely represented among member genes in the enriched KEGG pathways for a given cell type and (ii) differentially expressed (up-regulated) in the LMD interstitium at P < 0.05. Marker genes commonly identified in both snRNAseq and LMD regional transcriptomic sequencing, as well as genes unique to LMD, were annotated.

FISH assay

Fresh-frozen OCT-embedded human kidney tissue cross sections were obtained at a thickness of 10 μm. On the day of the experiment, the tissue sections were fixed with 4% paraformaldehyde (PFA) for 15 min at 4°C, and then the slides were dehydrated with a series of ethanol concentrations (50, 70, and 100%) each for 5 min at room temperature (RT). Following that the tissue sections were treated with hydrogen peroxide for 10 min at RT, they were then washed five times with 1X PBS. All tissue sections were pretreated with RNAscope protease IV (Advance Cell Diagnosis Inc., USA) for 30 min at RT; after the treatment period, the tissues were washed five times with 1X PBS and immediately proceeded to the in situ hybridization. RNA in situ hybridization was performed using RNAscope multiplex Fluorescent Reagent Kit v2 (Advance Cell Diagnosis Inc., USA) as per protocol. Probe sets that target the following human RNA targets Hs-TAS2R38 (cat. no. 405571), Hs-CELA3A-01 (cat. no. 886731), Hs-CRTAP (cat. no. 886721), Hs-COL6A1 (cat. no. 82461), Hs-KISS1 (cat. no. 507981), Hs-P3H1 (cat. no. 886741), Hs-ACTA2 (cat. no. 444771-C2), Hs-PECAM1-01 (cat. no. 487381-C2), Hs-PDGFRA (cat. no. 604481-C3), and Hs-TAGLN (cat. no. 498961-C3) were purchased from Advance Cell Diagnosis Inc. TSA Fluorescein Plus, Cyanine 3 Plus, and Cyanine 5 Plus Evaluation Kit (PerkinElmer Inc., MA, USA) were used as secondary probes for the detection of RNA signals. All slides were counterstained with DAPI, and coverslips were mounted using fluorescent mounting media (ProLong Gold Antifade Reagent, Life Technologies). The images were obtained using an LSM 800 confocal microscope (Carl Zeiss, Germany).

Immunofluorescence imaging

Reference tissue sections (50 μm) were obtained via cryosectioning and fixed in 4% PFA for 24 hours, followed by a PBS wash and blocking in 10% normal donkey serum for 2 hours, and then application of the primary antibodies: CELA3A (Sigma-Aldrich, cat. no. HPA028086) and CD31 (NovusBio, cat. no. NB600-562) for an overnight incubation. The tissues were then washed with PBST (PBS–Tween 20) throughout the day, and the following probes were applied: donkey anti-rabbit Alexa Fluor 555, donkey anti-mouse Alexa Fluor 633, Oregon Green phalloidin 488, and DAPI. After an overnight incubation, all sections were washed throughout the day, then fixed in 4% PFA for 15 min, washed again in PBS, and mounted on slides in Glass antifade mountant (Thermo Fisher Scientific, cat. no. P36980). All staining steps have been performed at RT.

A Leica SP8 confocal microscope was used to acquire the image stacks in four channels using 20× numerical aperture, 0.75× objective, and spacing of 1 μm. Image processing was performed using ImageJ.

Statistics

When not otherwise specified above, data were given as means ± SE. Groups were compared by two-tailed t test or analysis of variance (ANOVA) (for more than two groups). P values <0.05 were considered significant. An FDR was calculated after applying a Benjamini-Hochberg correction for all transcripts with expression measured in five or more reference samples (N = 24,387). All P values presented in the text are uncorrected unless otherwise stated; however, FDR-corrected levels of significance are provided in table S1 and denoted in each figure as required.

SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at http://advances.sciencemag.org/cgi/content/full/7/7/eabd3359/DC1

https://creativecommons.org/licenses/by-nc/4.0/

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.

REFERENCES AND NOTES

Acknowledgments: We would like to thank the investigators of the Kidney Precision Medicine Project (www.kpmp.org) for gracious support and advice. Funding: Support for this work was provided by the NIH/NIDDK K08DK107864 (M.T.E.), NIH/NIDDK UH3DK114923 (T.M.E.-A. and P.C.D.), R01DK099345 (T.A.S.), and an NIH-NIDDK-DK076169 Diacomp award (P.C.D. and T.M.E.-A.). Research reported in this manuscript was supported by the National Institute of Diabetes and Digestive and the Kidney Diseases (NIDDK) KPMP (www.kpmp.org), under award no. U2CDK114886. Author contributions: Conception, analysis, and interpretation of data: M.T.E., D.B., R.M.F., F.S., Y.-H.C., K.S.C., S.W., M.J.F., P.C.D., and T.M.E.-A.; drafting/revising the article: M.T.E., D.B., T.M.E.-A., and P.C.D.; providing intellectual content: C.L.P., K.J.K., K.W.D., B.H.R., T.H., T.A.S., and S.V.P.; final approval of the version to be published: all authors. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Data are archived in the GEO (GEO number GSE163603). Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article