Exploring the mechanism of BK polyomavirus-associated nephropathy through consensus gene network approach

Abstract

BK polyomavirus-associated nephropathy occurs in kidney transplant recipients under immunosuppressive treatment. BK polyomavirus is implicated in cancer development and invasion, and case reports of renal cell carcinoma and urothelial carcinoma possibly associated with BK polyomavirus has been reported. Further, it has been suggested that the immune responses of KT-related diseases could play a role in the pathogenesis and progression of renal cell carcinoma. Thus, we thought to examine the relationship between BK polyomavirus-associated nephropathy and renal cell carcinoma in terms of gene expression. To identify the common and specific immune responses involved in kidney transplantation-related diseases with a specific focus on BK polyomavirus-associated nephropathy, we performed consensus weighted gene co-expression network analysis on gene profile datasets of renal biopsy samples from different institutions. After the identification of gene modules and validation of the obtained network by immunohistochemistry of the marker across kidney transplantation-related diseases, the relationship between prognosis of renal cell carcinoma and modules was assessed. We included the data from 248 patients and identified the 14 gene clusters across the datasets. We revealed that one cluster related to the translation regulating process and DNA damage response was specifically upregulated in BK polyomavirus-associated nephropathy. There was a significant association between the expression value of hub genes of the identified cluster including those related to cGAS-STING pathway and DNA damage response, and the prognosis of renal cell carcinoma. The study suggested the potential link between kidney transplantation-related diseases, especially specific transcriptomic signature of BK polyomavirus associated nephropathy and renal cell carcinoma.

Fig 1. The results of the principal component analysis.

The results of the principal component analysis of the samples for each dataset are presented. The x and y-axis represent the principal component one and two respectively. The color indicates the disease category, and the 95% confidence ellipses are drawn for each category. IFTA, interstitial fibrosis and tubular atrophy; Normal, normal biopsy; BKPyVAN, BK polyomavirus-associated nephropathy; BKPyVB, normal biopsy findings with BKPyV DNAemia; AR, acute rejection; TCMR, T-cell mediated rejection; CAN, chronic allograft nephropathy.

Fig 2. The relationship between module eigengenes and clinical conditions.

(a) The relationship between clinical conditions and the module eigengenes is depicted for each dataset. The value and color of the cell of the heatmap represent coefficients of the linear models, and the asterisk mark indicates significance. Note that the range of the color value is specific to each dataset. The modules are sorted according to the number of contained genes. The color column represents the module color assigned by WGCNA. (b) The subset of significantly enriched pathways in specifically upregulated module in BK polyomavirus-associated nephropathy. The x-axis represents the number of genes that belong to the pathway. The color represents the values of the raw p-value.

Fig 3. The visualization of sub-network of the module related to BKPyVAN.

The sub-network of the top-5 scored hub gene and genes with the order of one from the gene within the corresponding module are visualized. The edges with the weight below the 99.9^th percentile and the nodes with zero degree were discarded beforehand. Additionally, only the edges with the weight above the 90^th percentile are included in the visualization. The color represents the consensus kME value. The edge width represents the strength of the connection between nodes. The node size represents the average normalized expression values of BK polyomavirus-associated nephropathy patients in GSE72925.

Fig 4. The summary of the survival analysis of the renal cell carcinoma.

The Kaplan-Meier plots of the overall survival and the representative expression value of each significant module using genes with the top-5 kME values. The p-values of log rank test are shown. The density plot represents the distribution of p-values calculated by the same method using randomly selected genes with the same number. The solid red line indicates the p-value of the gene signature, and the dashed line indicates p-value of 0.05.