Confirmation of the predictive function of cuproptosis-related gene FDX1 in clear cell renal carcinoma using qRT-PCR and western blotting

Expression of FDX1

The expression of FDX1 in tumors was significantly lower than that in adjacent tumor tissues, including BRCA, CHOL, COAD, KICH, KIRC, KIRP, LUAD, LUSC, PCPG, READ, and THCA (Figure 1A). In addition, the expression of FDX1 in tumors was lower than that in adjacent tumor samples from TCGA (KIRC), International Cancer Genome Consortium (ICGC) (RECA-EU), Gene Expression Omnibus (GEO) (GSE66272), and ArrayExpress (E-MTAB-3267) databases (Figure 1B, 1D–1F). We also obtained the same results when comparing the same patient in pairs (Figure 2A, 2C, 2D). Furthermore, we found that the expression of FDX1 in adjacent tumor specimens was significantly higher than that in tumors after performing qRT-PCR in our independent clinical database (NBUNH), regardless of whether it was a pair or discrete (Figure 1C, 2B). In addition, compared with HK-2 cells, the mRNA and protein expression levels of FDX1 were lower in 786-O and OS-RC-2 cells (Figure 3A, 3B).

Clinical analysis of FDX1

In this study, 532 patients with ccRCC were divided into high- and low-expression groups according to the median FDX1 expression level. Kaplan-Meier curves showed that the overall survival (OS) and progression free survival (PFS) in the low expression group were lower than those in the high expression group (P < 0.05, Figure 4A, 4B). However, there was no significant difference in the survival time in the NBUNH cohort (Figure 4C). In the E-MTAB-1980 dataset, when the best cutoff FDX1 expression level was selected, the OS in the low-expression group was significantly lower than that in the high-expression group (P < 0.01, Figure 4D). With the help of ROC (receiver operation curve) curves, we found that AUC (area under the curve) value of FDX1 was 0.966 (95% CI 0.943–0.982) in TCGA dataset (Figure 4E), 0.623 (95% CI 0.532–0.712) in the NBUNH cohort (Figure 4F), 0.979 (95% CI 0.948–0.999) in the ICGC dataset (Figure 4G), and 0.985 (95% CI 0.957–1.000) in the GEO dataset (Figure 4H).

The expression of FDX1 was lower in Stage IV than in Stage I in TCGA dataset (P < 0.05, Figure 5A), and the same was observed in grade 4 compared to grade 2 (P < 0.001, Figure 5B). There was no statistically significant difference between T1–2 and T3–4 (P > 0.05, Figure 5C). Our database (NBUNH) showed that the expression level of FDX1 in Stage II was higher than that in Stage I, G2 was higher than that in G1, and T2 was higher than that in T1 (all P < 0.05, Figure 5D–5F). The correlations between FDX1 expression levels and clinical characteristics from TCGA and NBUNH databases are shown in Table 1. Univariate and multivariate Cox analyses showed that FDX1 could be regarded as an independent prognostic indicator of OS in ccRCC (Figure 5G, 5H).

Function analyses of FDX1

There were 146 genes correlated to FDX1 considering the correlation coefficient >0.5 (Supplementary Table 1). The co-expression circle diagram shows the correlation between FDX1 and the other 11 genes with the largest absolute values of the correlation coefficient (Figure 6A). Gene Ontology (GO) categories included biological processes (BP), cellular components (CC), and molecular functions (MF) (Supplementary Table 2). We discovered that BP mainly contained cellular respiration, CC mainly included the mitochondrial inner membrane, and MF mainly contained the proton transmembrane (Figure 6B and 6C). The circular diagram depicts the top five GO functions (Figure 6D). Oxidative phosphorylation was closely associated with FDX1 expression in Kyoto Encyclopedia of Genes and Genomes (KEGG) (Figure 6E, 6F, Supplementary Table 3). Pathway graphs suggested that the genes correlated with the altered expression of oxidative phosphorylation (Figure 6G).

Altogether, 234 differentially expressed genes (DEGs) were screened according to the cutoff values (Supplementary Table 4). A heatmap illustrating the expression of the top 20 DEGs is shown in Figure 7A. GO and KEGG analyses were performed to explore their functions (Supplementary Tables 5 and 6). We mainly found monovalent inorganic cation homeostasis in BP, the apical part of the cell in CC, and the anion transmembrane in MF (Figure 7B, 7C). On the other hand, collecting duct acid secretion and synaptic vesicle cycle of differential genes were enriched in KEGG pathways (Figure 7E, 7F), which was analyzed in “c5.go.v7.4. symbols.gmt” (Supplementary Table 7). The alpha amino acid metabolic process was enhanced when FDX1 was expressed at different levels (Figure 7D). Oxidative phosphorylation and the peroxisome were main pathways after GSEA analysis via “c2.cp.kegg.v7.4. symbols.gmt” (Figure 7G and Supplementary Table 8).

FDX1 methylation, expression verification, and tumor mutational burden (TMB)

When mRNA is transcribed into a protein, it is modified by methylation in the transcription process. In the present study, FDX1 methylation was explored for its prognostic value in KIRC using MethSurv analysis. A DNA methylation heatmap illustrated that the highest FDX1 methylation level was in cg06674932 (Figure 8A). Overall, we identified eight CpGs in FDX1 that were significantly associated with ccRCC prognosis (Supplementary Table 9). TMB analysis demonstrated that FDX1 expression was negatively correlated with TMB levels (Spearman, R = −0.13, P = 0.019, Figure 8B). FDX1 was mainly expressed in the proximal and distal tubules of normal renal tissues but not in tumor tissues (Figure 8C).

Immune infiltration and tumor microenvironment (TME) analyses of FDX1

The CIBERSORT method explored 22 types of immune cells, of which five were found to be significantly different (P < 0.05, Figure 9A). In addition, correlation analysis showed that FDX1 was correlated with three types of immune cells (Figure 9B), suggesting that FDX1 could affect immune responses by regulating resting mast cells, resting NK cells, and regulatory T cells (Tregs). TME analysis indicated that stromal, immune, and ESTIMATE scores increased in the FDX1 low-expression group (P < 0.001; Figure 9C).

We evaluated the differences in common immune cells (ICs) between the FDX1 high- and low-expression groups, including programmed cell death 1 (PD1/PDCD1), programmed cell death ligand 1 (PDL1/CD274) and cytotoxic T lymphocyte antigen 4 (CTLA4). Our results showed that most ICs were upregulated (Figure 10A) in the low-expression group. Furthermore, a higher Tumor Immune Dysfunction and Exclusion (TIDE) score was obtained in the low-expression group, demonstrating stronger immune dysfunction and immune resistance (Figure 10B–10E).

Data collection and patient recruitment

RNA-seq and clinical information of patients with ccRCC were downloaded from TCGA database (KIRC; N = 72, T = 532; https://portal.gdc.cancer.gov/). In addition, the ICGC (RECA-EU; N = 45; T = 91; http://dcc.icgc.org), ArrayExpress (E-MTAB-1980; T = 101; E-MTAB-3267; N = 6, T = 53; https://www.ebi.ac.uk/arrayexpress/) and GEO (GSE66272; N = 25; T = 26; https://www.ncbi.nlm.nih.gov/geo/) databases were searched for external validation. All data were preliminarily processed using the “limma” R package. Details of the clinical information from each database are shown in Table 2.

This study was approved by the Ethics Committee of the NBUNH and written informed consent was obtained from all included patients. A cohort of 76 tissue samples was collected from patients with ccRCC at the NBUNH. This cohort included patients with primary ccRCC who underwent radical nephrectomy in the Department of Urology since 2010. One sample was not used because its RNA was degraded. All clinical characteristics were obtained from the electronic information system of the hospital. Details included the initial age at diagnosis, sex, stage, and grade. Follow-up data were collected via telephone and Ningbo residents’ health records. Detailed information on the 75 ccRCC patients is shown in Supplementary Table 10.

Western blot, qRT-PCR, and cell cultures

Cell lysates were harvested in RIPA buffer (Solarbio, Beijing, China) containing 1% PMSF protease inhibitor (Solarbio). Total protein concentration was calculated using a BCA protein assay kit (Beyotime, Beijing, China). Total protein (30 μg) samples were loaded and separated using 12% SDS-PAGE, transferred to PVDF membranes, blocked with 5% non-fat dry milk, and then incubated overnight with diluted primary antibody against FDX1 (Proteintech, Wuhan, China) or β-actin (Proteintech) at 4°C overnight. The blots were then washed with TBST, incubated with horseradish peroxidase-labeled secondary antibodies (Boster, Wuhan, China), and visualized using an enhanced chemiluminescence reagent.

Total RNA was extracted from the clinical samples and renal cancer cells using an RNA extraction kit (ServiceBio, Wuhan, China). The cDNA was synthesized according to the manufacturer’s instructions using a Servicebio^® RT First Strand cDNA Synthesis Kit (ServiceBio). The qRT-PCR used 2*SYBR Green qPCR Master Mix (ABclone, USA) for real-time quantitative PCR reaction. Overall, there were 40 cycles of denaturation at 95°C for 15 s and annealing/stretching at 60°C for 30 s. The following primers were used for qRT-PCR: FDX1, forward: 5′-CCACTTTATAAACCGTGATGGTG-3′; reverse: 5′-ACATGCACCAAAGCCATCAA-3′. GAPDH, forward: 5′-GGAAGCTTGTCATCAATGGAAATC-3′; reverse: 5′-TGATGACCCTTTTGGCTCCC-3′. GAPDH levels were used to standardize the data. The relative mRNA level of FDX1 was calculated by the 2^−ΔΔCt method.

The HK-2, OS-RC-2, and 786-O cell lines were purchased from the Cell Bank of the Chinese Academy of Sciences (Shanghai, China). HK-2 cells were cultured in DMEM (Hyclone, Logan, Utah, USA), and OS-RC-2 and 786-O cells were cultured in RMPI-1640 medium (Hyclone). All cells were incubated at 37°C in 5% CO₂ after supplementing the culture medium with 10% heat-inactivated fetal bovine serum (Hyclone, Auckland, NZ), 100 U/mL streptomycin, and 100 mg/mL penicillin (Hyclone, Logan, UT, USA).

Clinical survival analysis

TIMER 2.0 (http://timer.cistrome.org) was used to determine the expression levels of FDX1 in the tumors. The “ggplot2” and “ggpubr” R packages were applied to draw the box graph and compare the expression levels of FDX1 in tumor and adjacent normal samples. The R packages of “survival” and “survminer” were used to draw survival curves and analyze the different survival results in two groups with high- and low-expression levels according to the medium. The Receiver operating characteristic (ROC) curves were generated by the “pROC” R package, and the area under curve (AUC) values were calculated to evaluate the specificity and sensitivity of FDX1 in predicting benign and malignant tumors.

The clinical characteristics of the FDX1 high- and low-expression groups were compared in TCGA and NBUNH cohorts (details in Table 1). The clinicopathological parameters in high- and low-expression groups were compared using the “ggpubr” R package. Univariate and multivariate Cox regression analyses were performed using the Kaplan–Meier “survival” R package to assess the independence of FDX1 from other clinical factors.

Analysis of the function and pathways

We screened out genes correlated with FDX1 using R software, meeting the requirements of |Pearson correlation coefficient|>0.5 and P < 0.001. The R packages “circlize” and “corrplot” were used to visualize the co-expression results. The R packages “clusterProfiler”, “org.Hs.eg.db”, “dplyr”, “enrichplot”, “ggplot2”, “circlize”, “RColorBrewer”, “ComplexHeatmap”, “R.utils”, and “pathview” were applied for GO and KEGG analyses. GO circle graphs and KEGG path diagrams were constructed.

DEGs in groups with high and low expression levels of FDX1 were screened using R software, meeting the criteria |logFC| >2 and FDR < 0.05. The R package “limma” was used to identify differential genes, and then “pheatmap” visualized the top 20 DEGs. GO and KEGG analyses of the DEGs were performed using the R packages “clusterProfiler”, “org.Hs.eg.db”, “enrichplot”, and “ggplot2.” GSEA enrichment analysis was done in MSigDB gene sets through “c5.go.v7.4. symbols.gmt” and “c2.cp.kegg.v7.4. symbols.gmt” in the R software.

Methylation and TME analyses

MethSurv (https://biit.cs.ut.ee/methsurv/) was used to evaluate the prognostic value of FDX1 methylation in patients with ccRCC. The Protein Atlas Database (https://www.proteinatlas.org/) was used to display the protein expression level of FDX1. The TMB of FDX1 was calculated based on gene mutation data from TCGA.

The CIBERSORT algorithm revealed a relationship between FDX1 expression and 22 types of immune cells. The ESTIMATE algorithm was used to evaluate the immune microenvironment (ImmuneScore, robustness score, ESTIMATEScore, and TumorPurity) in groups with high and low FDX1 expression. The R packages “limma”, “estimate”, “e1071”, “reshape2”, “vioplot”, “ggExtra”, and “ggpubr” were used to complete the analyses. The box diagram and violin graph showed the results of immune cell infiltration and immune microenvironment. The correlation between FDX1 expression and 22 kinds of immune cells was illustrated by Lollipop graph.

The TIDE database was used to evaluate therapeutic effects. The expression matrix was uploaded to the website (http://tide.dfci.harvard.edu) to predict the possible immunotherapy effects with different expression levels of FDX1.

Statistical analyses

A t-test was used to analyze the differences between groups of variables with a normal distribution. Otherwise, the Mann-Whitney U test was applied. A chi-square test was used for quantitative comparisons. The Pearson correlation method was used to analyze the correlation between two different genes. All statistical analyses were performed using the R software (v4.1.1) (https://www.r-project.org/). P < 0.05 was considered to be statistically significant. We marked ^* in the results, where ^* represents P < 0.05, ^** represents P < 0.01, and ^*** represents P < 0.001.

Institutional review board statement

The study was conducted in accordance with the Declaration of Helsinki and approved by the Institutional Ethics Committee of Ningbo Urology and Nephrology Hospital.

Data availability statement

The data supporting the findings of this study are available from KIRC at https://portal.gdc.cancer.gov/, ICGC at http://dcc.icgc.org, ArrayExpress at https://www.ebi.ac.uk/arrayexpress/, and GEO at https://www.ncbi.nlm.nih.gov/geo/. The authors confirm that the data supporting the findings of this study are available in the article and supplementary material.