Abstract
Colorectal cancer is the third leading cause of cancer-related death in the United States. About 15% of colorectal cancers are associated with microsatellite instability (MSI) due to loss of function in the DNA mismatch repair pathway. This subgroup of patients has better survival rates and is more sensitive to immunotherapy. However, it remains unclear whether microsatellite stable (MSS) patients with colorectal cancer can be further stratified into subgroups with differential clinical characteristics. In this study, we analyzed The Cancer Genome Atlas data and found that Chr20q amplification is the most frequent copy number alteration that occurs specifically in colon (46%) and rectum (61%) cancer and is mutually exclusive with MSI. Importantly, MSS patients with Chr20q amplification (MSS-A) were associated with better recurrence-free survival compared with MSS patients without Chr20q amplification (MSS-N; P = 0.03). MSS-A tumors were associated with high level of chromosome instability and low immune infiltrations. In addition, MSS-A and MSS-N tumors were associated with somatic mutations in different driver genes, with high frequencies of mutated TP53 in MSS-A and mutated KRAS and BRAF in MSS-N. Our results suggest that MSS-A and MSS-N represent two subtypes of MSS colorectal cancer, and such stratification may be used to improve therapeutic treatment in an individualized manner.
This study shows that chromosome 20q amplification occurs predominately in microsatellite-stable colorectal cancer and defines a distinct subtype with good prognosis, high chromosomal instability, distinct mutation profiles, and low immune infiltrations.
Introduction
Colorectal cancer is the third leading cause of cancer-related death in the United States (1). Loss of function in DNA mismatch repair pathways is a common molecular aberration in colorectal cancers. Dysfunction of genes in this pathway reduces the ability of cells to repair DNA replication errors and thereby leads to a microsatellite instable (MSI) subtype of colorectal cancer (2). In all colorectal cancer cases, around 15% patients are MSI, which have higher somatic mutation burden and immune infiltration in the tumor microenvironment compared with their microsatellite stable (MSS) counterparts (3). Immune checkpoint genes such as CTLA4 and CD274 are more highly expressed in MSI than MSS patients (4, 5). Indeed, a number of antibody-based drugs targeting these immune checkpoint proteins such as pembrolizumab and nivolumab have been approved by FDA for treating patients with metastatic MSI colorectal cancer (6). Apart from high sensitivity to immunotherapy, MSI status itself is a good prognostic marker for patients with colorectal cancer subject to conventional treatment (7). MSI patients exhibit less clinical aggressiveness and longer survival time than MSS patients (8, 9).
The majority of colorectal cancers are MSS cases, for which the driver genomic alterations are unclear. It has been reported that some MSS colorectal cancers are associated with high level of chromosome instability (CIN; ref. 10). CIN is a common molecular aberration observed in various cancer types and plays critical roles during tumorigenesis (11). Tumors with high CIN are more likely to be from late stage and associated with poor prognosis (12). In colorectal cancer, CIN tumors are mostly MSS and have high copy number variation (CNV) burden (10). It has been reported that high CIN is associated with significantly shorter recurrence-free survival and overall survival (OS) rates in stage 2 to 3 colorectal cancers (13). High CIN is typically associated with high number of focal and arm-level CNVs, which in turn is correlated with markers of immune evasion and reduced response to immunotherapy (14).
In this study, we hypothesize that a significant fraction of MSS colorectal cancer cases is driven by certain CNV events. To test this hypothesis, we performed systematic analysis to identify arm-level CNVs that occurred with high frequency in colorectal cancers but not or with low frequency in other cancer types using The Cancer Genome Atlas (TCGA) data. We discovered that chromosome 20q (Chr20q) amplification predominantly occurred in MSS colorectal cancers, enabling us to stratify them into two groups: colorectal cancers with Chr20q amplification (MSS-A) or without (MSS-N). A series of comparisons indicated that MSS-A tumors was associated with significantly better prognosis, higher CIN level and lower immune infiltration compared with the MSS-N tumors. Moreover, Chr20q amplification status was also correlated with the sensitivity of colorectal cancer cell lines to a number of drugs. These results suggest that MSS-A and MSS-N represent two colorectal cancer subtypes with potential clinical values.
Materials and Methods
Datasets
Level 3 processed TCGA RNA-sequencing data for COAD and READ were downloaded from FireBrowse (http://firebrowse.org). These datasets contain expression profiles for 20,501 genes in 451 COAD samples and 94 READ samples. Somatic mutation of genes in COAD samples was determined on the basis of mutation annotation format (MAF) files, which were also downloaded from FireBrowse. Microsatellite status of 294 TCGA COAD patients have been determined by a computational method called MOSAIC using DNA-sequencing data and were downloaded from previous publication (15). MOSAIC determined microsatellite status by a parsimonious, weighted-tree classifier using gain of microsatellite alleles as input with an accuracy rate of 96.6%.
Processed microarray datasets with the accession IDs of GSE39582 (16), GSE41258 (17), and GSE87211 (18) were downloaded from the Gene Expression Omnibus (GEO) database (RRID: SCR_005012; ref. 19). Specifically, GSE39582 contains gene expression profiles for 566 colon cancer samples and 19 normal samples. Clinical parameters, overall survival information, recurrence-free survival information, and mutation status of the TP53, KRAS, and BRAF genes are also provided. TP53 mutations were determined by sequencing exons 4 to 9, and KRAS and BRAF mutations were determined by allelic discrimination assay (16). GSE41258 contains gene expression profiles for 186 primary colon tumor samples, 67 metastatic colorectal tumor samples (liver or lung metastasis), 2 microadenoma samples, 49 polyp samples, 74 normal tissues, and 12 cell lines. GSE87211 contains gene expression profiles for 203 rectum tumor samples and 160 normal mucosa samples. It also provides disease free and disease specific survival time as well as other clinical information. These datasets were processed into expression at the probe set level. We converted probe set expression into gene expression by referring to their corresponding platform annotation files. Some genes are represented by multiple probe sets and in this case the probe set with the maximum average intensity across all samples was selected to represent the gene.
Determine gain or loss of chromosomal bands
TCGA CNV segment files for 451 COAD and 165 READ samples were downloaded from FireBrowse (http://firebrowse.org). These files provide a list of genomic regions identified as significant CNV segments. For each segment, a segment mean was calculated to represent its copy number value. The segment mean values were determined by log2(copy-number/2) with 0 indicating normal copy number, and positive and negative values indicating amplification and deletion, respectively. We selected segments with absolute mean values >0.5 to obtain a confident list of CNV regions. Then we determined the copy number of all chromosome bands based on their overlap with selected segments. A chromosome band is considered gain or loss if it has >50% overlap with a segment. Bands on chromosome Y are excluded from this analysis. Since deletion event is more easily detected compared with amplification event, we normalized deletion frequency by multiplying ratio of total number of amplification events to deletion events.
Infer amplification of Chr20q from expression data
Since matched CNV data are not available for most cancer gene expression data, we adapted a previously proposed method (20) to infer Chr20q amplification status in tumor samples based on gene expression. When a chromosome band is amplified or deleted in a sample, we would expect genes in this band to be upregulated or downregulated. To this end, we selected a subset of Chr20q genes with expression highly correlated with their copy numbers based on TCGA gene expression profile and copy number data. Pearson correlation coefficients between expression levels and copy number of genes on Chr20q were calculated. 102 Genes with correlation coefficients above 0.6 are selected as the signature genes to infer Chr20q amplification status. Chr20q genes were downloaded from the C1 gene set of the MsigDB database (21).
We use gene set A to denote Chr20q signature genes and gene set B to denote all genes not located on Chr20q arm. Given the gene expression profile for a tumor sample, we compared the expression levels between the gene sets A and B using the Student t test. The resulting t-score indicates the upregulation or downregulation of gene set A relative to B:
where |$\overbar x$| denotes average value of expression levels, |${s^2}$| denotes variance and denotes number of genes in the gene set A or B. The t scores were calculated for all samples with the corresponding P values adjusted by using the Benjamini–Hochberg method to correct multiple tests. Adjusted P values less than 0.01 are considered as significant gain (t-score > 0) or loss (t-score < 0) of Chr20q arm.
Calculate immune cell infiltration
A deconvolution algorithm called binding association with sorted expression (BASE; ref. 22) was applied to infer the infiltration of different immune cell types based on tumor gene expression data. BASE is originally designed for inferring transcription factor activity and was later extended to infer immune cell infiltration in tumor samples by our lab (22, 23). This algorithm applies a rank-based method to examine the expression level of immune cell-specific genes in tumor samples and calculates infiltration scores to indicate immune cell infiltration levels. In this study, we focused on six major immune cell types, including memory B cells, naïve B cells, CD4+ T cells, CD8+ T cells, natural killer (NK) cells, and monocytes.
Statistical analysis
Univariate and multivariate Cox regression models were used to examine the association of Chr20q amplification status with patient prognosis. In the multivariant Cox regression model, based on the availability in a dataset the confounding variables were adjusted, including age, sex, tumor stages, tumor locations, and tumor invading depths. Specifically, patients over 70 years old were considered as old. Stage I to II were considered as low stage with III to IV as high stage. Kaplan–Meier plots were used to visualize survival difference between two patient groups (e.g., MSS-A vs. MSS-N). P values were calculated by using the log-rank test.
Survival analysis was performed using the R package “survival” and “survminer”. Specifically, “survfit” and “ggsurvplot” functions were used for Kaplan–Meier plots. “survdiff” function was used to compare the difference between survival curves. “coxph” function was used to build Cox proportional hazard models. Kaplan–Meier curves and forest plot were used to visualize the Cox regression results.
The P values of other analysis were calculated as described in the corresponding sections.
Drug sensitivity data
The Drug sensitivity dataset was downloaded from Genomics of Drug Sensitivity in Cancer (GDSC) database (https://www.cancerrxgene.org, RRID: SCR_011956; ref. 24). From this database, we retrieved the sensitivity data of 22 colon cancer cell lines to a total of 138 drugs. The sensitivity is represented as natural logarithm of the half maximal inhibitory concentration. Microsatellite stability status and gene expression data were also available from this database. Chr20q amplification status of these cell lines was inferred from gene expression data using the method described above. Statistical comparison of drug sensitivity between two cell line groups (MSS-A and MSS-N) was performed by one-sided Wilcoxon rank-sum test.
Differential expression analysis
The R package “limma” was used to identify differentially expressed genes between two sample groups. Since there are more COAD samples than READ samples, we applied different thresholds when select differential expressed genes. Specifically, we adopted threshold as P = 1e−20 for COAD differential expressed genes and P = 1e−7 for READ differential expressed genes. Moreover, due to the higher fold change between tumor and normal samples compared with that between different subtypes, different threshold for log2-fold change were applied, with 0.8 for tumor versus normal and 0.4 for amplification versus nonamplification samples.
Gene essentialities
DepMap is an ongoing project that systematically identifies genetic and pharmacologic dependencies in hundreds of cancer cell lines (25). Gene dependency, represented by probability of knocking out a gene has a real depletion effect, for 19 colorectal cancer cell lines were obtained from datasets of CRISPR (Avana) Public 20Q2 (26).
Results
Chr20q amplification occurs in colorectal cancers with high frequency
The overall workflow is shown in Supplementary Fig. S1. Through pan-cancer analysis on TCGA data, we identified Chr20q amplification as the most frequent chromosome band aberration in colorectal cancers. On the basis of Chr20q amplification and microsatellite status, colorectal cancer cohorts can be further stratified into three subtypes: MSI, MSS-A, and MSS-N. Then we subsequentially investigated Chr20q amplification's association with prognostic value, CIN, mutation profile, immune infiltration, and drug sensitivity. Finally, through integration of diverse data sources, potential cancer driver genes on Chr20q were determined.
To start with, we systematically examined the CNV profiles of 33 cancer types from TCGA to identify chromosome bands that are amplified or deleted with high frequency in colorectal cancer but not in other cancer types. Of all human chromosome bands, we observed that Chr20q1, Chr20q2, and Chr20q3 had the highest amplification frequency in colon adenocarcinoma (COAD) and rectum adenocarcinoma (READ; Fig. 1A). Specifically, 41% of Chr20q1, 44% of Chr20q2, and 41% of Chr20q3 were amplified in COAD, whereas 54%, 58%, and 54% were amplified in READ. In the 293 COAD samples, 56 are known to be MSI tumors. After excluding MSI samples, even higher amplification frequencies were observed for the three chromosome bands (48% for Chr20q1, 51% for Chr20q2, and 48% for Chr20q3). Because only three READ samples were reported to be MSI, the MSS-specific amplification frequencies were not calculated. The amplification frequencies of Chr20q1–3 are much higher than all the other chromosome band amplification or deletion events (Fig. 1A). In fact, the next most frequent amplification is on Chr13q13 with a frequency of 33%, and the most frequent deletion is on Chr18q22 with a frequency of 33%. Since Chr20q1, Chr20q2, and Chr20q3 are the only three bands of the chromosome arm Chr20q and their amplification in COAD and READ samples are highly correlated, we merge them into a single arm-level amplification event of Chr20q. Moreover, our results also suggested the frequent amplification of Chr20q were specific to COAD and READ cancer types: none of the other cancer types has an amplification frequency greater than 40% (Fig. 1B). Given the high frequency of Chr20q amplification compared with the other amplification/deletion events, we postulated that it might play critical roles during the tumorigenesis of colorectal cancer and focus on it in this study.
Chr20q is frequently amplified in colon and rectum cancers. A, Amplification/deletion frequencies of chromosomal bands across 33 cancer types. Red and blue denote amplification and deletion frequency, respectively. B, Chr20q is most frequently amplified in colon and rectum cancers. C, Chr20q amplification rates in different stages of colon and rectum cancer. D, Chr20q amplification rates in MSI and MSS samples. E, A combination of MSI and Chr20q amplification status stratifies colon cancer samples into three subtypes. F, The majority of Chr20q genes are upregulated in MSS-A versus MSS-N samples. G, Distribution of correlation coefficients between expression level and copy number for Chr20q genes. Genes with correlation coefficients over 0.6 were selected to form a signature. H, The t scores calculated by comparing Chr20q signature genes with non-Chr20q genes discriminate MSS-A from MSS-N samples. P value was calculated by using two-sided Wilcoxon rank-sum test. I, ROC for classifying MSS-A versus MSS-N samples based on inferred Chr20q amplification scores. Plots A to C are based on the TCGA data for 33 cancer types. All the other plots are based on the TCGA COAD data unless otherwise specified.
Chr20q is frequently amplified in colon and rectum cancers. A, Amplification/deletion frequencies of chromosomal bands across 33 cancer types. Red and blue denote amplification and deletion frequency, respectively. B, Chr20q is most frequently amplified in colon and rectum cancers. C, Chr20q amplification rates in different stages of colon and rectum cancer. D, Chr20q amplification rates in MSI and MSS samples. E, A combination of MSI and Chr20q amplification status stratifies colon cancer samples into three subtypes. F, The majority of Chr20q genes are upregulated in MSS-A versus MSS-N samples. G, Distribution of correlation coefficients between expression level and copy number for Chr20q genes. Genes with correlation coefficients over 0.6 were selected to form a signature. H, The t scores calculated by comparing Chr20q signature genes with non-Chr20q genes discriminate MSS-A from MSS-N samples. P value was calculated by using two-sided Wilcoxon rank-sum test. I, ROC for classifying MSS-A versus MSS-N samples based on inferred Chr20q amplification scores. Plots A to C are based on the TCGA data for 33 cancer types. All the other plots are based on the TCGA COAD data unless otherwise specified.
We examined the frequencies of Chr20q amplification in TCGA COAD samples at different stages. We found that Chr20q had comparable amplification frequency at different stages with slightly increased frequency in later stages (Fig. 1C).
Chr20q amplification occurs predominantly in MSS colon cancer
We compared Chr20q amplification frequencies between MSI and MSS tumor samples in the TCGA COAD dataset. We found that Chr20q amplification presented in 2 of the 56 MSI samples (4%), compared with 128 of the 237 MSS samples (54%). The frequency was significantly different with a P value of 7e−14 (Fig. 1D). This disproportional distribution suggests that Chr20q amplification may play important roles particularly in MSS colon cancer. On the basis Chr20q amplification status, we propose to stratify MSS tumors into two groups: MSS tumors with Chr20q amplification (MSS-A) and MSS tumors without Chr20q amplification (MSS-N). In combination with MSI status, the TCGA COAD samples were divided into three groups with 128 MSS-A, 109 MSS-N, and 56 MSI samples, respectively (Fig. 1E).
A few high-quality gene expression datasets have been generated in previous studies providing expression profiles for a large number of colorectal cancer samples. In contrast, the CNV data for colorectal cancer samples are quite limited. To overcome this problem, we applied a previously proposed method (20) to determine Chr20q amplification status based on gene expression data. As shown in Fig. 1F, a large fraction of Chr20q genes showed significantly higher expression levels in MSS-A than in MSS-N colorectal cancer samples, suggesting that Chr20q amplification status can be inferred on the basis of the expression of genes from this region. Specifically, we selected 102 genes that are located on Chr20q and display high correlation between their expression levels and copy numbers according to the TCGA COAD data (Fig. 1G). These genes form a signature that can be used to infer Chr20q amplification status by comparing their gene expression with non-Chr20q genes. The resultant t scores and P values indicate the likelihood of tumor sample to have Chr20q amplification event. As shown, MSS-A samples displayed significantly higher t scores than MSS-N samples (Fig. 1H). To evaluate the performance, we applied the t scores to discriminate MSS-A versus MSS-N samples and achieved an AUC of 0.96 in COAD (Fig. 1I). The performance was further validated in READ data by an AUC of 0.91 (Fig. 1I). In practice, we predict samples with t-score >0 and adjusted P value (FDR) <0.01 as MSS-A samples and the others as MSS-N samples.
We then applied this signature to the GSE41258 dataset (17), which contains expression profiles for 186 primary colorectal tumor samples, 67 metastatic colorectal tumor samples, and 49 polyp samples. According to the inferred Chr20q amplification status, 10 of 49 polyp samples, 102 of 186 all primary colorectal tumor samples, 85 of 132 MSS primary colorectal tumor samples, and 43 of 67 metastatic colorectal tumor samples have Chr20q amplification events (Supplementary Fig. S2). The Chr20q amplification frequencies are similar between primary and metastatic samples but significantly different between polyps and cancerous samples (P = 1e−6).
Chr20q amplification is associated with better prognosis in MSS colon cancer
We next examined the prognostic value of the Chr20q amplification status in the GSE39582 data, which contains expression profile for 519 colon cancer samples (Table S1). Since the CNV data are not provided, we inferred the Chr20q amplification status using our gene signature and divided all patients into three groups: 209 MSS-A, 235 MSS-N, and 75 MSI samples. Consistent with previous studies, MSI patients had significant longer recurrence-free survival time than MSS samples (P = 0.008, Fig. 2A; refs. 8, 9). More importantly, MSS-A patients showed better recurrence-free survival rate compared with MSS-N patients (P = 0.03, Fig. 2A). Result from multivariant Cox regression model indicated that the association of Chr20q amplification with prognosis remained significant after considering established clinical variables including age, sex, tumor stage, and location. As shown, MSS-N patients are associated with a hazard ratio of 1.55 compared with MSS-A patients (P = 0.03, Fig. 2B).
Association of Chr20q amplification status with prognosis in colon cancer. A, Kaplan–Meier plots for MSI, MSS-A, and MSS-N patient groups. B, Forest plot for visualizing results from the multivariate Cox regression model. C–F, Survival difference between the MSS-A and MSS-N groups in patient subsets with mutated TP53 (C), wild-type TP53 (D), distal colon cancer (E), and proximal colon cancer (F), respectively.
Association of Chr20q amplification status with prognosis in colon cancer. A, Kaplan–Meier plots for MSI, MSS-A, and MSS-N patient groups. B, Forest plot for visualizing results from the multivariate Cox regression model. C–F, Survival difference between the MSS-A and MSS-N groups in patient subsets with mutated TP53 (C), wild-type TP53 (D), distal colon cancer (E), and proximal colon cancer (F), respectively.
As an important tumor suppressor, the mutation status of the TP53 gene is of great prognostic value in colorectal cancer, and therefore we further stratified patients based on their TP53 mutation status. In TP53 mutated patients, Chr20q amplification was significantly associated with prognosis with MSS-A patients having longer recurrence free survival time than MSS-N patients (P = 6e−4, Fig. 2C). In contrast, no significant survival difference was observed between MSS-A and MSS-N patients in TP53 wild-type patients (Fig. 2D).
Similarly, we performed patient stratification based on tumor location, which represented another important clinical factor in colon cancer (27). In patients with distal colon tumor, MSS-A patients were found to have significantly longer recurrence free survival time than MSS-N patients (P = 0.002, Fig. 2E). In contrast, there was no significant survival differences between MSS-A and MSS-N in patients with proximal colon tumor (Fig. 2F). Taking together, our results suggest that Chr20q amplification status is a useful prognostic factor and provides independent prognostic values in addition to well established clinical variables.
MSS-A and MSS-N colon tumors have differential chromatin instability levels and mutation profiles
To elucidate potential mechanisms by which Chr20q amplification affects patient prognosis in colorectal cancer, we compared MSS-A with MSS-N patients in several important molecular features. First, we examined CIN, which has been reported to be a poor prognostic marker in COAD (13). The CIN levels of TCGA COAD samples have been reported in a previous study (28), which used fraction of altered genomic regions (fraction altered) and aneuploidy score to quantify CIN levels. In addition, a homologous recombination deficiency (HRD) score was calculated to measure the deficiency of the homologous recombination DNA repair pathway, which is the major cause of increased CIN levels in tumors. Our results indicated that MSS-A tumors were associated with significantly increased CIN levels compared with MSS-N tumors as revealed by both fraction altered (P = 8e−11, Fig. 3A) and aneuploidy scores (P = 8e−10, Fig. 3B). Consistently, significantly higher HRD scores were observed in MSS-A than in MSS-N tumors (P = 5e−6, Fig. 3C). It should also be noted that MSI tumors showed much lower CIN levels (fraction altered: P = 4e−23, aneuploidy score: P = 5e−21, Fig. 3A and B) and HRD scores (P = 1e−19, Fig. 3C) than MSS tumors.
Association of Chr20q amplification with genomic aberrations. A and B, MSI, MSS-A, and MSS-N subtypes have differential CIN levels represented as fraction of altered genomic regions (A) and aneuploidy score (B). C, Deficiency of the homologous recombination (HR) DNA repair pathway, represented as HR deficiency (HRD) scores, in MSI, MSS-A, and MSS-N subtypes. D, MSI, MSS-A, and MSS-N subtypes vary in tumor mutation burdens, represented as the total number of nonsynonymous somatic mutations. E, Top 25 most differentially mutated genes between MSS-A and MSS-N. F–H, MSS-A and MSS-N subtypes have significant differential mutation frequencies in TP53 (F), KRAS (G), and BRAF (H). Plots A–E are based on the TCGA COAD data, and F to H are based on the GSE39582 data. P values were calculated from two-sided Wilcoxon rank-sum test (A–D) or Fisher exact test (E–H).
Association of Chr20q amplification with genomic aberrations. A and B, MSI, MSS-A, and MSS-N subtypes have differential CIN levels represented as fraction of altered genomic regions (A) and aneuploidy score (B). C, Deficiency of the homologous recombination (HR) DNA repair pathway, represented as HR deficiency (HRD) scores, in MSI, MSS-A, and MSS-N subtypes. D, MSI, MSS-A, and MSS-N subtypes vary in tumor mutation burdens, represented as the total number of nonsynonymous somatic mutations. E, Top 25 most differentially mutated genes between MSS-A and MSS-N. F–H, MSS-A and MSS-N subtypes have significant differential mutation frequencies in TP53 (F), KRAS (G), and BRAF (H). Plots A–E are based on the TCGA COAD data, and F to H are based on the GSE39582 data. P values were calculated from two-sided Wilcoxon rank-sum test (A–D) or Fisher exact test (E–H).
Next, we compared the tumor mutation burdens between MSS-A and MSS-N colon tumors in the TCGA data. Tumor mutation burden is calculated as the total number of nonsynonymous mutations in a tumor sample. As shown in Fig. 3D, no significant difference in tumor mutation burdens was observed between MSS-A and MSS-N tumors, whereas as expected, MSI tumors showed significantly higher tumor mutation burdens than MSS tumors. Following that, we identified all genes with significantly different mutation frequency between MSS-A and MSS-N tumors. In Fig. 3E, we showed the top 25 most differently mutated genes. As shown, TP53 had significantly higher mutation frequency in MSS-A tumors, whereas KRAS and BRAF had significantly higher mutation frequency in MSS-N tumors. These results were further validated from the GSE39582 dataset, in which the mutation status of TP53, KRAS, and BRAF had been determined (Fig. 3F–H).
Chr20q amplification is associated with lower immune infiltration in MSS colon cancer
Previous studies reported the association between aneuploidy and immune infiltration (14). Having shown the difference between MSS-A and MSS-N in aneuploidy scores, we then compared their immune infiltration. First, we compared the infiltration scores of leukocytes and lymphocytes in TCGA COAD samples, which were calculated by computational analysis based on a previous study (28). MSS-N tumors were found to have significantly higher infiltration levels of both leukocytes (P = 1e−6, Fig. 4A) and lymphocytes (P = 4e−7, Fig. 4B) than MSS-A tumors. As expected, MSI tumors had much higher infiltrations levels compared with MSS tumors. Specifically, comparison between MSI and MSS-N tumors showed that the leukocyte fraction was lower in MSS-N than MSI with moderate significance (P = 0.003), and lymphocyte infiltration was similar between MSS-N and MSI (P > 0.1). These results suggested that MSS-N tumors are more similar to MSI tumors than to MSS-A tumors in terms of leukocyte and lymphocyte infiltration.
Association of Chr20q amplification with immune infiltration. A and B, Differential infiltration levels of leukocytes (A) and lymphocytes (B) among MSI, MSS-A, and MSS-N subtypes. C and D, Differential infiltration levels of naïve B cells, memory B cells, CD8+ T cells, and monocytes among the three subtypes in the TCGA COAD data (C) and the GSE39582 data (D). E, Relative expression levels of immune-related genes among the three subtypes. Plots A, B, C, and E are based on the TCGA COAD data. P values were calculated by two-sided Wilcoxon rank-sum test.
Association of Chr20q amplification with immune infiltration. A and B, Differential infiltration levels of leukocytes (A) and lymphocytes (B) among MSI, MSS-A, and MSS-N subtypes. C and D, Differential infiltration levels of naïve B cells, memory B cells, CD8+ T cells, and monocytes among the three subtypes in the TCGA COAD data (C) and the GSE39582 data (D). E, Relative expression levels of immune-related genes among the three subtypes. Plots A, B, C, and E are based on the TCGA COAD data. P values were calculated by two-sided Wilcoxon rank-sum test.
Next, we examined which specific immune cell types contribute to the difference between MSS-A and MSS-N tumors. Specifically, by using a previous developed deconvolution algorithm we calculated and compared six major immune cell types, naïve B cell, memory B cell, CD4+ T cell, CD8+ T cell, NK cell, and monocyte, based on TCGA COAD gene expression data. Our results indicated that the infiltration levels of naïve B cell, memory B cell, CD8+ T cell were significantly higher in MSS-N than in MSS-A tumors in TCGA COAD samples (Fig. 4C). These results were further validated in the GSE39582 dataset for Naïve B cell and CD8+ T cell (Fig. 4D). For memory B cell, the same trend was observed but did not reach significance threshold.
In addition, we examined the expression pattern of 70 immune-related genes in MSI, MSS-A, and MSS-N tumors. Included in this gene list are major immune stimulatory and inhibitory genes. As shown in Fig. 4E, these genes tend to have higher expression levels in the MSI tumors than in the MSS tumors. This is consistent with previous studies and supports that MSI tumors have a hot immune environment. Within the MSS tumors, MSS-N tumors displayed higher expression levels in immune-related genes compared with MSS-A tumors.
Taking together, these results suggest that in colon cancer MSS tumors tend to have cold immunity compared with MSI tumors. There is significant difference between MSS-A and MSS-N tumors with the later having significantly higher immune cell infiltration and immune gene expression.
The effect of Chr20q amplification in rectum cancer
As shown in Fig. 1A, the frequency of Chr20q amplification in rectum cancer is even higher than that in colon cancer. We have been focused on colon cancer in the previous sections. In this section, we extend our analyses to rectum cancer to examine the association of Chr20q amplification with prognosis, genomic features, and immune cell infiltration. Unlike colon cancer, only a small subset of rectum cancers is associated with MSI. Specifically, out of the 95 TCGA READ samples, only three are MSI tumors. As such, we simply divided rectum tumors into two groups (Chr20q-A and Chr20q-N) without considering MSI status.
First, we performed survival analysis using the GSE87211 dataset, which contained gene expression profiles for a total of 203 rectum tumor samples (Table S1). On the basis of inferred Chr20q amplification status, we separate these samples into a Chr20q-A group (amplification-positive) and a Chr20q-N group (amplification-negative) with 80 and 123 samples. We found that Chr20q-A patients had significantly longer disease-free survival time than Chr20q-N patients (P = 0.05, Fig. 5A), which is consistent with results obtained in colon cancer. After adjusting several confounding variables including age, gender, and invading depth, Chr20q amplification status remained significant according to the multivariant Cox regression model (Fig. 5B).
Chr20q amplification status is associated with prognosis in rectum cancer. A, Chr20q-A patients have significantly longer survival than Chr20q-N patients. B, Forest plot for results from multivariate Cox regression model. Note that Chr20q amplification status is significant after adjusting clinical variables. C–E, CIN characterized by aneuploidy score (C), fraction altered (D), and HRD score (E) between Chr20q-A and Chr20q-N patients in TCGA READ dataset. F and G, Immune infiltration levels of leukocytes (F) and lymphocytes (G) between Chr20q-A and Chr20q-N patients. H, Memory B cell, CD8+ T cell, and nature killer cell were found to be infiltrated in Chr20q-A and Chr20q-N differently at statistical significance. Plots A and B are based on the GSE87211 data, and plots C to H are based on TCGA READ data. All P values for box plots are calculated from two-sided Wilcoxon rank-sum test.
Chr20q amplification status is associated with prognosis in rectum cancer. A, Chr20q-A patients have significantly longer survival than Chr20q-N patients. B, Forest plot for results from multivariate Cox regression model. Note that Chr20q amplification status is significant after adjusting clinical variables. C–E, CIN characterized by aneuploidy score (C), fraction altered (D), and HRD score (E) between Chr20q-A and Chr20q-N patients in TCGA READ dataset. F and G, Immune infiltration levels of leukocytes (F) and lymphocytes (G) between Chr20q-A and Chr20q-N patients. H, Memory B cell, CD8+ T cell, and nature killer cell were found to be infiltrated in Chr20q-A and Chr20q-N differently at statistical significance. Plots A and B are based on the GSE87211 data, and plots C to H are based on TCGA READ data. All P values for box plots are calculated from two-sided Wilcoxon rank-sum test.
Second, we examined the correlation between Chr20q amplification and the CIN level and mutation profiles in TCGA READ samples. Consistent with observations in colon cancer, our results indicated that Chr20q-A rectum tumors tend to have higher CIN levels compared with Chr20q-N patients (aneuploidy score: P = 2e−4, fraction altered: P = 0.005, Fig. 5C and D) and more deficient homologous recombination DNA repair pathway (HRD score: P = 0.01, Fig. 5E). TP53 mutation was enriched in rectum tumors with Chr20q amplification (Chr20q-A; Supplementary Fig. S3A, P = 0.04). For mutation frequency of other genes including KRAS (P > 0.1) and BRAF (P = 0.09), the same trend as those in colon cancer was observed but did not reach significance threshold (Supplementary Figs. S3B and S3C). In addition, similar levels of tumor mutation burden were observed between Chr20q-A and Chr20q-N tumors (Supplementary Fig. S3D).
Finally, we compared the immune infiltrations between Chr20q-A and Chr20q-N tumors in rectum cancer. We found that Chr20q-N rectum tumors had higher fraction of leukocytes (P = 0.03, Fig. 5F). The same trend was observed between Chr20q-N and Chr20q-A samples, but the difference did not reach significance threshold, presumably due to smaller sample size of the READ dataset (Fig. 5G). More specifically, the infiltration levels of memory B cells, CD8+ T cells, and NK cells were lower in Chr20q-A rectum tumors (Fig. 5H).
Association of Chr20q amplification with drug sensitivity in colorectal cancer cell lines
Given the potential roles of Chr20q amplification during tumorigenesis in colorectal cancer, we then examined whether its status had impact on drug sensitivity of colorectal cancer cell line. Specifically, we selected 22 MSS colorectal cancer cell lines from the GDSC dataset (24) and inferred their Chr20q amplification status (Supplementary Table S2). Then, we compared the sensitivity of cells to a total of 138 drugs between the MSS-A and MSS-N colon cancer cell lines. Tipifarnib, Mitomycin-C, and AUY922 were identified as the top 3 significant drugs that were more effective in MSS-N than in MSS-A cell lines (Fig. 6A). Tipifarnib is a farnesyltransferase inhibitor targeting Ras signaling pathway (29). KRAS and BRAF with high mutation frequency in MSS-N are both in this pathway. Mitomycin-C is a DNA cross linker, thereby inhibit transcription in tumor cells. A possible reason underlining different efficacy of mitomycin-C is high CIN level in MSS-A subtypes.
Drug efficacy difference between MSS-A and MSS-N cell lines. A, Top three drugs more effective in MSS-N cell line. B, Top three drugs more effective in MSS-A cell lines. Drug efficacies are represented by natural logarithm of the IC50 [ln(IC50)]. All P values for box plots were calculated from one-sided Wilcoxon rank-sum test.
Drug efficacy difference between MSS-A and MSS-N cell lines. A, Top three drugs more effective in MSS-N cell line. B, Top three drugs more effective in MSS-A cell lines. Drug efficacies are represented by natural logarithm of the IC50 [ln(IC50)]. All P values for box plots were calculated from one-sided Wilcoxon rank-sum test.
On the other hand, EHT-1864, AP-24534, and NSC87877 were identified as the top 3 significant drugs more effective in MSS-A than MSS-N cell lines (Fig. 6B). Among them, AP-24534 is an inhibitor of FGFR, a therapeutic target in many cancer types (30). FGFR involves in similar pathways of EGFR, high frequency of EGFR mutation may render MSS-A subtypes sensitivity to AP-24534.
Potential cancer driver genes on Chr20q
The Chr20q arm contains a total of 324 genes. Twenty-five genes with zero expression levels were excluded from this analysis. To identify the candidate driver genes from this arm, we (i) calculate gene-specific amplification frequency in TCGA COAD and READ samples, (ii) compared gene expression between tumor and adjacent normal tissues in COAD and READ, (iii) compared gene expression between MSS-A and MSS-N COAD samples and between Chr20q-A and Chr20q-N READ samples, (iv) determined the association of gene expression with patient prognosis in colon (GSE39582) and rectum cancer (GSE87211), and (v) examined the essentiality of genes in colorectal cancer cell line as reported in the DepMap database (26).
All genes located on Chr20q show similar amplification frequency in both COAD and READ, indicating their amplifications are caused by a shared arm-level amplification event of Chr20q. Out of the 299 genes, 29 and 38 genes are significantly upregulated in tumor than adjacent normal samples in COAD (adjusted P < 1e−20 and log2-fold change >0.8) and READ (adjusted P < 1e−7 and log2-fold change >0.8), respectively; 52 are significantly upregulated (P < 1e−20 and log2-fold change >0.4) in MSS-A than MSS-N COAD subtypes; 78 are significantly upregulated in Chr20q-A than Chr20q-N READ subtypes (P < 1e−7 and log2-fold change >0.4); 79 are significantly associated with prognosis (P < 0.05) in colon cancer with 7 being hazardous (HR>1) and 72 being protective (HR<1); 29 are significantly associated with prognosis (P < 0.05) in rectum cancer with 5 being hazardous and 24 being protective. By integrating these results, we selected 37 genes and visualized them in Fig. 7. In addition, essentialities of these 37 genes across 19 colorectal cancer cell lines downloaded from DepMap database were visualized in Fig. 7 as well as in ref. 26.
Potential cancer driver genes on Chr20q arm. Amplification frequencies of 299 genes on Chr20q arm of colon and rectum cancer are shown in top two curves, respectively. Genes are ordered by their positions on Chr20q arm along x-axis. Probability of genes affecting cell viability indicated by dependency of selected 37 genes across 19 colorectal cancer cell lines from DepMap are visualized in middle heatmap. Comparison/association analyses of selected 37 genes are visualized in the bottom heatmap. All significant P values were divided into three categories as low, medium, and high significance. Red gene names indicate possible candidate driver genes. These genes were selected on the basis of both significant essentiality from DepMap result and significance in at least four of the six comparison/association analyses.
Potential cancer driver genes on Chr20q arm. Amplification frequencies of 299 genes on Chr20q arm of colon and rectum cancer are shown in top two curves, respectively. Genes are ordered by their positions on Chr20q arm along x-axis. Probability of genes affecting cell viability indicated by dependency of selected 37 genes across 19 colorectal cancer cell lines from DepMap are visualized in middle heatmap. Comparison/association analyses of selected 37 genes are visualized in the bottom heatmap. All significant P values were divided into three categories as low, medium, and high significance. Red gene names indicate possible candidate driver genes. These genes were selected on the basis of both significant essentiality from DepMap result and significance in at least four of the six comparison/association analyses.
From these genes, we further identified 10 genes that showed significant essentiality according to the DepMap result (26) and were significant in at least four of the six comparison/association analyses described above. These genes are TPX2, PDRG1, AHCY, RPN2, TP53RK, CSE1L, MOCS3, AURKA, RAE1, and PSMA7. TPX2 and AURKA were previously reported to promote amplification of Chr20q synergically (31) and drive tumorigenesis in colorectal cancer (32). PDRG1 and TP53RK are p53-related genes that regulate cellular cycle and apoptosis (33, 34). CSE1L, RPN2, and PSMA7 have also been reported to promote proliferation, cellular migration, and metastasis (35–37).
Discussion
Precision therapy is now advancing cancer treatments, especially for colorectal cancer, the third leading cause of cancer-related death in United States. Cancer genotyping is an important way to divide patients into subgroups with distinct molecular profiles, thereby provides accurate prognostic prediction and personalized treatments. In this study, we identified exclusive high amplification frequency of Chr20q arm in MSS patients with colorectal cancer. Stratifying patients with colorectal cancer based on Chr20q amplification, we found Chr20q amplification accurately predict recurrence-free survival rates within both colon and rectum cancer patients. Furthermore, we observed Chr20q amplification is associated with distinct chromosomal stability, mutation profiles, and immune infiltration levels, which is correlated with a previous publication (38). To examine clinical impacts of Chr20q amplification, we compared 138 drugs' efficacies between MSS-A and MSS-N subtypes and identified top three useful drugs for each subtype. Altogether, these results indicate the ability of Chr20q amplification in stratifying colorectal cancers.
In addition, we observed twenty percent of polyps have Chr20q amplification as well (Supplementary Fig. S2). This suggests Chr20q occur at early stage and could be potentially helpful in distinguishing benign and malignant polyps. Moreover, to identify potential cancer driver genes on Chr20q, we integrate different analysis and DepMap results to select 10 most possible genes (TPX2, PDRG1, AHCY, RPN2, TP53RK, CSE1L, MOCS3, AURKA, RAE1, and PSMA7), which are further proved to be involved in tumorigenesis of colorectal cancer according to previous publications (31–37). Ptashkin and colleagues has reported the SRC gene as one of the driver gene according to the high correlation between DNA copy number and expression level, but it was not included in our list (38). Although we verified that the DNA copy number of SRC is correlated with its expression level in both TCGA COAD (R = 0.60) and READ (R = 0.60) cohorts, we found that SRC was not significantly differential expressed between tumor and normal tissues in both COAD (log2-fold change = 0.02, P > 0.1) and READ (log2-fold change = 0.62, P = 0.002). In addition, SRC is not associated with prognosis in both COAD (P > 0.1) and READ (P > 0.1) and has low dependency scores across the 19 colorectal cancer cell lines in DepMap dataset with the highest score is 0.4. Therefore, the SRC gene reported to be potential driver genes by Ptashkin and colleagues was not included in our list (38).
This study has the following limitations. First, the inferred chr20q amplification status was used to determine its correlation with prognosis. This result would be consolidated in colon cancer data with directly measured copy number data and matched survival information. Second, the data were analyzed retrospectively from public databases, more work is needed to validate the results prospectively in patients with colorectal cancer. Third, the drug sensitivity analysis was based on a relatively small number of COAD MSS cell lines. The difference between cell lines and tumors might further confound our results. Large-scale clinical trial data will be helpful to verify the drugs identified in our cell line-based analysis.
It is interesting that we observed higher immune infiltration levels were associated with MSS-N patients who have poor survival rate. Specifically, CD8+ T cell, usually correlated with good prognosis (39), was found to be highly infiltrated in Chr20q nonamplification samples from both colon and rectum cancer as well. Recent studies illustrated heterogeneity of several important immune cells including monocytes and CD4+ T cells (40, 41). Subsets of these cells have different, even opposing roles in tumor microenvironment. For example, two main subsets of CD4+ T cells are regulator T cells and helper T cells (41). Regulator T cells can suppress immune responses and helper T cells can polarize immune responses. It is possible that higher fraction of regulator T cells exist in tumor microenvironments of MSS-N patients, which explains the association of CD4+ T cell with their shorter poor survival. Therefore, more nuanced research will be needed to verify this finding.
In summary, we have identified Chr20q amplification as a distinct genomic alteration to stratify and predict recurrence-free survival in patients with colorectal cancer. This discovery gives new insights into cancer detection and precision therapy in colorectal cancers.
Authors' Disclosures
No disclosures were reported.
Authors' Contributions
B. Zhang: Formal analysis, validation, investigation, visualization, writing–original draft. K. Yao: Writing–review and editing. E. Zhou: Writing–review and editing. L. Zhang: Writing–review and editing. C. Cheng: Supervision, funding acquisition, project administration, writing–review and editing.
Acknowledgments
This work was supported by the Cancer Prevention Research Institute of Texas (CPRIT; RR180061 to C. Cheng) and the NCI of the NIH (1R21CA227996 to C. Cheng). C. Cheng is a CPRIT Scholar in Cancer Research.
The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.