A Two-DNA Methylation Signature to Improve Prognosis Prediction of Clear Cell Renal Cell Carcinoma.

PURPOSE: Effective biomarkers and models are needed to improve the prognostic prospects of clear cell renal cell carcinoma (ccRCC). The purpose of this work was to identify DNA methylation biomarkers and to evaluate the utility of DNA methylation analysis for ccRCC prognosis.
MATERIALS AND METHODS: An overview of genome-wide methylation of ccRCC tissues derived from The Cancer Genome Atlas (TCGA) database was download for analysis. DNA methylation signatures were identified using Cox regression methods. The potential clinical significance of methylation biomarkers acting as a novel prognostic markers was analyzed using the Kaplan-Meier method and receiver operating characteristic (ROC) curves.
RESULTS: This study analyzed data for 215 patients with information on 23171 DNA methylation sites and identified a two-DNA methylation signature (cg18034859, cg24199834) with the help of a step-wise multivariable Cox regression model. The area under the curve of ROCs for the two-DNA methylation signature was 0.819. The study samples were stratified into low- and high-risk classifications based on an optimal threshold, and the two groups showed markedly different survival rates. Moreover, the two-DNA methylation marker was suitable for patients of varying ages, sex, stages (I and IV), and histologic grade (G2).
CONCLUSION: The two-DNA methylation signature was deemed to be a potential novel prognostic biomarker of use in increasing the accuracy of predicting overall survival of ccRCC patients. © Copyright: Yonsei University College of Medicine 2019.

INTRODUCTION

Clear cell renal cell carcinoma (ccRCC) is the most potentially lethal subtype of renal cell carcinoma (RCC) and accounts for 94% of metastatic RCCs.123 Early diagnosis of ccRCC is difficult, and by the time it is identified, it typically has already metastasized in 30–35% of patients.4 In early stage patients (stages I and II), most ccRCCs are treated by resection via surgery, with survival rates over 5 years greater than 70%, while as for patients in stages III (regional spread) and IV (metastatic), which have a poor prognosis, survival rates at 5 years are only 50% and 10%, respectively.567 Therefore, predictive methods with high specificity and sensitivity may improve the prognosis of patients.

Important targets in cancer research, cancer initiation and proliferation are regulated by epigenetic modifications.8 As one of the basic mechanism of epigenetic modification, DNA methylation plays as a pivotal role in the processes of chromatin structure regulation and gene regulation.9 Abnormal DNA methylation was reported to show a significant relationship with the occurrence and prognosis of various cancers.10 Dubrowinskaja, et al.11 demonstrated that neurofilament heavy (NEFH)-specific hypermethylation in RCC is relevant to overall survival rates among patients undergoing targeted therapy. Mir-874, which is downregulated in breast cancer, is correlated with the prognosis of breast cancer.12 Zhang, et al.13 demonstrated that decreased of mir-874 expression is mediated by DNA methylation of the promoter region and that DNA methylation is directly related to the prognosis of breast cancer. Furthermore, DNA methylation of specific genes has already been used as biomarkers in early diagnosis and prognosis,14,15 such as aldehyde dehydrogenase 1 family member A3 (ALDH1A3) promoter hypermethylation in glioma CpG island methylator phenotype-primary glioblastoma16 and hypomethylation of cyclin dependent kinase inhibitor 2A (CDKN2A) in colon cancer.17 Moreover, for patients with stage I non-small cell lung cancer, a shorter relapse-free survival rate has been found to be associated with higher methylation levels of CpG sites within five genes (HIST1H4F, PCDHGB6, NPBWR1, ALX1, and HOXA9).18

It was reported that ccRCC typically shows features of epigenetic alterations and methylated CpGs accompanying progress of disease,1920 although a specific prognostic model has not been proposed. This study analyzed whole-genome methylation profiles for ccRCC in The Cancer Genome Atlas (TCGA) database to identify potential biomarkers for predicting survival times for ccRCC patients.

MATERIALS AND METHODS

DNA methylation information for ccRCC tissues in the TCGA dataset

The DNA methylation information of ccRCC (platform: Illumina Human Methylation 27; downloaded in March 27, 2019; depending on GRCh38, the genomic coordinates of the CpGs were gained) and relevant clinical data were obtained from the TCGA database (https://cancergenome.nih.gov/). No further approval from the Ethics Committee was requested due to the data being obtained from the TCGA dataset. We studied how DNA methylation levels impacted survival in ccRCC using dates reported in survival records of the patients. β values were used to express DNA methylation levels, and M/(M+U) was computed. In M/(M+U), M is defined as the signal from methylated beads, while U is defined as the signal from unmethylated beads at targeted CpG sites. Ultimately, 215 samples comprised 96 ccRCCs of grade 2, 82 of grade 3, and 28 of grade 421 and 113 of stage I, 27 of stage II, 49 of stage III, and 25 of stage IV.22 In total, 23171 DNA methylation sites were selected in this work, and samples with relevant clinical records were obtained from the TCGA dataset. According to the TCGA series number, these 215 samples were assigned into a training dataset (60%) and a validation dataset (40%): one was used for discriminating and designing prognostic markers, while the other was applied to prove the precision of the biomarkers in survival prediction.

Prognostic identification and choosing CpG sites

Detection of methylation signatures was performed by means of the computing environment R with Survival package (R version 3.5.1; Revolution Analytics, Mountain View, CA, USA). First, methylation markers (p<0.05) found to be closely related to patient survival were defined as candidate markers through univariate Cox proportional hazard analysis. Then, robust likelihood-based survival analysis modeling and multivariable Cox regression analysis was applied to further screen for CpG sites.23

Establishment and validation of risk scores

We established risk scores for prognosis-related CpG sites by means of estimation of the regression coefficients in the multivariable Cox regression analysis of the training set. p values were characterized using the Wilcoxon rank-sum test, followed by false discovery rate correction. Then, we ranked the samples according to the risk score. We drew cross-validated time-dependent receiver operating characteristic (ROC) curves to verify the prediction accuracy of this model and to develop an effective optimal threshold. Based on this threshold, the samples were classified into high- and low-risk groups. Thereafter, cumulative survival time was calculated by a nonparametric statistic (Kaplan-Meier estimator), with a log-rank test (Mantel-Cox). Also, differences between the two groups in overall survival time were compared. Furthermore, the risk score system was evaluated by fitting in the checking set and the complete set. Also, by using different regrouping methods, we predicted the performance of the DNA methylation signature.

RESULTS

Collection of data from TCGA

The study was performed as indicated in Fig. 1. In total, 23171 DNA methylation sites were collected for 215 patient samples from the TCGA DNA methylation and clinical dataset. The median age of the patients was 60 years (range, 33–86 years), and the median survival time was 2386 days. The patients were randomly assigned to a training set (n=129, used to identify key CpG sites) and to a testing set (n=86, used to verify the methylation signature) based on status (alive/dead). The clinical data, including age, sex, ethnicity, stage, and status, are listed in the Table 1.

Fig. 1

Flow-chart of the study. The order of analyses to develop the risk score model and to validate the efficiency of the signature to predict prognostic outcomes. TCGA, The Cancer Genome Atlas; ccRCC, clear cell renal cell carcinoma.

Table 1

Clinicopathological Characteristics of ccRCC Patients from TCGA

Characteristics	Groups	Patients
		Total (n=215)	Training dataset (n=129)	Validation dataset (n=86)
Age (yr)	Median	60	59	60
	Range	33–86	33–83	33–86
	≤60	106	66 (62.26)	40 (37.74)
	>60	108	62 (57.41)	46 (42.59)
	Unknown	1	1 (100.00)	0 (0.00)
Stage	I	109	65 (59.63)	44 (40.37)
	II	30	16 (53.33)	14 (46.67)
	III	51	31 (60.78)	20 (39.22)
	IV	24	17 (70.83)	7 (29.17)
	Unknown	1	0 (0.00)	1 (100.00)
Histologic	G2	94	53 (56.38)	41 (43.62)
	G3	89	57 (64.04)	32 (35.96)
	G4	25	15 (60.00)	10 (40.00)
	Others	7	4 (57.14)	3 (42.86)
Gender	Male	145	83 (57.24)	62 (42.76)
	Female	70	46 (65.71)	24 (34.29)
Status	Alive	146	86 (58.90)	60 (41.10)
	Dead	69	43 (62.32)	26 (37.68)

ccRCC, clear cell renal cell carcinoma; TCGA, The Cancer Genome Atlas.

Values are presented as n (%) unless otherwise indicated.

Screening and identification of DNA methylation markers in the training dataset

Methylation levels were used as variables in the training dataset to identify DNA methylation markers related with ccRCC, and based thereon, univariate Cox proportional hazard regression analysis and robust likelihood-based survival analysis were performed. The first 19 significant prognosis-related CpG sites were chosen (p<0.001, Supplementary Table 1, only online). Next, the 19 CpG sites of prognostic significance were put into multivariable Cox regression survival analysis to further screen CpG sites in order to identify an optimum prognostic model for predicting ccRCC. Finally, two methylation sites (cg18034859, cg24199834) in the hazard ratio (HR) model were identified; their risk coefficients are listed in Table 2. cg24199834 was indicated as a risk factor with a HR of 1.03, while the methylation level of cg24199834 was inversely proportional to the risk of death, with an HR of 0.78, indicating that cg18034859 was a protective factor. Accordingly, we developed the following risk score formula: Risk score=-0.2455×β value of cg18034859+0.0304×β value of cg24199834. The formula demonstrated that the hypermethylation levels of cg24199834 are associated with higher risk, while the hypomethylation levels of cg18034859 are related to lower risk.

Table 2

Multivariable Survival Analysis

CpG sites	Gene	Position	Coef*	Exp (coef)*	p value†	95% CI*
cg18034859	MYLK2	chr20:29870349	-0.2455	0.7823	<0.001	0.6992–0.8752
cg24199834	POU4F2	chr4:147779576	0.0304	1.0309	<0.001	1.0159–1.0461

CI, confidence interval.

*In multivariable Cox regression analysis; †p value was calculated through the Wilcoxon rank-sum test followed by false discovery rate adjustment for multiple corrections.

Evaluation of the predictive performance of the two-DNA methylation signature

After screening, 0.54 was identified as the optimal threshold cut-off point, exhibiting the best performance (Fig. 2A). The area under the curve of ROCs (AUC) for the two-DNA methylation signature was 0.819, indicating a favorable prognostic value in prediction of patient survival. Then, the training dataset was divided into high-risk (n=87) and low-risk groups (n=42). Kaplan-Meier survival analysis was applied to compare the survival time of the two groups. The results indicated that patients in the high-risk group had shorter survival time than the low-risk group. A significantly worse prognosis was embodied in the high-risk group (p<0.001, Fig. 2B). To further validate the prognostic ability of the two-DNA methylation signature, we divided the validation dataset into high-risk (n=58) and low-risk groups (n=28). Similar to the above, the high-risk group had significantly shorter survival time in Kaplan-Meier survival analysis (Fig. 2B). Fig. 3A shows the risk score distribution and the scatter plot of patient states, with higher risk corresponding to more deaths. Heat mapping depicts the methylation levels of cg18034859 and cg24199834 between the two groups for the entire dataset. Significantly higher methylation levels for cg24199834 were confirmed, while the opposite was true for cg18034859 (Fig. 3B).

Fig. 2

Evaluation of the predictive performance of the two-DNA methylation signature. (A) Receiver operating characteristic analysis of the sensitivity and specificity for survival time by the two-DNA methylation signature in the training dataset. The black dot represents the optimal validation dataset using the two-DNA methylation signature. (B) The Kaplan-Meier analysis was used to visualize the survival probability for the low-risk versus high-risk groups of patients based on the optimal cut-off point. Rows represent survival time (days), and columns represent survival rate. TPR, true positive rate; FPR, false positive rate; AUC, area under the curve of receiver operating characteristics.

Fig. 3

Risk score analysis of the training set. The distribution of two-DNA methylation based risk core, patient survival, and methylation levels of two CpGs were analyzed in the training set (n=129). (A) Two-DNA methylation signature risk score distribution. ‘Time’ means ‘survival time’. (B) Heat-map of the DNA methylation profiles. Rows represent CpG sites, and columns represent patients.

Predictive performance of the two-DNA methylation signature

Studies have reported that age,24 sex,25 pathological stage,26 and histologic grade etc.27 are associated with prognostic survival. Accordingly, all patient samples were regrouped according to clinicopathological characteristics, and the two-DNA methylation signature was validated separately in the different groups. First, the effect of the two-DNA methylation signature among the different age groups of patients (based on survival status, “alive” used the final follow-up time; otherwise the time of death): ≤60 (n=106, 49.53%) and >60 (n=108, 50.47%) years were analyzed (Fig. 4A). From Kaplan-Meier curves, we found that the overall survival of patients in different risk groups differed between the two age cohorts (p<0.01), and AUC values were 0.860 in patients younger than 60 years old and 0.666 in patients older than 60 years old (Fig. 4B), exhibiting the independent relationship between the two-DNA methylation signature and age. Ricketts, et al.25 found that the spread of ccRCC among the sexes is distinctly uneven. Sex may affect mutation spectra and hold implications for prognosis. A summary analysis of age, sex, pathological stage, and histologic grade is listed in Table 3. The two risk groups both in male and female patients (p<0.01) reflected obviously distinct survival times, and the AUCs for male and female sex were 0.772 and 0.832, respectively. Taking the sample number into account, we verified the predictive value of the two-DNA methylation signature in G2 (n=96) and G3 (n=82), with correlated AUC values of 0.877 and 0.628. About one quarter of the patients showed distant metastasis at initial diagnosis, and our results showed that overall survival was significantly different between the risk groups in pathological stage IV patients (p=0.01). The AUCs in stages I, II, III, and IV cohorts were 0.714, 0.686, 0.721, and 0.577, respectively. From the above results, the two-DNA methylation signature was proven to be an independent applicable value for predicting patient survival.

Fig. 4

Stratified analysis. Kaplan-Meier and receiver operating characteristic (ROC) analyses of patients with clear cell renal cell carcinoma of different age. Grouping was based on age at initial diagnosis: ≤60 (n=106, 49.53%), >60 (n=108, 50.47%). (A) Kaplan-Meier analysis with two-sided log-rank test was performed to estimate the differences in survival time between the low-risk and high-risk patients. (B) ROC curves of the two-DNA methylation signature were used to demonstrate the sensitivity and specificity in predicting the survival time of patients. TPR, true positive rate; FPR, false positive rate; AUC, area under the curve of receiver operating characteristics.

Table 3

Results of Kaplan-Meier and ROC Analyses Based on the Different Regrouping Methods

Regrouping factors	Group	Sample size	p value	AUC
Age (yr)	≤60	106	0.001	0.860
	>60	108	0.003	0.666
Stage	I	113	0.006	0.714
	II	27	0.290	0.686
	III	49	0.100	0.721
	IV	25	0.012	0.577
Histologic grade*	G2	96	<0.001	0.877
	G3	82	0.145	0.628
Gender	Female	76	<0.001	0.832
	Male	139	0.005	0.772

ROC, receiver operating characteristic; AUC, area under the curve of ROCs.

*G1 and GX (the degree of differentiation is unknown) were excluded in this study as these tumors may exhibit different biological behaviors. G4 was also excluded due to too few samples.

DISCUSSION

A lack of effective and reliable prognostic biomarkers or models remain the focus of ccRCC research. Studies have demonstrated that molecular markers, including DNA methylation, can be used to predict the prognosis of ccRCC.28293031 Compared to using individual DNA methylation factors as predictors, combinations of DNA methylation could achieve higher sensitivity and specificity.32 In the present study, we generated a risk score model based on a two-DNA methylation signature significantly related to ccRCC prognosis in Cox regression and ROC analysis. Median survival time was used to divide high- and low-risks groups, which has more clinical practicability than median risk score. Calculating patient risk scores could be helpful to predict the probability of patients reaching median survival time. Time-dependent ROC curves showed that our two-DNA methylation signature has stable and good performance in prognosis prediction. The model could distinguish different risk groups with totally varied survival time, which was verified in training and validation datasets. Independent grouping tests showed that, the two-DNA methylation signature possessed independent predictive value from age and sex, and was more suitable to forecast survival time for patients in stage I and G2.

Compared to other acquainted prognostic markers, it is found that the two-DNA methylation signature possesses noticeable higher sensitivity and specificity in predicting ccRCC prognosis. Zhan, et al.28 reported a five-gene prognostic model to predict prognosis in ccRCC, with an AUC of 0.783. The sixlncRNA signature identified in a previous study33 was described as a prognostic biomarker of ccRCC, with an AUC of 0.696. In this study, low- and high-risk groups could be well distinguished by way of the two-DNA methylation signature. Also, the two-DNA methylation signature performed well in associated log-rank tests, confirming its independent predictive value for patient survival.

Exploring another database (http://genome.ucsc.edu/), we found cg18034859 and cg2419834 are located in the promoter regions of myosin light chain kinase 2 (MYLK2) and POU class 4 homeobox 2 (POU4F2) genes, respectively, indicating the potential influence of its methylation on gene expression. POU4F2 is a member of the POU-domain transcription factor family, which is involved in determining the fate of retinal ganglion cells (RGCs) and initiating RGC differentiation.34 The methylation levels of POU4F2 were quantified by quantitative methylation-specific PCR using urine sediment DNA from another study, and the researchers generated a two-DNA methylation signature [combination of POU4F2 and protocadherin 17 (PCDH17)] for bladder cancer. The model showed high sensitivity and specificity of 90.00% and 93.96%, respectively, reflecting the ability of detecting bladder cancer effectively among pathologically varied sample groups.35 Moreover, urinary levels of POU4F2 hypermethylation were considered to be having a bearing on non-muscle invasive bladder cancer recurrence, which is in good association with our results for a high DNA methylation level of cg2419834 at POU4F2.36 In addition, gene expression levels of POU4F2 in ccRCC were not detectable based on Gene Expression Profiling Interactive Analysis (GEPIA, http://gepia.cancer-pku.cn/index.html) and TCGA database (data not shown). Meanwhile, MYLK2 encodes a calcium/calmodulin-dependent serine/threonine kinase and is a risk factor for the development of acute lung injuries.37 Soung, et al.38 reported that MYLK2 does not play a significant role in cancer pathogenesis because the kinase domain of MYLK2 is rarely mutated in common human carcinomas. The gene expression data of MYLK2 were downloaded from TCGA, while the clinical data of ccRCC were the same. We merged and ranked samples according to the expression levels of MYLK2 and divided them into high-expression (n=89) and low-expression groups (n=89). Kaplan-Meier survival analysis methods exhibited a difference in survival time between high-expression and low-expression groups: the survival time of patients in the high-expression group was significantly shorter (p=0.028) (Supplementary Fig. 1, only online). The results indicated that MYLK2 may be an independent prognostic marker of survival for patients with ccRCC.

In conclusion, we discovered that the survival time of patients is greatly affected by a two-DNA methylation signature, which was detected from genome-whole analysis of DNA methylation data for 215 samples with ccRCC, and the significance for patients of varying ages, sex, stages, and histologic grades was confirmed. The two-DNA methylation signature could play an important role as a prognostic marker for forecasting the survival time of patients with ccRCC. In the future, research may want to focus on a functional mechanism for the two-DNA methylation marker and its potential association in carcinogenesis.