Biological pathways and comparison with biopsy signals and cellular origin of peripheral blood transcriptomic profiles during kidney allograft pathology.

Kidney transplant injury processes are associated with molecular changes in kidney tissue, primarily related to immune cell activation and infiltration. How these processes are reflected in the circulating immune cells, whose activation is targeted by strong immunosuppressants, is poorly understood. To study this, we analyzed the molecular alterations in 384 peripheral blood samples from four European transplant centers, taken at the time of a kidney allograft biopsy, selected for their phenotype, using RNA-sequencing. In peripheral blood, differentially expressed genes in 136 rejection and 248 no rejection samples demonstrated upregulation of glucocorticoid receptor and nucleotide oligomerization domain-like receptor signaling pathways. Pathways enriched in antibody-mediated rejection (ABMR) were strongly immune-specific, whereas pathways enriched in T cell-mediated rejection were less immune related. In polyomavirus infection, upregulation of mitochondrial dysfunction and interferon signaling pathways was seen. Next, we integrated the blood results with transcriptomics of 224 kidney allograft biopsies which showed consistently upregulated genes per phenotype in both blood and biopsy. In single-cell RNASeq (scRNASeq) analysis of seven kidney allograft biopsies, the consistently overexpressed genes in ABMR were mostly expressed by infiltrating leukocytes in the allograft. Similarly, in peripheral blood scRNASeq analysis, these genes were overexpressed in ABMR in immune cell subtypes. Furthermore, overexpression of these genes in ABMR was confirmed in independent cohorts in blood and biopsy. Thus, our results highlight the immune activation pathways in peripheral blood leukocytes at the time of kidney allograft pathology, despite the use of current strong immunosuppressants, and provide a framework for future therapeutic interventions.

Despite excellent short-term outcomes of kidney transplantation, allograft survival in the long term is still hampered by transplant-specific diseases, most importantly allograft rejection and polyomavirus-associated nephropathy.1, 2, 3, 4 In kidney allograft rejection, graft-infiltrating cells are activated primarily in lymphoid organs before traveling and infiltrating the graft, despite the use of strong immunosuppressants. Apart from high-dose steroids and T cell–depleting agents such as thymoglobulin, there is no sufficiently efficacious therapy for rejection, partly because the primary pathophysiological processes behind this immune escape are poorly understood. In polyomavirus-associated nephropathy, acquired immunodeficiency permits reactivation of latent polyomavirus in the urinary tract epithelium and infection of kidney tubular epithelial cells, with direct viral toxic effects on the cells and activation of an antiviral host response, leading to tubulo-interstitial inflammatory infiltrates in the graft.

In the past decade, transcriptional studies in kidney transplant biopsies have improved our pathophysiological insights into rejection. As immune cell infiltrates in the graft tissue originate from immune activation outside the kidneys (i.e., activation in lymphoid tissue), one hypothesis is that molecular changes seen in kidney allografts at the time of rejection or virus-associated nephropathy are also present in circulating immune cells. We recently illustrated the potential of identifying biologically relevant genes and pathways in peripheral blood at the time of antibody-mediated rejection.

In this study, we used RNA-sequencing for whole transcriptome screening in a set of peripheral blood samples selected for rejection and polyomavirus nephropathy, and we compare this to the transcriptomics of heterogeneous control samples without these diseases. We hypothesized that the transcriptional changes in peripheral blood can unveil the key pathways activated in circulating peripheral blood immune cells in kidney transplant recipients, despite the use of powerful immunosuppressants, facilitate their noninvasive diagnosis, and provide potential targets for new therapies.

Methods

Study design, patient population, and sample collection

Peripheral blood samples were selected from an existing biobank prospectively collected as part of the BIOMARGIN (BIOMarkers of Renal Graft Injuries, www.biomargin.eu, ClinicalTrials.gov, NCT02832661) and ROCKET (Reclassification using OmiCs integration in KidnEy Transplantation) multicenter studies performed in 4 European transplant centers (Hôpital Necker, Paris, France; University Hospitals Leuven, Leuven, Belgium; Medizinische Hochschule, Hannover, Germany, and Centre Hospitalier Universitaire Limoges, Limoges, France). Samples were prospectively and consecutively collected at the time of kidney allograft biopsies, between June 2011 and August 2016. In the 4 clinical centers, protocol kidney allograft biopsies were performed at 3, 12, and sometimes 24 months after transplantation, according to local center practice, in addition to clinically indicated biopsies (biopsies at the time of graft dysfunction). All biopsies were used for histologic diagnosis according to the Banff 2019 classification, whereas in some of the concomitant biopsy samples (n = 95), additional molecular analysis was performed. Single-end, 50–base pair (bp) sequencing was performed after RNA extraction and quality control of the peripheral blood samples, at Genomics Core facility, KU Leuven, Leuven, Belgium, using an Illumina HiSeq 4000 device, after library preparation using the Lexogen QuantSeq mRNA kit according to manufacturer instructions. Microarray gene expression data of 224 renal allograft biopsy samples (including 95 matched samples) from the same 4 European centers were available for data integration. Details on sample collection, RNA extraction, and the microarray platform used were published previously (Gene Expression Omnibus database under the accession number GSE147089). Previously published human single-cell data from 7 kidney allograft biopsies were used to determine the corresponding cell types expressing the antibody-mediated rejection (ABMR) genes observed in both blood and biopsy. The associated raw counts or matrices were downloaded from the Gene Expression Omnibus (GEO, GSE145927, https://www.ncbi.nlm.nih.gov/geo) and the Kidney Precision Medicine Project (https://atlas.kpmp.org/repository). Cell enrichment analysis was performed using the xCell online tool, which is a computational method used to investigate changes in cell distribution based on bulk transcriptomic profiles. We validated the observed ABMR signals obtained in this study in 2 independent publicly available datasets—1 set of 403 biopsy samples (GSE36059), and 1 set of 16 blood samples (GSE120649), with available rejection phenotyping. To further explore the expression of the ABMR signals in the peripheral blood leukocytes at the time of ABMR, scRNAseq analysis was performed on 12 peripheral blood samples from 6 ABMR and 6 stable patients (4 with biopsy, 2 without biopsy) from Nantes University Hospital, Nantes, France.

Additional information on patient population, sample collection, RNA extraction, histologic diagnosis, RNA-sequencing and differential gene expression, integration with biopsy microarray data, validation in publicly available datasets, and single-cell RNASeq experiments can be found in the Supplementary Methods.

Statistical analysis

Variables with normal distribution are displayed as mean ± SD. Median and range are given for variables not normally distributed. For variance analysis of xCell enrichment scores in different groups, a 2-way analysis of variance was used. To correct for multiple testing, a false discovery rate (FDR) correction was applied during differential gene expression and pathway analysis, and for comparison of xCell enrichment scores. A 2-tailed FDR-corrected P value of <0.05 was considered significant for differentially expressed genes and xCell enrichment scores, and an FDR-corrected P value of < 0.20 was considered significant for Ingenuity Pathway Analysis (IPA) pathway analyses. For all other analyses, a 2-tailed nominal P value of < 0.05 was considered statistically significant. R (version 4.0, R Development Core Team), SAS (version 9.4; SAS Institute), and GraphPad Prism (version 9; GraphPad Software) were used for data analysis and presentation.

Results

Patient and biopsy characteristics

From 365 patients, 384 peripheral blood samples taken at the time of a concomitant allograft biopsy for genome-wide expression analysis were included in a case–control study design, as follows: 136 rejection cases (n = 86 with histology of ABMR [ABMRh], and n = 68 with T cell–mediated rejection [TCMR]; of which n = 18 had concomitant ABMRh and TCMR in the biopsy); 21 polyoma-virus associated nephropathy (PVAN) cases; and 227 heterogeneous control biopsies (absence of rejection and PVAN). Patient demographics and clinical characteristics are provided in Table 1; detailed histologic characteristics are given in Supplementary Table S1. Details on the biopsy samples used for microarray gene expression (n = 224) have been provided previously.

Table 1

Characteristics of the patients and biopsies

Variable	Value
Transplant characteristics (N = 365 patients)
Recipient age at transplantation, yr	49.3 ± 15.5
Recipient age at time of biopsy, yr	51.9 ± 14.6
Recipient sex (male/female)	217/148 (59.4/40.6)
Repeat transplantation	61 (16.7)
Recipient ethnicity (European/Asian/African/other)a	309/6/11/36 (84.7/1.6/3.0/9.9)
Donor age, yra	49.5 ± 16.5
Donor sex (male/female)a	169/182 (46.3/49.9)
Deceased/living donora	279/82 (76.4/22.5)
Heart-beating/non–heart-beating donor	247/32 (67.7/8.8)
Cold ischemia time, ha	12.9 ± 8.5
Biopsy characteristics (N = 384 biopsies)
Indication/protocol biopsy	141/243 (36.7/63.3)
Time after transplantation, d	361 (4 – 12,564)
Time after transplantation, yr
< 1	203 (52.9)
> 1	181 (47.1)
MDRD eGFR, (ml/min per 1.73 m2) at biopsy	41.9 ± 18.3
Proteinuria (g/g creatinine) at biopsya	0.16 (0.0–8.0)
Immunosuppression at biopsy
Cyclosporine	38 (9.9)
Tacrolimus	323 (84.1)
Mycophenolate	54 (85.9)
Azathioprine	7 (1.8)
mTOR inhibitor	43 (11.2)
Corticosteroids	353 (91.9)
Otherb	14 (3.6)
Histologic diagnosis
Control biopsy (no rejection, no PVAN)	227 (59.1)
T-cell mediated rejection (including mixed)
Borderline changes	42 (10.9)
Grade 1 or 2	26 (6.8)
Histology of antibody-mediated rejection (incl. mixed)	86 (22.4)
Mixed rejection	18 (4.7)
Interstitial fibrosis/tubular atrophy grade
0	200 (52.1)
1	70 (18.2)
2	70 (18.2)
3	44 (11.5)
Polyomavirus-associated nephropathy	21 (5.5)
Polyoma viremiaa	60 (15.6)

CIT, cold ischemia time; eGFR, estimated glomerular filtration rate; MDRD, modification of diet in renal disease; mTOR, Mechanistic target of rapamycin; PVAN, polyoma-virus associated nephropathy.

Values are mean (median) ± SD (minimum – maximum), or n (%).

Donor age missing in n=12; CIT missing in n = 14; proteinuria missing in n = 5; donor sex missing in n = 14; ethnicity unknown in n = 3; deceased/living unknown in n = 4; polyoma viremia missing in n = 57.

Other: belatacept (n = 5) / eculizumab (n = 9).

Differential gene expression, pathway enrichment analysis, and upstream regulators

In the blood samples from patients with any rejection (n = 136) versus no rejection (n = 248), ZEB2, DGKH, PATL2, IL18R1, and SLAMF7 were the top 5 most significant differentially expressed genes (DEGs; Figure 1a; Supplementary Table S2). Pathway enrichment analysis on the 319 differentially expressed genes (nominal P < 0.005) revealed glucocorticoid receptor signaling, hepatic fibrosis/hepatic stellate cell activation, and Th1 and Th2 activation pathways as the top 3 canonical pathways (Figure 2a). Additionally, in competitive set enrichment analysis using gene set enrichment analysis (GSEA), the nucleotide-binding oligomerization domain (NOD)-like receptor signaling pathway was the most upregulated (normalized enrichment score [NES] 2.65, P = 0.003), followed by the Janus kinase signal transducer and activator of transcription (Jak-STAT) signaling pathway (NES 2.21, P = 0.003) and Th17 cell differentiation (NES 2.19, P = 0.003; Supplementary Figure S1). The top 5 most contributing genes to the NOD-like receptor signaling pathway were GBP5, STAT1, GBP1, GBP2, and CASP5. Predicted upstream factor analysis uncovered interferon regulatory factors and T cell activating genes as the most activated upstream regulators (Supplementary Table S3; Supplementary Figure S2). In a sensitivity analysis comparing any rejection versus pristine biopsies (without any lesions; n = 110), similar significant DEGs were found (top 10 DEGs are shown in Supplementary Table S4). When looking at differential expression for clinical (n = 78) versus subclinical rejection (n = 58), no DEG showed a significant difference, with adjusted P value < 0.05.

Figure 1

Differentially expressed genes in the peripheral blood for different phenotype comparisons. Red line indicates the false discovery rate (FDR) 0.05 significance level. The most significant features are denoted by name. ABMRh, histology of antibody-mediated rejection; HLA–DSA, human leukocyte antigen–donor-specific antibody; PVAN, polyoma-virus associated nephropathy; TCMR, T cell–mediated rejection.

Figure 2

Top 10 canonical pathways for the different phenotype comparisons obtained using Ingenuity Pathway Analysis. The orange line denotes the false discovery rate (FDR) P value 0.2 level of significance. The top 5 contributing features are shown per pathway. ABMRh, histology of antibody-mediated rejection; ARE, antioxidant response element; CXCR, C-X-C C-X-C chemokine receptor; fMLP, N-Formylmethionyl-leucyl-phenylalanine; GM-CSF, Granulocyte-macrophage colony-stimulating factor; HLA–DSA, human leukocyte antigen–donor-specific antibody; IL, interleukin; ILK, Integrin-linked kinase; iNOS, Inducible nitric oxide synthase; LXR, Liver X Receptor; NF, nuclear factor; MSP-RON, macrophage-stimulating protein - recepteur d'origine nantais; NFAT, Nuclear factor of activated T cells; PKR, protein kinase R; PPAR, Peroxisome proliferator-activated receptor; PVAN, polyoma-virus associated nephropathy; RXR, retinoid X receptor; STAT, signal transducer and activator of transcription; TCMR, T cell–mediated rejection; Th, T helper; VDR, vitamin D receptor.

Comparing ABMRh (n = 86) versus no ABMRh (n = 298), GBP5, STAT1, GBP1, GBP4, and PATL2 were the top 5 most significant DEGs (Figure 1b; Supplementary Table S5). Primary immunodeficiency signaling, protein ubiquitination, and antigen presenting pathways were identified as the top 3 canonical pathways (Figure 2b). Competitive set enrichment analysis with GSEA again identified the NOD-like receptor signaling pathway as the top upregulated pathway (NES 2.71, P = 0.008; Supplementary Figure S3). Interferon-regulatory factors were the most activated upstream regulators (Supplementary Table S6; Supplementary Figure S4). In a sensitivity analysis comparing ABMRh versus pristine biopsies, similar significant DEGs were found (top 10 DEGs are shown in Supplementary Table S7).

Comparing TCMR (n = 68) versus no TCMR (n = 316), IFI27 (down), DGKH, RSAD2 (down), IFI44L (down), and IL18R1 were the top 5 most significant DEGs (Figure 1c; Supplementary Table S8). Hypercytokinemia/hyperchemokinemia in pathogenesis of influenza (down-regulated), interferon signaling (down-regulated), and interleukin (IL)-10 signaling were identified as the top 3 canonical pathways (Figure 2c). Competitive set enrichment analysis using GSEA identified endocytosis as the most upregulated enriched pathways (NES 1.78, P = 0.03), and viral infections were significantly downregulated (Supplementary Figure S5). Interferon regulatory factors were strong inhibitory upstream regulators (Supplementary Table S9; Supplementary Figure S6). In a sensitivity analysis comparing TCMR versus pristine biopsies, few DEGs were found, with only DGKH, MIAT, NKG7, and ZEB2 with FDR < 0.05 (top 10 DEGs are shown in Supplementary Table S10).

When comparing pure ABMRh (n = 68) versus pure TCMR (n = 50; excluding mixed rejection cases), STAT1, PARP14, RSAD2, APOL6, and GBP5 were the most significant DEGs, all upregulated in the ABMRh group. Genes that were more upregulated in TCMR included B and T cell–related genes (CD79A, BLNK, LSP1, TCL1A, and IGHM; Supplementary Table S11).

When we compared human leukocyte antigen–donor-specific antibody (HLA–DSA)–positive samples (n = 100) versus HLA–DSA negative samples (n = 284), CD27 (down), MT-ND1 (down), SEC24A, SLC2A1, and BIRC2 were the top 5 significant DEGs (Supplementary Table S12; Figure 1d). Nuclear factor of activated T cells (NFAT)–signaling pathways were identified as the top canonical pathways (Figure 2d). The interferon-α group, Krüppel-like Factor 3 (KLF3), IL-5, X-box binding protein 1 (XBP1), and interferon-γ (IFN-γ) were the most strongly activated upstream regulators (Supplementary Table S13; Supplementary Figure S7).

Comparing PVAN (n = 21) versus no PVAN (n = 363, including rejection), MT-ND1, MT-ND2, IFI27, MT-ND4, and MT-ATP6 were the top 5 significant DEGs (Supplementary Table S14; Figure 1e). Oxidative phosphorylation and mitochondrial dysfunction were identified as the top canonical pathways (Figure 2e). GSEA demonstrated extracellular matrix–receptor interaction, glycerolipid metabolism, and antigen processing and presentation as the top 3 upregulated pathways (Supplementary Figure S8). INFL1 and DAP3 (involved in mediating IFN-γ–induced cell death) were the most activated upstream regulators (Supplementary Table S15; Supplementary Figure S9). In a sensitivity analysis comparing PVAN versus pristine biopsies, similar significant DEGs were found (top 10 DEGs are shown in Supplementary Table S16).

For polyomavirus viremia (with or without nephropathy; n = 60), compared to absence of viremia (n = 267), LY6E, IFI27, MT-ND2, MT-ND4, and IFI44L were the top 5 significant DEGs (Supplementary Table S17; Figure 1f). Hypercytokinemia/hyperchemokinemia in the pathogenesis of influenza, oxidative phosphorylation, and interferon signaling were identified as the top 3 upregulated canonical pathways (Figure 2f). Interferons were the most activated upstream regulators (Supplementary Table S18; Supplementary Figure S10). In a sensitivity analysis comparing viremia versus control (n = 174, excluding rejection), similar significant DEGs were found (the top 10 are shown in Supplementary Table S19). No peripheral blood transcriptomic differences were found for polyoma viremia with versus without nephropathy after FDR adjustment.

Integration with biopsy microarray data

Next, we compared, per rejection phenotype, the DEGs from blood with the DEGs from 224 biopsies. Integrative meta-analysis revealed that the top 10 overexpressed genes consistent across blood and biopsy samples in any rejection, versus no rejection, were GBP5, CCL4, C1QA, FCGR1B, CRTAM, GBP1, NKG7, SLAMF7, KLRD1, and STAT1 (Figure 3). These genes were consistent across the different ranking methods (Supplementary Table S20). Similarly, the top 10 differentially expressed genes in ABMRh versus no ABMRh, in blood and biopsy samples, were GBP5, CCL4, GBP1, C1QA, CRTAM, FCGR1B, C1QB, GBP4, AIM2, and SLAMF7 (Figure 3; Supplementary Table S20). In TCMR versus no TCMR, C1QA, VSIG4, LCN2, KCNJ2, IL18R1, CCL4, FKBP5, C1QB, CST7, and BCAT1 were the top 10 overexpressed genes when integrating blood and biopsy differential expression data (Figure 3; Supplementary Table S20). Given that only 3 PVAN cases were present in the biopsy dataset, no integration was performed for the PVAN phenotype. The top consistently enriched upregulated pathways from GSEA across the different rejection phenotypes in blood and biopsies are represented in Supplementary Table S21.

Figure 3

The top 10 overexpressed features across biopsy microarray data and peripheral blood RNAseq data for any rejection, histology of antibody-mediated rejection (ABMRh), and T cell–mediated rejection (TCMR). The scatter plot shows the distribution of the individual genes, supplemented with a violin plot demonstrating overall distribution with included median and interquartile range.

Cellular origin of consistently overexpressed genes in the kidney allograft

To assess the cellular origin of these consistently overexpressed genes in the kidney allograft, we used publicly available single-cell data from kidney allograft biopsies. Briefly, scRNASeq was performed on 2 kidney allograft biopsies with ABMR, and 5 healthy surveillance biopsies. After quality control and filtering was completed, 33,216 cells were detected, and unsupervised clustering revealed 14 clusters corresponding to the main kidney cell subtypes but also infiltrating immune cells (Supplementary Figure S11A–C). We evaluated the expression of genes corresponding to the top 10 consistently overexpressed genes in ABMR from both blood and biopsy (see above, GBP5, CCL4, GBP1, C1QA, CRTAM, FCGR1B, C1QB, GBP4, AIM2, and SLAMF7) across the kidney structural cells and infiltrated immune cells (Supplementary Figure S11D). The genes of interest were mainly expressed by the infiltrated leukocytes in the kidney—that is, the antigen-presenting cells and lymphocytes.

Cell-type enrichment analysis

We found that the peripheral blood signals were consistent with the biopsy signals at the time of rejection. Moreover, these differentially expressed molecules were expressed by the leukocyte pool in the allograft. Therefore, we hypothesized that transcriptomic changes in the peripheral blood reflect changes in circulating immune cells. To further investigate this possibility, we conducted a cell enrichment analysis using xCell. Cell enrichment scores did not discriminate among the different phenotypes (Figure 4). A comparison of the different immune cell types across the different phenotypes showed few significant differences between enrichment scores after FDR adjustment (Supplementary Figure S12). In mixed rejection versus control, there was an overrepresentation of Th2 cells and an underrepresentation of Th1 cells. In PVAN versus control, basophils and Th1 cells were overrepresented, and Th2 cells were underrepresented (Supplementary Figure S12). In contrast, cell enrichment scores based on the biopsy transcriptomic data did discriminate rejection phenotypes from no rejection (Figure 4). In the biopsy data, the most significantly overrepresented cell types in the rejection phenotypes were basophils, CD4+ and CD8+ central and effector memory T cells, dendritic cells, monocytes, NK cells, and Th2 cells. In contrast, naïve CD4+ and CD8+ T cells, and T regulatory cells were significantly underrepresented in all rejection phenotypes (Supplementary Figure S12). The significant differences between pure ABMRh and pure TCMR included the overrepresentation of NK cells in ABMRh and the overrepresentation of B-cells, CD4+ T cells (including CD4+ memory T cells), and CD8+ central memory T cells in TCMR (FDR P < 0.05).

Figure 4

Cell enrichment scores in the different phenotypes of kidney allograft pathology as obtained from blood (upper panel) and biopsy transcriptomic data. Heatmaps of the cell enrichment scores, ordered by hierarchical clustering. In the blood samples, cell enrichment scores did not discriminate among the different phenotypes, whereas in the biopsy samples, discrimination of rejection versus no-rejection phenotypes based on their cell enrichment scores in immune cells is seen. ABMRh, histology of antibody-mediated rejection; CD, cluster of differentiation; DC, dendritic cell; HLA–DSA, human leukocyte antigen–donor-specific antibody; NK, natural killer cell; NKT, natural killer T cell; PVAN, polyoma-virus associated nephropathy; Tcm, Tem, TCMR, T cell–mediated rejection; Th, T helper cell.

Validation of ABMR signals in independent cohorts

Finally, we validated observed ABMR signals from this study in independent cohorts: one large biopsy microarray cohort (n = 403; GSE36059) and a small but well-phenotyped blood RNAseq cohort (n = 16; GSE120649). We started from the top 10 consistently overexpressed genes from both blood and biopsy in ABMR—that is, GBP5, CCL4, GBP1, C1QA, CRTAM, FCGR1B, C1QB, GBP4, AIM2, and SLAMF7. The 10 genes of interest were overexpressed in ABMR in both cohorts (Supplementary Figure S13). All 10 genes were differentially expressed at the adjusted P value level of < 0.05 in the biopsy cohort, whereas in the smaller blood dataset, not all 10 genes reached significance, probably owing to low statistical power (Supplementary Table S22).

To further explore the expression of these genes in the peripheral blood leukocytes at the time of ABMR, scRNAseq analysis was performed on 12 peripheral blood samples (6 ABMR and 6 stable patients). After quality control and filtering was completed, along with removal of 2 clusters containing only doublets, 69,127 cells were detected. Unsupervised clustering revealed 15 clusters corresponding to the main myeloid and lymphoid cells and granulocytes/platelets (Figure 5a–c). We evaluated the expression of genes corresponding to the top 10 consistently overexpressed genes in ABMR from both blood and biopsy across the peripheral blood cells (Figure 5d). The genes of interest were mainly expressed by the monocytes, NK cells, T cells, and B cells from the ABMR condition. When performing differential expression for ABMR versus no ABMR per cell type, significantly differentially expressed genes in ABMR in nonclassical monocytes included C1QA, C1QB, GBP5, GBP1, GBP4, and FCGR1B; in T cells, GBP5, GBP1, and GBP4; and in B cells, AIM2 (Figure 5e).

Figure 5

Overview of the single-cell RNA-sequencing analysis on 12 peripheral blood samples with and without antibody-mediated rejection (ABMR) to determine the cellular origin of the ABMR signals in the blood. Briefly, scRNAseq was performed on 6 peripheral blood samples from kidney transplant recipients with a concomitant diagnosis of ABMR, and 6 stable patients without ABMR. After quality control and filtering were completed, along with removal of 2 clusters containing only doublets, 69,127 cells were detected. (a–c) Unsupervised clustering revealed 15 clusters corresponding to the main myeloid and lymphoid cells and granulocytes/platelets. (d) We evaluated the expression of genes corresponding to the top 10 consistently overexpressed genes in ABMR from both blood and biopsy (GBP5, CCL4, GBP1, C1QA, CRTAM, FCGR1B, C1QB, GBP4, AIM2, and SLAMF7) across the peripheral blood cells. (e) When performing differential expression for ABMR versus no ABMR per cell type, significantly differentially expressed genes in ABMR in nonclassical monocytes included C1QA, C1QB, GBP5, GBP1, GBP4, and FCGR1B; in T cells, GBP5, GBP1, and GBP4; and in B cells, AIM2. CD, cluster of differentiation; cDC, conventional dendritic cell; IL, interleukin; MT, mitochondrial transcripts; NK, natural killer cell; NKT, natural killer T cell; PBMC, peripheral blood mononuclear cell; pDC, plasmacytoid dendritic cell.

Discussion

In this study, we report the transcriptomic changes in peripheral blood at the time of kidney transplant rejection and polyomavirus-associated nephropathy. The gene expression changes observed in blood partly mirror those observed in kidney transplant biopsies, although the latter are more pronounced, reflecting immune cell infiltration in addition to activation of resident cells. The differentially expressed genes in blood of cases with rejection reflect immune activation pathways, with interferons as the most likely upstream regulators, which are very sensitive to corticosteroid therapy but nonetheless seem to escape our current immunosuppressive armamentarium. Analysis per rejection subtype reveals that it is ABMRh that associates with these activated immune pathways, whereas the TCMR signals are nonspecific. Polyomavirus nephropathy also associates with transcriptional changes in interferon-regulated genes, in addition to mitochondrial gene overexpression. Although the composition of the immune cell infiltration in kidney transplant biopsies is being charted,, and overlapping rejection genes from both blood and biopsy are expressed in the leukocytes present in the kidney allograft biopsy, our transcriptomic analyses do not suggest major changes in the immune cell composition in peripheral blood. Overall, we identified the key genes and pathways related to immune activation after kidney transplantation, despite the use of strong immunosuppressants.

We included a heterogeneous control group, containing varying degrees of fibrosis, with inclusion of low-grade inflammatory scores not reaching the threshold for definition of rejection. This makeup reflects the heterogeneity of samples encountered in biopsies taken in clinical practice, rather than highly selected clear-cut samples, and attenuates confounding factors (such as kidney dysfunction and time after transplantation) in the comparisons. Also, the gene expression seems to be related to the underlying process, independent of the clinical presentation, as illustrated by the absence of differences in clinical versus subclinical rejection. Similarly, we compared samples with the phenotype versus without the phenotype, including other diagnoses, as this is most representative for samples encountered in clinical practice, yet it could have impacted the statistical power to detect signals in more severely inflamed cases, for example, of TCMR. Therefore, we also included sensitivity analyses comparing each phenotype against pristine controls.

To our knowledge, this is the largest discovery set of biopsy-paired peripheral blood samples analyzed with RNA-sequencing, studying both rejection phenotypes and polyomavirus infection. The enrichment of these phenotypes allowed for sufficient power to detect robust and biologically relevant signals. The strong enrichment of immune pathways observed in rejection cases in our study contrasts with a previous RNAseq study on 37 peripheral blood samples, which primarily identified noncoding genes and non–immune pathway enrichment in cases with ABMR. The reason for the difference in results between the present study and this previous study could lie in the larger power of our study, with enrichment of the phenotypes of interest leading to stronger detectable signals. Nevertheless, even with this enriched and large cohort, fold changes are quite small in the blood DEGs, likely related to the heterogeneity in phenotypes and the weaker signals in blood.

The NOD-like receptor signaling pathway was one of the most upregulated and enriched pathways identified in peripheral blood of rejection cases, more specifically for ABMRh, as we also observed in kidney transplant biopsies. The guanylate binding protein (GBP) genes and interferon-related genes are the strongest contributors to the upregulation of this pathway. In our previous transcriptomic study using gene expression microarray in peripheral blood samples, the GBP genes were also identified as the top discriminating genes for ABMRh, and were proven to be valid in a noninvasive signature for this phenotype. NOD-like receptors are important intracellular receptors that initiate innate immune recognition but can also activate adaptive immune responses by identifying pathogen-associated molecular patterns. They are directly linked to the nuclear factor (NF)-κB pathway and the inflammasome, through which they lead to the release of various proinflammatory cytokines. These pattern-recognition receptors can recognize non–self-components and respond to various forms of cell stress. Activation of these receptors by missing of “self” has been described in mice; however, in humans, the role of these cytosolic pattern-recognition receptors in allograft rejection has not been clarified. We also identified other pattern-recognition receptor- like toll-like receptors and retinoic acid inducible gene I (RIG-I)-like helicases among the top enriched canonical pathways in rejection cases.

Next, we studied the transcriptomic differences in patients with versus without HLA–DSA, irrespective of histologic diagnosis. Clinical studies have demonstrated worse outcome in patients with HLA–DSA, even in the absence of histologic lesions in the biopsy, whereas transcriptomic studies observed increased intrarenal ABMR-associated transcripts in these patients., Here, we observed peripheral blood upregulation of calcineurin-NFAT–associated genes such as CALM1 (Calmodulin 1), CHP1 (calcineurin homologous protein), and PPP3R1 (Calcineurin B, Type I) in patients with HLA–DSA. The strong upregulation of this lymphocyte-activating pathway in patients with HLA–DSA might result in insufficient therapeutic action of the calcineurin inhibitor therapy, given to >90% of patients in this study. This potential imbalance between the activated lymphocyte-activating pathway and our targeted immunosuppressive therapy could suggest that higher doses of calcineurin inhibitors or additional therapy are needed in patients with HLA–DSA, to prevent the deleterious lymphocyte activation associated with HLA–DSA. Further confirmation of the relevance of this pathway in patients with HLA–DSA is needed from other transcriptomic studies, as well as clinical studies investigating the impact of therapy (and nonadherence?) on this pathway.

One of the major upstream regulators predicted in these cases, XBP1, is required for the transcription of a subset of class II major histocompatibility genes. XBP1 expression is controlled by the cytokine IL-4 and the antibody Immunoglobulin Heavy Constant Mu (IGHM), and in turn, XBP1 controls the expression of IL-6, which promotes plasma cell growth and production of immunoglobulins in B lymphocytes. XBP1 has been implicated in human cancers and metabolic diseases, for which potential therapeutics have been described.,

Although ABMRh is associated with clear immune activation in peripheral blood, this was not the case for TCMR. In TCMR, many of the significant differentially expressed genes were downregulated, reflected by downregulation of the enriched pathways, such as hypercytokinemia and interferon signaling, and rather unspecific upregulated pathways. The absence of biologically relevant upregulated pathways in peripheral blood in TCMR was also seen in previous whole transcriptome studies., The reason these signals are weaker and unspecific in TCMR is not clear, other than being due to attenuation by T cell–targeted immunosuppression. Also, the heterogeneity of this group, comprising borderline changes, tubulo-interstitial inflammation, and intimal arteritis, might contribute to these less-specific findings. A majority of cases in the TCMR group were borderline changes (threshold t > 0 and i > 0), reflecting also the lower prevalence of more severe forms of TCMR in current clinical practice, again related to the use of T cell–targeted immunosuppressive drugs that largely suppress the clinical phenotype of TCMR. In comparing pure ABMRh with pure TCMR cases, we found overexpression of some B cell transcripts in TCMR. We could speculate that these B-cell transcripts are related to the modulating role of B cells for T-cell activation through antigen presentation, costimulation, or production of cytokines, independently from antibody production., Possibly, the B-cell activation from this interaction during TCMR then leads to plasma cell formation, antibody production, and eventually ABMR, which is mirrored in the clinical observations of arising de novo DSA and ABMR after previous TCMR episodes.24, 25, 26

From the integrative analysis with rejection signals from the biopsy, we extracted consistently upregulated genes and pathways in both blood and biopsy per rejection phenotype. Given that the biopsy signals are more local to the injury compared to the blood signals, we believe this integration further strengthens and uncovers robust and specific signals. The overexpression of these integrated signals in ABMR was confirmed in external datasets, both from biopsy and blood and from different platforms (microarray, RNASeq, and small conditional [sc]RNASeq).

The infiltration and composition of immune cells in kidney transplant biopsies at the time of rejection have been reported before, at the transcriptomic level, and using multiplex imaging.,,27, 28, 29 Whether the immune cell composition in peripheral blood changes in parallel was unclear. From the scRNASeq analysis, we found that the genes that were overexpressed in both blood and biopsies from ABMRh cases were found to correspond to expression in the leukocytes present in the kidney allograft. Next, using immune cell enrichment analysis in biopsy samples, we confirmed and expanded on previous observations, showing that the most significantly overrepresented cell types in rejection were basophils, CD4+ and CD8+ central and effector memory T cells, dendritic cells, monocytes, NK cells, and Th2 cells. In ABMRh cases, we observed overrepresentation of NK cells, whereas in TCMR, there was overrepresentation of B cells, CD4+ T cells (including CD4+ memory T cells), and CD8+ central memory T cells. Naïve CD4+ and CD8+ T cells, and T regulatory cells, were significantly underrepresented in all rejection phenotypes. These solid observations strengthen our approach of using immune cell enrichment analysis in biopsy specimens; yet, this approach did not yield clear differences in immune cell enrichment scores inferred from peripheral blood transcriptomic data. The intrarenal immune cell composition thus does not seem to be mirrored by alterations in peripheral blood cell populations. This discrepancy could be due to confounders influencing the cell populations in the peripheral blood, such as systemic infections. Also, in contrast to the graft tissue, these cells are always present in the peripheral blood and do not infiltrate solely upon inflammation. Moreover, others have previously suggested that there are compartment-specific differences between the peripheral blood leukocytes in the circulation and the subset of leukocytes that are activated and recruited to the allograft at the time of rejection. Finally, xCell reports abundances in cell types to allow comparisons between conditions but is not fully equipped to study the activation changes of these cells. As an illustration of these activation changes, in scRNASeq analysis of peripheral blood, we found differentially expressed genes in peripheral leukocyte cell types in ABMR versus no ABMR, suggesting changes in expression and activation state within the same cell type upon rejection.

The presence of polyomavirus viremia, with or without evidence of PVAN, was associated with clear upregulation of interferon-regulated genes, and remarkably, mitochondrial genes. Related to these mitochondrial genes, oxidative phosphorylation and mitochondrial dysfunction were identified among the top upregulated canonical pathways. Recently, the polyomavirus agnoprotein was shown to target mitochondria and modulate their functions. The BK polyomavirus disrupts the mitochondrial network and membrane potential when expressing the 66aa-long agnoprotein in the late viral life cycle. This agnoprotein impairs nuclear interferon regulatory factor 3 (IRF3)-translocation, impairs interferon-beta expression, and promotes mitophagy, thereby allowing the polyomavirus to evade innate immune sensing. This mechanism could explain mitochondrial injury in tubular cells, which might be reflected in the peripheral blood but will need confirmation through gene expression analyses on kidney allograft tissue.

Our study has several limitations. For instance, our population was almost universally treated with a calcineurin inhibitor–based immunosuppressive regimen, and the majority of our population was of Caucasian ethnicity, which might limit generalizability. Moreover, the diagnostic criteria for rejection after kidney transplantation are a topic of active discussion, making our reference standard of histology imperfect. Also, this study consists of selected cases of rejection phenotypes, PVAN, and control cases from per-protocol and for-cause biopsies, not representing the real-life prevalence. No information on HLA–DSA specificity was available in the centralized database, nor did we have available information on previous rejection treatments. Downstream analysis and validation were confined to ABMR signals based on the more robust signals in this phenotype and the availability of this phenotype in other datasets.

In conclusion, using peripheral blood RNASeq data, we uncover the immune activation in the peripheral blood leukocytes at the time of kidney allograft pathology, bypassing our current strong immunosuppressant armamentarium, which mirrors the molecular changes in the kidney allograft and provides a framework for future therapeutic interventions.

Disclosure

All the authors declared no competing interests.

Data Statement

RNA-sequencing data have been deposited at the National Institutes of Health Gene Expression Omnibus at http://www.ncbi.nlm.nih.gov/geo under the series accession number GSE175718 together with clinical metadata. Single-cell RNA-sequencing data from peripheral blood samples have been deposited at the ArrayExpress database (http://www.ebi.ac.uk/arrayexpress) under the series accession number E-MTAB-11450. Other individual deidentified participant data that underlie the results reported in this article (text, tables, figures, and appendices) can be made available on a collaborative basis following institutional review board approval. Requests should be directed to the corresponding author.