Introduction

Colorectal cancer is the third leading type of cancer and second leading cause of cancer-related deaths worldwide. According to the World Health Organization (WHO), the worldwide age-adjusted incidence rate of colorectal cancer was 19.5 per 100,000, and the age-adjusted mortality rate of colorectal cancer was 9.0 per 100,000 in 2020.¹ Although a decline in the incidence and mortality of colorectal cancer has been observed in the United States since the mid-1980s, it is still the fourth leading type of cancer.² With an annual age-adjusted incidence of 36.3 per 100,000 and age-adjusted mortality of 12.8 per 100,000, colorectal cancer accounts for 8% of all new cancer cases and 8.7% of all cancer-related deaths in the United States.³ The 5-year relative survival rate for colorectal cancer patients in the United States is 63.7%.⁴ Despite the overall decline in the colorectal cancer incidence and mortality, the colorectal cancer incidence among people younger than 50 years has increased every year by 2% from 2012 to 2016, and mortality among patients younger than 55 years has increased by 1% per year from 2008 to2017.³

Various modifiable and non-modifiable risk factors increase the risk of colorectal cancer development. Modifiable risk factors include behavioral factors such as overweight or obesity, physical inactivity, certain types of diets, smoking, and heavy alcohol use.^5,6 On the other hand, non-modifiable risk factors include, older age, past history or family history of colorectal polyps or colorectal cancer, history of inflammatory bowel disease, some inherited syndromes, and type 2 diabetes.⁶ In addition to these risk factors, recent studies have identified the roles of several genetic factors in the etiology and prognosis of colorectal cancer.

Genetic expression is an essential factor for studying colorectal cancer.⁷ Previous studies have identified more than 300 genes associated with colorectal cancer.^8–11 In addition to acquired gene mutations, certain types of colorectal cancer are associated with inherited gene mutations. Familial adenomatous polyposis (FAP) and hereditary nonpolyposis colorectal cancer (HNPCC), also known as Lynch syndrome, are the two most common inherited conditions that increase the risk of colon cancer. FAP syndrome is associated with inherited gene mutation in APC gene.¹² Lynch syndrome is caused by inherited gene mutations in MSH2, MSH6, PMS2, and MLH1 genes.¹¹ Peutz-Jeghers syndrome, another inherited condition caused by a mutation in the STK11 gene, is also associated with colorectal cancer.¹¹ FAP patients inevitably develop colorectal cancer, as 95% of them develop adenomas by 35 years of age.¹³ Similarly, individuals with Lynch syndrome have cumulative lifetime risk of 78% for colorectal cancer.¹⁴

The genomic profiling of cancers has a significant impact on cancer treatment, cancer risk prediction, and risk reduction. Recent studies have used gene expression to predict cancer survival and have developed several survival prediction models to estimate the survival of patients with different types of cancer.¹⁵ The multi-omics approach has contributed to the identification of various cancer-specific mutations, epigenetic alterations, and cellular responses to chemotherapy.¹⁶ It provides a comprehensive understanding of the disease by studying the molecular intricacies at multiple levels. These levels –omics – include the genome, transcriptome, metabolome, proteome, exome, and the microbiome.¹⁷ Using a multi-omics approach aids in determining genomic alterations that are specific to different types of cancers.¹⁸ Additionally, there is evidence that microRNAs (miRNAs) play vital roles in the carcinogenesis and development of colorectal cancer.¹⁹ Therefore, miRNA biomarkers play an important role in early diagnosis, treatment, recurrence, and survival prediction of colorectal cancer.²⁰

Previous research has examined the overall burden of colorectal cancer by estimating its prevalence, incidence, and mortality in different populations. Recent studies have identified variations at the molecular level, including genetic factors, associated with colorectal cancer.²¹ However, there is a paucity of studies using gene expression data to evaluate the prognosis of colorectal cancer. Furthermore, survival prediction models based on multiple factors, including gene expression, have not been studied in patients with colorectal cancer. Polygenic risk score analysis is an important approach to assess the risk of cancer at any time in life and to easily identify individuals who are at high risk of cancer-related morbidity and mortality.^22,23 This study aimed to examine the effects of various factors on colorectal cancer survival, as well as, determine different factors including genetic variations that predict colorectal cancer survival and to develop a polygenic survival scoring system for colorectal cancer. The Cancer Genome Atlas (TCGA) data were analyzed to identify the genetic, clinical, and demographic predictors of colorectal cancer survival. Based on the identified predictors, a survival prediction model and polygenic survival score were developed and validated.

Material and Methods

Data

Data for this study were obtained from the Cancer Genome Atlas data portal (TCGA), known as the Web Portal for Genomic Data Commons (GDC), and the National Cancer Institute Genomic Data Commons.²⁴ TCGA provides a comprehensive profiling of more than 30 cancers, including clinical and genomic data. Two datasets were obtained from TCGA: The Cancer Genome Atlas Colon Adenocarcinoma (TCGA-COAD) and the Cancer Genome Atlas Rectum Adenocarcinoma (TCGA-READ). Each TCGA-COAD and TCGA-READ dataset contained approximately 500 specimens.^25–27 In this study, we used clinical data and genomic measurements, specifically mRNA expression.

Study Design and Sample Size

This cohort study included 459 patients from the TCGA-COAD dataset, which was used to identify the predictors of overall survival among colorectal cancer patients, develop a survival prediction model, and develop a polygenic survival scoring system. A total of 170 patients from the TCGA-READ dataset were used to validate the survival prediction model and the polygenic survival scoring system.

Study Variables

Overall Survival

The primary outcome of interest was the overall survival (OS) of patients with colorectal cancer. This was defined as the time interval in years from the date of cancer diagnosis to the date of death. Accordingly, patients who did not experience the event (ie, death) were censored on the date of their last known contact.

Demographic Variables

Information related to the patient age, sex, and race was also included. Age was used as a categorical variable and dichotomized into less than 65 and ≥ 65 years. Due to the small sample size of some racial groups, we categorized race into three groups: White, Black, and Others.

Clinical Variables

Clinical characteristics, such as pathological tumor stage and histological type, were included. The stages were classified into four groups: Stage I contained I and IA; Stage II contained IIA, IIB, and IIC; Stage III contained III, IIIA, IIB, and IIIC; and Stage IV, IV, IVA, and IVB.

Gene Expressions

A total of 205 genes previously associated with colorectal cancer were included in the analysis.^8–11,28,29 The Database for Annotation, Visualization, and Integrated Discovery (DAVID) bioinformatics gene ID conversion tool was used to convert official gene symbols to Ensembl gene ID.³⁰ The mRNA gene expression data were standardized with a mean of 0 and standard deviation of 1 to achieve a normal distribution.

Statistical Analysis

Descriptive statistics were calculated for all the study variables. To examine the differences between TCGA-COAD and TCGA-READ datasets, the chi-square test was used for categorical variables, whereas the Log rank test was used for follow-up time. Additionally, the median overall survival times for TCGA-COAD and TCGA-READ cohorts were calculated using the Kaplan-Meier method.

Identification of Survival Predictors and Survival Prediction Model

Univariate analyses were performed to assess the effect of each study variable, including age, sex, race, cancer stage, and 205 gene expressions on OS using Cox proportional hazards (CoxPH) models. The hazard ratio (HR) from the Cox proportional hazards model was used to measure effect size. Variables that had a significant univariate association with OS were considered for multivariate Cox proportional hazards regression analysis. For the survival prediction model, gene expression was converted to z-scores using the sample mean and standard deviation from the TCGA-COAD data. Multivariate Cox proportional hazards regression models were used to develop the survival prediction model. As a first step, variables included in the survival prediction model were determined based on stepwise selection with a significance level of 0.99 to enter the model and 0.995 to stay in the model, which resulted in a full model consisting of all potential predictors.³¹ Subsequently, the corrected Akaike Information Criterion (AIC_c) was used to identify a set of (k) top AIC-optimized models with the lowest 3 ± models to obtain more candidates and to increase robustness when including more predictors around the smallest AIC_c. The variables in the (k) set of models were included in the final model. The receiver’s Harrell’s C-index, also known as the concordance index, was used to evaluate the performance of the survival prediction model.³²

Calculation of Polygenic Survival Score

The polygenic survival score calculation was based on the scoring procedure proposed by Fu et al to predict ovarian cancer survival.³¹ The regression coefficients obtained from the survival prediction model (multivariable Cox proportional hazards model) were used to calculate the polygenic survival score (PSS). The following formula summarizes the PSS calculation:

where PSS_j, is the PSS for the jth patient for survival from colorectal cancer, zx_i is the standardized z-score of the ith gene expression predictor, w_i is the ith categorical (demographic/clinical) predictor, and β_i is the corresponding regression coefficient obtained from the multivariate Cox proportional hazards model.³¹ The PSS is a continuous measure of colorectal cancer survival, with a lower score predicting longer survival.³¹ The overall accuracy of the PSS in predicting survival was examined using a Receiver Operating Characteristic (ROC) curve. The area under the curve (AUC) and 95% confidence intervals (CI) were computed to evaluate the performance of PSS for one, three, and five-year survival among patients with colon cancer. Furthermore, the median PSS risk score was used as the threshold to classify colon cancer patients into high risk group (PSS greater than the median) and low risk group (PSS less than the median). Kaplan-Meier plot was generated for the two PSS based on risk groups. Median survival was also obtained for each risk group.

Validation of the Survival Prediction Model and Polygenic Survival Score

The TCGA-READ dataset was used to validate the survival prediction model, and PSS was developed based on TCGA-COAD. Predictors of survival among patients with colon cancer, age, stage, and the ten genes identified based on TCGA-COAD analysis were used to examine survival among patients with rectal cancer. The predictor variables were standardized by converting them to z-scores using mean and standard deviation. A multivariate Cox proportional hazards model was used to examine its association with overall survival. The performance of this survival prediction model was further assessed using the Harrell’s concordance index (C_H-index). The accuracy of PSS in predicting survival among rectum cancer was evaluated by ROC analysis and AUC. The Kaplan-Meier method was used to estimate the median survival time for the PSS based on the high and low risk groups.

Data were analyzed using SAS v9.4 and R RStudio v1.12.1335 with packages including: SummarizedExperiment, DESeq2, and TCGAbiolinks. Statistical significance was set at an alpha level of 0.05.

Results

The majority of patients in both datasets were 65 years or older (62% of TCGA-COAD patients and 54% of TCGA-READ patients) and males. There was no significant difference between TCGA-COAD and TCGA-READ cohorts for age (p= 0.052), gender (p= 0.80), cancer stage (p = 0.3284), and histological type of cancer (p=0.0541). Similarly, the median overall survival among TCGA-COAD (7 years) and TCGA-READ (5 years) patients was not significantly different (p= 0.211) (Table 1 and Figure 1).

Table 1 Descriptive Statistic of TCGA-COAD and TCGA-READ Patients Characteristics

Figure 1 Kaplan–Meier Survival Plot For TCGA-COAD and TCGA-READ.

Identification of Survival Predictors and Survival Prediction Model

The results of Cox proportional hazards regression analysis showed that age, cancer stage, and expression of 27 genes had a univariate association with colon cancer survival. Further analysis of these variables with a multivariable Cox proportional hazards regression model using stepwise and AIC_c selection methods retained age, cancer stage, KLF5, ETS2, SMAD2, TRAP1, HIGD1A, AXIN2, MMP1, TOP2A, and HSPE1 gene expression as predictors of colon cancer survival. (Table 2) The risk of death among older cancer patients was more than twice that among those aged < 65 years (HR: 2.39, 95% CI: 1.46, 3.93) when adjusted for cancer stage and 9 gene expressions. Similarly, increasing cancer stage had a significantly worsening impact on survival when adjusted for other covariates (Table 2). In terms of gene expression, except for HSPD1, the expression of all other genes reduced the risk of death in patients with colon cancer. TRAP1 expression was associated with the highest reduction in the risk of death by 37.4% (HR: 0.63, 95% CI: 0.46, 0.85), whereas TOP2A expression had the lowest reduction in the risk of death (HR: 0.86, 95% CI: 0.66, 1.12) when adjusted for other covariates.

Table 2 Multivariate Cox Proportional Hazards Model for Survival Predictors of Colon Cancer (TCGA-COAD)

Performance Evaluation of Survival Prediction Model

Predictive accuracy of the survival prediction model was assessed with the ROC analysis. The model had high predictive accuracy based on the AUC results. For year 1 the overall accuracy of the model in predicting survival was 0.75 (95% CI: 0.65, 0.84); for year 3 it was 0.81 (95% CI: 0.74–0.88), and for year 5, the AUC was 0.78 (95% CI: 0.68–0.88). Figure 2 illustrates ROC curves for the model at year 1, 3, and 5. The overall predictive performance of the final Cox proportional hazards model (C_H-index) was 0.76 (95% CI: 0.70, 0.82).

Figure 2 Time-dependent-Roc Curves Estimated for TCGA-COAD Patients. (A) ROC at year 1. (B) ROC at year 3. (C) ROC at year 5.

Polygenic Survival Score

Parameter estimates (β) for the variables and the values of the corresponding variables were used to calculate the polygenetic survival score (PSS). The mean of the polygenic survival score for the colon cancer patients was 0.44 (SD= 1.13), and the median was 0.38 (min= −2.93, max= 3.84). The median survival score was used as the cutoff value to categorize patients into high and low risk groups. Those who had a polygenic survival score greater than 0.38 were categorized as the high risk group. The results of Figure 3 showed that patients in the low risk group had a significantly longer median survival time than those in the high risk group (p <0.0001). Three-year survival probability for the low risk group was 91.62%, whereas that for the high risk group was 62.71%. Similarly, the 5-year survival rate was higher in the low risk group than that in the high risk group (82.36% vs 45.16%).

Figure 3 Kaplan–Meier Plot for overall survival by polygenic survival score categories among colon cancer patients (TCGA-COAD).

External Validation

Survival Prediction Model

The predictive performance of the survival model was assessed using the data obtained from patients with rectal cancer (TCGA-READ). The results of the ROC analysis showed that, as time increased from year 1 to year 5, the AUC of the model increased. Thus, the survival prediction model became more accurate as duration increased. The accuracy of the model in predicting survival was highest for year 5 (AUC: 0.82, 95% CI: 0.63–1.00), followed by year 3 (AUC: 0.80, 95% CI: 0.67, 0.93), and year 1 (AUC: 0.75, 95% CI: 0.58, 0.92). ROC curves (Figure 4) showed high sensitivity and specificity of the model in predicting survival at different years. The overall adequacy of the prediction model assessed by C_H-index was 0.80 (95% CI: 0.71, 0.88).

Figure 4 Time-dependent-Roc Curves Estimated for TCGA-READ Patients. (A) ROC at year 1. (B) ROC at year 3. (C) ROC at year 5.

Polygenic Survival Score

Parameter estimates obtained from TCGA-COAD were multiplied by the corresponding variables from TCGA-READ and the values were summed to obtain the PSS for each patient. The mean of the PSS for TCGA-READ was −0.24 (SD= 1.13) and the median was −0.40 (range: −2.64 to 3.89). The median risk score from TCGA-COAD (0.38) was used as the cutoff value to categorize patients in TCGA-READ into high and low risk groups. The median survival time was significantly higher in the low risk group than in the high risk group (p= 0.0012) (Figure 5). Three-year survival probability for the low risk group was 90%, while that for the high risk group was 84%. Similarly, the 5-year survival rate was higher in the low risk group than that in the high risk group (74% vs 11%).

Figure 5 Kaplan–Meier Plot for overall survival by polygenic survival score categories among rectum cancer patients (TCGA-READ).

Discussion

The findings from this study indicate that age and cancer stage are associated with survival in patients with colorectal cancer. The risk of death was twice as high among older patients (≥ 65 years) than that among younger patients. Recent research based on TCGA found that colorectal cancer patients aged 65 years or older had 2.61 times higher risk of death than those aged less than 65 years (95% CI: 1.59–4.28).¹⁰ Findings from the NY State Cancer Registry and Statewide Planning and Research Cooperative System cohort concluded that older patients had a higher risk of colon cancer-specific deaths (65–74 year: HR= 1.56, 95% CI: 1.26–2.00; ≥ 75 years: HR= 2.57, 95% CI: 2.29–6.05).³³

The findings of the current study suggest that there is a gradient in the association between cancer stage and risk of death from colon cancer. Compared to stage I, the hazard of death was highest among stage IV patients, followed by stage III and stage II patients. Although stage I patients had significantly higher survival rates than stage III and IV patients, no significant difference was observed in stage II patients. Our results regarding Stage II survival contradict the findings of previous research studies.^34–37 Other studies have shown that stage I colorectal cancer patients have a significantly better prognosis in terms of survival compared to all advanced stages, including stage II. A retrospective study of colorectal cancer patients found that stage II, III, and IV patients had significantly lower colorectal cancer survival than stage I patients: stage II (HR= 1.40, 95% CI: 1.09–1.80), stage III (HR= 2.32, 95% CI: 1.76–3.04), and stage IV (HR= 6.462, 95% CI: 4.50–9.24).³⁶ Similarly, another retrospective study concluded that all advanced stages were associated with significantly lower survival rates than stage I.³⁷ The differences in our findings may be attributed to the differences in patient populations between these studies.

Previous studies identified more than 200 genes as risk factors for colorectal cancer. The current study used previously identified genes to examine their roles in the survival of patients with colorectal cancer. The results showed that nine genes, KLF5, ETS2, SMAD2, TRAP1, HIGD1A, AXIN2, MMP1, TOP2A, and HSPE1, increased survival, whereas one gene, HSPD1, decreased survival among patients with colorectal cancer. Previous studies have shown that increased KLF5 expression is linked to poor prognosis in patients with colorectal cancer by promoting cancer stem cell-like properties and cell proliferation.³⁸ A study examined the role of ETS2 in colorectal cancer and found that the expression of ETS2 in human colon cancer cell lines is linked to increased cell survival and proliferation.³⁹

Similarly, our findings of a significant association between SMAD2, TRAP1, HIGD1A, and AXIN2 gene expression and increased survival among colorectal cancer patients are in accordance with the findings of previous studies of colorectal cancer. A previous study showed that there was a considerable difference in the overall survival of colorectal patients who demonstrated low or high expression of SMAD2/3 in normal tissues.⁴⁰ Colorectal cancer patients with increased TRAP1 expression have aggressive lymph node metastasis, multidrug tolerance, and shorter overall median survival.^41,42 Another study reported that HIGD1A inhibited tumor growth but facilitated tumor cell survival in vivo in different types of cancer including colorectal cancer.⁴³ There are inconsistencies in the findings of previous studies examining the association between AXIN2 and colorectal cancer.⁴⁴ A few studies have shown a protective effect of AXIN2 in colorectal cancer, whereas other studies have found a poor prognosis among colorectal patients with mutations in AXIN2 gene.^44–46

Although this study found a significant role of MMP-1 in colorectal cancer survival, there are inconsistent findings regarding the effect of MMP-1 on colorectal cancer tissue from in-vitro studies.^47,48 A previous study determined that MMP-1 transcription rates are associated with poorer results in colorectal cancer.^47,48 However, another study indicated that high MMP-1 expression is associated with longer overall cancer tissue survival.⁴⁸ The function of HSPD1 in colorectal cancers remains unknown, particularly in adenocarcinoma.⁴⁹ A cohort study showed that high levels of HSPD1 protein are associated with tumorigenesis and poor prognosis in patients with colorectal cancer, suggesting that HSPD1 causes tumor metastasis.⁴⁹ Recent studies have shown that TOP2A plays a major role in tumorigenesis in many malignancies, such as prostate, breast, and ovarian cancers. However, there is limited evidence regarding the role of TOP2A in colorectal cancer.⁵⁰ HSPs are generally overexpressed in various tumor types. Increased expression of HSPs is usually associated with a poor prognosis and increased resistance to therapies.⁵¹ However, the association between HSPE1, TOP2A, and colorectal cancer survival rates has not yet been established.

The survival prediction model developed in this study had adequate accuracy, with an AUC ranged from 0.74 to 0.81. Although the accuracy of the model increased from year 1 to year 3 and decreased from year 3 to year 5, the difference between years was not statistically significant.

The overall accuracy of our model was adequate, as suggested by previous studies.^52,53 We further validated the survival prediction model by using data from patients with rectal cancer. Both the colon and rectum form parts of the large intestine with no anatomical or physiological differences. Additionally, colon and rectal cancers are extremely similar at cellular and molecular levels. As colon and rectal cancers share significant similarities, data from patients with rectal cancer were used for external validation of the survival prediction model. The results of the external validity analysis showed a similar accuracy of the survival prediction model among rectal cancer patients to that observed for colon cancer patients.

Conclusion

The polygenic survival score (PSS) is an important measure for estimating survival based on risk factors, including changes in gene expression that impact colorectal cancer. To obtain polygenic survival scores, standardized z-scores of gene expression and parameter estimates from the survival prediction model were used. This approach ensured that the risk associated with each variable was measured on a uniform scale for all the patients. Thus, the estimation process improved the accuracy of the scoring system. In addition to evaluating PSS as a continuous measure of survival, we evaluated its accuracy by classifying patients into PSS based on low and high risk groups. We used the median PSS score as the cutoff value for this classification. The median was selected because it is a central parameter that is unaffected by extreme values. Therefore, it provides a better evaluation of PSS-based categories and their utility in clinical and research settings. Findings showed that the PSS based low risk group had a significantly higher median survival than the high risk group.

This study has several strengths. A comprehensive approach was used to include gene expression and a thorough literature search was conducted to include all possible known genes associated with colorectal cancer. The TCGA dataset, which provides information related to a wide range of genetic data, was used for this study. Therefore, the effects of more than 200 genes on colorectal cancer survival were examined. Potential confounding and effect modifications of the covariates were assessed. Another strength of this study was the validation of the survival prediction model and PSS using an external dataset. Despite its strengths, this study has several limitations. This study did not include the colorectal cancer specific mortality rates. Owing to the secondary analysis of the existing data, information about some important risk factors, such as behavioral factors, clinical characteristics, and environmental factors, was not available. Hence, the role of these important factors as potential confounders or effect modifiers has not yet been assessed. Similarly, TCGA data did not include geographic location, including the country of residence of the patients; thus, the generalizability of our findings could not be determined. The majority of patients were Caucasian, which might have introduced a selection bias and limited the generalizability of the study findings. This study did not use the staging system recommended by the American Joint Commission on Cancer. Due to the small sample size for some cancer stages, we grouped them into four broader stages. This approach addresses the issue of low statistical power at the cost of a more in-depth understanding of the effects of each cancer stage on survival.

This study is the first step towards developing a comprehensive risk scoring system that includes a wide range of genetic and other risk factors associated with survival among patients with colorectal cancer. The survival prediction model and PSS developed in this study can serve as important tools in colorectal cancer research. Furthermore, the PSS has clinical relevance as a scoring system for patients. The findings of this study will direct future research to utilize a multi-omics approach to study colorectal cancer survival.