Purpose:

Irinotecan/5-fluorouracil (5-FU; FOLFIRI) or oxaliplatin/5-FU (FOLFOX), combined with bevacizumab or cetuximab, are approved, first-line treatments for metastatic colorectal cancer (mCRC). We aimed at identifying germline variants associated with survival in patients with mCRC treated with these regimens in Cancer and Leukemia Group B/SWOG 80405.

Experimental Design:

Patients with mCRC receiving either FOLFOX or FOLFIRI were randomized to either cetuximab or bevacizumab. DNA from peripheral blood was genotyped for approximately 700,000 SNPs. The association between SNPs and overall survival (OS) was tested in 613 patients of genetically estimated European ancestry using Cox proportional hazards models.

Results:

The four most significant SNPs associated with OS were three haplotypic SNPs between microsomal glutathione S-transferase 1 (MGST1) and LIM domain only 3 (LMO3, representative HR, 1.56; P = 1.30 × 10−6), and rs11644916 in AXIN1 (HR, 1.39, P = 4.26 × 10−6). AXIN1 is a well-established tumor suppressor gene in colorectal cancer, and rs11644916 (G>A) conferred shorter OS. Median OS for patients with the AA, AG, or GG genotypes was 18.4, 25.6, or 36.4 months, respectively. In 90 patients with stage IV colorectal cancer from The Cancer Genome Atlas (TCGA), rs11649255 in AXIN1 [in almost complete linkage disequilibrium (LD) with rs11644916], was associated with shorter OS (HR, 2.24, P = 0.0096). Using rs11648673 in AXIN1 (in very high LD with rs11644916 and with functional evidence), luciferase activity in three colorectal cancer cell lines was reduced.

Conclusions:

This is the first large genome-wide association study ever conducted in patients with mCRC treated with first-line standard treatment in a randomized phase III trial. A common SNP in AXIN1 conferred worse OS and the effect was replicated in TCGA. Further studies in colorectal cancer experimental models are required.

Translational Relevance

This article demonstrates, for the first time, that common germline variation in AXIN1 in patients with metastatic colorectal cancer (mCRC) may affect patient survival. This article also provides new results about the clinical role of the entire genomic variation in patients. If confirmed in additional studies, these results can be used for reaching decisions to obtain optimized treatment in patients with mCRC.

Cancer and Leukemia Group B (CALGB)/SWOG 80405 was designed to test whether cetuximab or bevacizumab or both added to 5-fluorouracil (5-FU)/leucovorin with either oxaliplatin (FOLFOX) or irinotecan (FOLFIRI) lead to superior overall survival (OS) as first-line therapy in advanced or metastatic colorectal cancer (mCRC). CALGB is now a part of the Alliance for Clinical Trials in Oncology (Alliance). In the primary cohort patients, such as patients with KRAS wild-type (WT) for codon 12/13 and treated with either cetuximab or bevacizumab, there was no statistically significant difference in OS between the bevacizumab arm and the cetuximab arm (1). Over the years from the start of the trial, the combinations of those biologics with chemotherapy became the current standard-of-care treatment in the first-line setting of mCRC.

As part of the protocol, prospective collection of germline DNA from peripheral blood had been obtained from patients who provided informed consent. The analysis of the common germline variation of patients with cancer can allow discoveries of new biological determinants of outcome. Although most of the attention for biomarker discovery efforts in oncology has been focused on the somatic genome, many hallmarks of cancer lie outside of the cancer genome. They include angiogenesis, inflammation, adaptive immunity, and others (2, 3).

We hypothesized that novel germline variants that drive tumor progression and/or resistance to therapy might affect patient outcome. Outcome data and DNA from peripheral blood collected from this relatively large study might provide an unprecedented opportunity for identifying germline DNA markers that affect OS. In this study, the germline DNA of patients was scanned using high-density SNP chips that detect hundreds of thousands of SNPs. We used this approach to identify novel genes associated with OS in patients with mCRC. In addition, we used The Cancer Genome Atlas (TCGA) data and in vitro experiments to corroborate the clinical findings from CALGB/SWOG 80405.

CALGB/SWOG 80405

The trial was designed to compare either cetuximab, bevacizumab, or cetuximab/bevacizumab combined with FOLFOX or FOLFIRI as first-line treatment of advanced and mCRC. Within 3 years, interim analysis of outcome data led to a restriction of eligibility to patients with KRAS WT (codon 12 and 13) and closure of the dual antibody arm. A revised two-arm trial (chemotherapy plus either cetuximab or bevacizumab) became the “primary cohort.” Randomization was stratified by chemotherapy, prior adjuvant chemotherapy, and prior pelvic radiation. Eligibility criteria can be found in the original article (1).

In this genomic study, the primary endpoint of OS was defined as time of randomization until death due to any cause. Patients without reported deaths were censored at their last known follow-up. The median follow-up was 60.4 months [95% confidence interval (CI), 56.3–64.3].

Table 1 reports demographics and clinical characteristics of patients in this study, who belong to a subset of the primary cohort. Institutional review board approval was obtained for this study and all participating patients provided written informed consent to this analysis. Patient registration, data collection, and data analysis were performed by the Alliance Statistics and Data Center. Data quality was ensured by careful review of data by Alliance Statistics and Data Center staff and the study chair of the clinical trial. This study was conducted in accordance with the ethical guidelines from the Declaration of Helsinki.

Table 1.

Clinical characteristics and demographics of patients genotyped from the primary cohort of CALGB/SWOG 80405.

GWASPrimary cohort
(N = 613)(N = 1,137)
Age (years) 
 Median 59.2 59.1 
 Range (21.8–84.3) (20.8–89.5) 
Arm 
 Chemotherapy/bevacizumab 312 (50.9%) 559 (49.2%) 
 Chemotherapy/cetuximab 301 (49.1%) 578 (50.8%) 
Chemotherapy 
 FOLFOX 459 (74.9%) 835 (73.4%) 
 FOLFIRI 154 (25.1%) 302 (26.6%) 
Prior adjuvant chemotherapy 
 No 534 (87.1%) 977 (85.9%) 
 Yes 79 (12.9%) 160 (14.1%) 
Prior pelvic radiation 
 No 554 (90.4%) 1,035 (91.0%) 
 Yes 59 (9.6%) 102 (9.0%) 
Gender 
 Male 380 (62.0%) 697 (61.3%) 
 Female 233 (38.0%) 440 (38.7%) 
Race 
 More than one race 
 Unknown 2 (0.3%) 24 (2.1%) 
 White 611 (99 8%) 934 (82.1%) 
 African American 0 (0.0%) 129 (11.3%) 
 Asian 0 (0.0%) 35 (3.1%) 
 Native Hawaiian or Pacific Islander 0 (0.0%) 3 (0.3%) 
 American Indian or Alaska Native 0 (0.0%) 6 (0.5%) 
 Not reported 0 (0.0%) 1 (0.1%) 
ECOG PS 
 0 363 (59.2%) 657 (57.8%) 
 1 250 (40.8%) 478 (42.0%) 
 2 0 (0.0%) 2 (0.2%) 
Number of metastatic sites 
 Missing 18 42 
 1 344 (57.8%) 635 (55.8%) 
 2 190 (31.9%) 355 (31.2%) 
 3+ 61 (10.3%) 105 (9.2%) 
Tumor location 
 Left 375 (61.2%) 689 (60.6%) 
 Right 154 (25.1%) 280 (24.6%) 
 Transverse 36 (5.9%) 62 (5.5%) 
 Multiple 0 (0.0%) 1 (0.1%) 
 Unknown 48 (7.8%) 105 (9.2%) 
Liver metastases only 
 Missing 
 No 367 (60.3%) 683 (60.1%) 
 Yes 242 (39.7%) 445 (39.1%) 
OS (months) 
 Median (95% Cl) 29.5 (26.7–32.6) 29.4 (27.6–31.4) 
PFS (months) 
 Median (95% Cl) 10.6 (9.7–11.3) 10.6 (9.8–11.1) 
GWASPrimary cohort
(N = 613)(N = 1,137)
Age (years) 
 Median 59.2 59.1 
 Range (21.8–84.3) (20.8–89.5) 
Arm 
 Chemotherapy/bevacizumab 312 (50.9%) 559 (49.2%) 
 Chemotherapy/cetuximab 301 (49.1%) 578 (50.8%) 
Chemotherapy 
 FOLFOX 459 (74.9%) 835 (73.4%) 
 FOLFIRI 154 (25.1%) 302 (26.6%) 
Prior adjuvant chemotherapy 
 No 534 (87.1%) 977 (85.9%) 
 Yes 79 (12.9%) 160 (14.1%) 
Prior pelvic radiation 
 No 554 (90.4%) 1,035 (91.0%) 
 Yes 59 (9.6%) 102 (9.0%) 
Gender 
 Male 380 (62.0%) 697 (61.3%) 
 Female 233 (38.0%) 440 (38.7%) 
Race 
 More than one race 
 Unknown 2 (0.3%) 24 (2.1%) 
 White 611 (99 8%) 934 (82.1%) 
 African American 0 (0.0%) 129 (11.3%) 
 Asian 0 (0.0%) 35 (3.1%) 
 Native Hawaiian or Pacific Islander 0 (0.0%) 3 (0.3%) 
 American Indian or Alaska Native 0 (0.0%) 6 (0.5%) 
 Not reported 0 (0.0%) 1 (0.1%) 
ECOG PS 
 0 363 (59.2%) 657 (57.8%) 
 1 250 (40.8%) 478 (42.0%) 
 2 0 (0.0%) 2 (0.2%) 
Number of metastatic sites 
 Missing 18 42 
 1 344 (57.8%) 635 (55.8%) 
 2 190 (31.9%) 355 (31.2%) 
 3+ 61 (10.3%) 105 (9.2%) 
Tumor location 
 Left 375 (61.2%) 689 (60.6%) 
 Right 154 (25.1%) 280 (24.6%) 
 Transverse 36 (5.9%) 62 (5.5%) 
 Multiple 0 (0.0%) 1 (0.1%) 
 Unknown 48 (7.8%) 105 (9.2%) 
Liver metastases only 
 Missing 
 No 367 (60.3%) 683 (60.1%) 
 Yes 242 (39.7%) 445 (39.1%) 
OS (months) 
 Median (95% Cl) 29.5 (26.7–32.6) 29.4 (27.6–31.4) 
PFS (months) 
 Median (95% Cl) 10.6 (9.7–11.3) 10.6 (9.8–11.1) 

Abbreviation: ECOG PS, Eastern Cooperative Oncology Group performance status.

Genotyping data

DNA was extracted from peripheral blood. The first batch of genotyping was performed on the Illumina HumanOmniExpress-12v1 platform at the Riken Institute (Tokyo, Japan), and included 731,412 SNPs. The second batch of genotyping was performed on the Illumina HumanOmniExpress-8v1 and included 964,193 SNPs. A total of 719,461 SNPs from HapMap from batch 1 were also on the chip from batch 2. A total of 1,178 of these were excluded for having mismatched annotation between the two platforms. Among the 718,283 SNPs in common between the two platforms, 70,108 SNPs were excluded because of call rates less than 99%, 107,761 SNPs with minor allele frequencies (MAFs) less than 0.05 were removed, and 77 SNPs with strong evidence for departure from Hardy–Weinberg equilibrium (P < 10−8) were also removed. Among the 540,337 SNPs passing the filter, 540,021 SNPs were autosomal, and were used in the association analyses. Patients whose genotypic data were identified as unintentional duplicates, or whose genotypically estimated sex did not match their self-reported gender were also excluded. Genetic ancestral origin of patients was estimated using principal components analysis prior to filtering the variants (4). This resulted in 613 patients of European ancestry in the primary cohort (Supplementary Fig. S1).

Bioinformatic analysis of SNP function

Putative functional effects of the four most significant SNPs associated with OS and variants in high linkage disequilibrium (LD, r2 > 0.8 in Europeans] were investigated using SNPnexus (5), which integrates resources from ENCODE and the Roadmap Epigenomics Consortium to map and annotate noncoding variants onto different classes of regulatory regions and predict their functional impact using eight noncoding variant scoring algorithms and computational methods. Variants were also examined for annotation of putative function using RegulomeDB (https://www.regulomedb.org/regulome-search/) and HaploReg v4.1 (www.broadinstitute.org/mammals/haploreg/haploreg.php).

Statistical analysis

Our primary analysis determined the association between SNPs and OS in both arms combined, in patients of genetically estimated European ancestry in the primary cohort (N = 613). Patients with genetically estimated ancestries other than European were not included in this study. The Cox score (log-rank) test was used, and the analysis was conducted under an additive genetic model. To assess potential confounding of baseline risk factors, we conducted an adjusted analysis within the framework of a multivariable additive log-linear Cox proportional hazards model (6). These covariates were age (log10 transformed), sex, and tumor location (left vs. right/transverse). In an alternative adjusted analysis, tumor location was replaced by BRAF V600E mutation status, microsatellite instability status (MSI), and consensus molecular subtype (CMS). For the unadjusted analyses, the HR of the genetic effect and a 95% CI were reported, along with the P value for the score test. For the adjusted analyses, the Wald P value of the SNP effect was reported instead. As an exploratory analysis, an interaction model was used to test for differences in SNP effect between the two treatment arms (cetuximab and bevacizumab). The interaction model was adjusted for age (log10 transformed), sex, and tumor location (left vs. right/transverse) and stratified by chemotherapy. The Wald P value of the interaction term was reported, along with the HR and 95% CI of the genetic effect in each treatment arm (estimated separately). The interaction model was repeated using progression-free survival (PFS, measured as reported previously; ref. 7) as the outcome. The interaction between SNP and chemotherapy was also tested, stratified by treatment arm.

As secondary analyses, selected SNPs from the genome-wide analysis were tested for independence from tumor location (left vs. right/transverse), BRAF V600E, and MSI using generalized linear models, and independence from CMS using the χ2 test. An analysis of deviance was also conducted to determine whether tumor location contributed additional information to the survival model beyond that of each SNP alone. Analyses involving covariates were limited to the subset of patients whose data were available. The P values and CIs were not adjusted for multiple comparisons, as they were used as a feature selection for additional testing in TCGA and in vitro experiments. LD r2 values were reported on the basis of the European population from the 1000 Genomes Project, unless otherwise noted.

A quantile–quantile plot was generated to illustrate the distribution of P values from the primary analysis (Supplementary Fig. S2). Kaplan–Meier and LocusZoom (8) plots were generated to illustrate the association of selected SNPs with survival. The R software environment for statistical computing and graphics (9) and the survival (10) extension package were used to conduct the statistical analyses and generate the figures. The analyses were carried out with adherence to the principles of reproducible analysis using the knitr package (11) for generation of dynamic reports and Duke's GitLab (https://gitlab.oit.duke.edu) for source code management. The code for replicating the statistical analysis was made accessible through a public source code repository (https://gitlab.oit.duke.edu/dcibioinformatics/pubs/calgb-80405-gwas).

TCGA colorectal cancer replication analysis

We used OS and germline genetic data from patients with TCGA who had stage IV colon adenocarcinoma and rectum adenocarcinoma (n = 90; ref. 12). rs11644916 in AXIN1 was not genotyped in TCGA. Hence, we used a proxy, such as rs11649255 (A>G) in AXIN1, which is in almost complete LD with rs11644916 (r2 = 0.97). The association between rs11649255 and OS was assessed using a Cox proportional hazards model. Because there were only three patients that were homozygous for the rs11649255 variant allele, a resampling was performed using 100,000 permutations without replacement. rs10772941 was also tested for association with OS, as a proxy for rs10846370 [r2 = 0.87, intergenic between microsomal glutathione S-transferase 1 (MGST1) and LIM domain only 3 (LMO3)].

Cell culture, plasmids, and luciferase assay

HCT116, SW480, and RKO colorectal cancer cell lines were obtained from the ATCC. HCT116 and SW480 cells were grown in McCoy's 5A (Mediatech) supplemented with 10% FBS (Omega Scientific, Inc.) and 1% penicillin/streptomycin, and incubated at 37°C and 5% CO2. RKO cells were grown in DMEM (Mediatech) supplemented with 10% FBS (Omega Scientific, Inc.) and 1% penicillin/streptomycin, and incubated at 37°C and in 5% CO2.

To test the effects of the AXIN1 variant on transcriptional efficiency, DNA fragments corresponding to the candidate enhancer region in AXIN1 (chr16:327363-328800, hg38) were amplified by PCR from normal human genomic DNA using the following forward and reverse primers, respectively: 5′-GGGCAAGAGAATGAGTCACC-3′ and 5′-CATGAGTGTTTGGGTTCCTG-3′. PCR fragments were subcloned into the SacII restriction enzyme site in both directions, upstream of a thymidine kinase (TK) minimal promoter firefly luciferase vector obtained courtesy of Dr. G.A. Coetzee (Van Andel Research Institute, Grand Rapids, MI). Fragments were cloned using the In-Fusion HD Cloning Kit (Takara). rs11648673 (in LD with rs11644916, r2 = 0.87) was used instead of rs11644916 in the assays because of stronger functional prediction based upon the bioinformatic analysis (see Results section). The AXIN1 functional variant was created using site-directed mutagenesis. Plasmid clones were sequenced by Sanger sequencing to confirm the presence of the candidate variants and the absence of any PCR amplification–induced mutations. A region of Chr8q24 previously shown to have no activity in any of the cell lines served as the negative control, and a region of Chr8q24 that was previously shown to have enhancer activity in all of the cell lines served as the positive control. For enhancer assays, SW480 (10 × 104 cells/well), HCT116, and RKO cells (6 × 104 cells/well) were seeded into 96-well plates. Cells were cotransfected with reporter plasmids and constitutively active pRL-TK Renilla luciferase Plasmid (Promega) using Lipofectamine 2000 Reagent (Life Technologies) according to the manufacturer's instructions. After 24 hours, cells were harvested, and extracts were assayed for luciferase activity using the Dual-Luciferase Reporter Assay System (Promega) according to the manufacturer's instructions and measured using a Tecan Infinite F200Pro microplate reader. The ratio of normalized luminescence from the experimental sample to the negative control reporter was calculated for each sample and defined as the relative luciferase activity. Luciferase activity was tested in a minimum of three independent clones for each allele. The data are presented as mean ± SD of four independent transfection experiments. Statistical differences in activity between the two alleles were calculated using a two-sided Student paired t test.

Gene editing by CRISPR

To delete an intronic fragment surrounding rs11648673 in AXIN1, RNA sequences were designed using Benchling (oligonucleotide sequences in Supplementary Table S1). Complementary oligonucleotides for each target sequence were incubated at 37°C for 30 minutes, followed by 95°C for 5 minutes, and annealed by decreasing 5°C/minute to 25°C using a PCR Thermocycler (Eppendorf). The annealed fragments were then ligated into the Cas9-T2A-GFP-Puro plasmid linearized by the BsaI restriction enzyme site. The plasmids were transformed into DH5α cells and were grown on ampicillin plates at 37°C overnight. Bacterial colonies were selected and expanded. Plasmid DNA was extracted using the HiSpeed Plasmid Maxi Kit (Qiagen). The guide RNA (gRNA) inserts were verified through Sanger sequencing. HCT-116 cells were transfected with Cas-T2A-GFP-Puro plasmid containing gRNA insert of interest using Lipofectamine 2000 to generate the AXIN1-knockout cell lines. HCT-116 cells were electroporated with Cas-T2A-GFP-Puro plasmid containing gRNA insert of interest and donor oligo using the Nucleofector (Lonza) to generate a deletion spanning the AXIN1 intronic SNP, rs11648673. Cells were grown in a T25 flask for 48 hours prior to single-cell sorting by dilution in 96-well plates. Cells were expanded and deletion spanning rs11648673 was verified through PCR and Sanger sequencing. The AXIN1 knockouts (Supplementary Fig. S3) were verified by complete abrogation of protein through Western blot analysis, while AXIN1 with 32-bp intronic deletion (Supplementary Fig. S3) preserved protein expression (results not shown).

TCF/LEF reporter assay

To test the effect of the AXIN1-mutant cell lines on WNT signaling activity, HCT-116 cells were plated at the density of 10,000 cells/well in a 96-well plate 5 hours prior to infection of 10 μL of Cignal Lenti TCF/LEF Reporter (Qiagen). After 24 hours at 37°C, the media were replaced and a final concentration of 10 μmol/L CT99021 was added to the cells. A Dual-Luciferase Reporter Assay (Promega) was performed after 24 hours and luminescence was read on a GloMax (Promega). Firefly luciferase was normalized to Renilla luminescence. The data are presented as mean ± SD of three independent experiments. Statistical differences between constructs were calculated using a two-sided Student paired t test.

Genomic analysis of OS

The top three ranked SNPs according to their corresponding unadjusted P values for association with OS were rs7300446, rs12426370, and rs10846370. They are intergenic between MGST1 and LMO3 and belong to a haplotype block (r2 > 0.85 among them). The next highly ranked SNP was rs11644916, located in intron 4 of AXIN1 (Table 2; for a LocusZoom plot of these variants, see Supplementary Fig. S4). Among these three genes, that is, MGST1, LMO3, and AXIN1, we selected AXIN1 for further evaluations and mechanistic studies, as the biology of AXIN1 is directly involved in the pathogenesis of colorectal cancer and progression of disease (13–15). AXIN1 is a cytoplasmic protein that negatively regulates the WNT signaling pathway by forming a β-catenin destruction complex with APC and GSK-3β. The WNT signaling pathway is highly dysregulated in colorectal cancer, making the SNP in AXIN1 (rs11644916, G>A) a relevant candidate to investigate further.

Table 2.

Top 10 SNPs ranked by P value for association with OS in CALGB/SWOG 80405 of the primary cohort.

RSIDChrGeneAllelesFeatureMAFHR95% LCL95% UCLP
rs7300446 12 MGST1/LMO3 A>G Intergenic 0.15 1.56 1.30 1.87 1.30 × 10−6 
rs12426370 12 MGST1/LMO3 A>G Intergenic 0.15 1.55 1.29 1.86 1.76 × 10−6 
rs10846370 12 MGST1/LMO3 A>G Intergenic 0.15 1.50 1.27 1.78 2.29 × 10−6 
rs11644916 16 AXIN1 G>A Intronic 0.30 1.40 1.20 1.59 4.26 × 10−6 
rs2461035 STMN2/HEY1 A>G Intergenic 0.33 0.72 0.62 0.83 4.65 × 10−6 
rs16985284 KCNS3/RDH14 A>C Intergenic 0.12 1.56 1.29 1.89 5.12 × 10−6 
rs11246159 11 ANO9 A>G Intronic 0.17 1.47 1.24 1.75 6.38 × 10−6 
rs11649255 16 AXIN1 A>G Intronic 0.29 1.38 1.20 1.58 6.83 × 10−6 
rs898838 MSGN1/KCNS3 A>G Intergenic 0.45 0.74 0.65 0.85 7.51 × 10−6 
rs2467754 STMN2/HEY1 A>G Intergenic 0.33 0.72 0.63 0.84 8.42 × 10−6 
RSIDChrGeneAllelesFeatureMAFHR95% LCL95% UCLP
rs7300446 12 MGST1/LMO3 A>G Intergenic 0.15 1.56 1.30 1.87 1.30 × 10−6 
rs12426370 12 MGST1/LMO3 A>G Intergenic 0.15 1.55 1.29 1.86 1.76 × 10−6 
rs10846370 12 MGST1/LMO3 A>G Intergenic 0.15 1.50 1.27 1.78 2.29 × 10−6 
rs11644916 16 AXIN1 G>A Intronic 0.30 1.40 1.20 1.59 4.26 × 10−6 
rs2461035 STMN2/HEY1 A>G Intergenic 0.33 0.72 0.62 0.83 4.65 × 10−6 
rs16985284 KCNS3/RDH14 A>C Intergenic 0.12 1.56 1.29 1.89 5.12 × 10−6 
rs11246159 11 ANO9 A>G Intronic 0.17 1.47 1.24 1.75 6.38 × 10−6 
rs11649255 16 AXIN1 A>G Intronic 0.29 1.38 1.20 1.58 6.83 × 10−6 
rs898838 MSGN1/KCNS3 A>G Intergenic 0.45 0.74 0.65 0.85 7.51 × 10−6 
rs2467754 STMN2/HEY1 A>G Intergenic 0.33 0.72 0.63 0.84 8.42 × 10−6 

Abbreviations: Chr, chromosome; LCL, lower confidence limit; UCL, upper confidence limit.

Association between AXIN1 rs11644916 and OS in CALGB/SWOG 80405

The A allele of rs11644916 (G>A) was associated with shorter OS with both treatment arms combined (P = 4.26 × 10−6; HR, 1.39; 95% CI 1.20–1.59; Fig. 1A). The median OS for patients with the AA, AG, or GG genotypes was 18.4 (95% CI, 14.2–27.6), 25.6 (95% CI, 23.6–30.4), or 36.4 (95% CI, 32.9–40.6) months, respectively.

Figure 1.

Kaplan–Meier plots of OS for rs11644916 (G>A) in AXIN1 in CALGB/SWOG 80405 (A) and for rs11644916 (A>G) in AXIN1 in TCGA (B).

Figure 1.

Kaplan–Meier plots of OS for rs11644916 (G>A) in AXIN1 in CALGB/SWOG 80405 (A) and for rs11644916 (A>G) in AXIN1 in TCGA (B).

Close modal

In a model adjusted for covariates and stratified by chemotherapy, there was no evidence of an interaction between rs11644916 and treatment (bevacizumab vs. cetuximab) in the association with OS (P = 0.23). HRs estimated separately in the bevacizumab and cetuximab arms were 1.52 (95% CI, 1.24–1.85) and 1.16 (95% CI, 0.94–1.44), respectively (Fig. 2). For the PFS outcome, there was no evidence of interaction between the SNP and treatment arm (bevacizumab vs. cetuximab; P = 0.87). Similarly, no evidence of interaction between the SNP and chemotherapy (FOLFOX vs. FOLFIRI) was observed (P = 0.68).

Figure 2.

Kaplan–Meier plots of OS for rs11644916 by arm (bevacizumab vs. cetuximab) in CALGB/SWOG 80405.

Figure 2.

Kaplan–Meier plots of OS for rs11644916 by arm (bevacizumab vs. cetuximab) in CALGB/SWOG 80405.

Close modal

Association between AXIN1 rs11644916 and OS in TCGA

In TCGA, because of the difference in the genotyping platform used, AXIN1 rs11649255 (A>G) was used for analysis with OS instead of rs11644916 (G>A). The two variants were almost in complete LD (r2 = 0.99). In 90 patients with stage IV colorectal cancer from TCGA, the AXIN1 rs11649255 (A>G) was associated with shorter survival (HR, 2.24; 95% CI, 1.20–4.19, P = 0.0096; one-degree of freedom; Fig. 1B). The direction of the observed effect was concordant with that of rs11644916 (G>A) in CALGB/SWOG 80405 (P = 4.26 × 10−6; HR = 1.39; 95% CI, 1.20–1.59).

The median OS for patients with the AA, AG, and GG genotypes was 44.4 (95% CI, 36.5–not yet reached), 22.7 (95% CI, 13.4–44.4), and 9.7 (95% CI, 5.1–not yet reached) months. Given that only 3 patients had GG genotype, a resampling permutation test was used as a sanity check for the asymptotic P value (P = 0.011 of the permutation test). For a direct comparison, the association of AXIN1 rs11649255 with OS was tested in CALGB/SWOG 80405, with similar results (P = 6.83 × 10−6; HR, 1.38; 95% CI, 1.20–1.58; Table 2).

Although we selected the AXIN1 variant for replication, we also evaluated the association of rs10772941 with OS as a proxy for rs10846370 (intergenic between MGST1 and LMO3, r2 = 0.87; Table 2). In the same 90 patients with stage IV colorectal cancer (which included only two variant homozygotes and 23 heterozygotes), there was no evidence that rs10772941 was associated with OS (P = 0.9802; HR, 1.01; 95% CI, 0.54–1.87).

Tumor location analysis

As expected, patients with left-sided tumors were associated with longer OS than right/transverse-sided tumors (P = 1.37 × 10−7; HR, 0.59; 95% CI, 0.49–0.72). Both AXIN1 rs11644916 (associated with OS) and rs11648673 (functional variant selected by bioinformatics) were associated with an increased odds of having the tumor on the right/transverse side of the colon (P = 1.85 × 10−3 and 2.21 × 10−3, respectively; Fig. 3A and B). Both SNPs were associated with worse OS in a model adjusted for age, sex, and tumor location (P = 1.80 × 10−4; HR, 1.32; 95% CI, 1.14–1.53 and P = 3.28 × 10−3; HR, 1.25; 95% CI, 1.08–1.46, respectively; Fig. 3C and D).

Figure 3.

Tumor location by genotype for rs11644916 (A) and rs11648673 (B) in AXIN1. Kaplan–Meier plots of OS by tumor location and rs11644916 (C) and rs11648673 (D).

Figure 3.

Tumor location by genotype for rs11644916 (A) and rs11648673 (B) in AXIN1. Kaplan–Meier plots of OS by tumor location and rs11644916 (C) and rs11648673 (D).

Close modal

Supplementary Table S2 shows the results of the top 10 SNPs, ranked by P value from the unadjusted analysis of OS in CALGB/SWOG 80405, adjusted for age (log10 transformed), sex, and tumor location (left vs. right/transverse). Supplementary Table S3 shows the results for the same SNPs, as Supplementary Table S2, when tumor location was replaced in the model by BRAF V600E, MSI, and CMS, for the 289 patients for which those data were available.

Bioinformatics assessment of SNPs associated with OS

For the three intergenic SNPs between MGST1 and LMO3, none of them appeared to be located in regulatory regions or have enhancer functions (Supplementary Table S4). Similarly, AXIN1 rs11644916, located in intron 4, did not appear to be located in a regulatory region or affect transcription factor binding. The analysis of HaploReg for other SNPs in high LD with AXIN1 rs11644916 selected rs11648673 (in intron 2 of AXIN1) as the functional SNP (LD r2 = 0.81). This SNP is located in a region containing several regulatory elements, and RegulomeDB indicated that it is more likely to affect transcription factor binding (score = 1f; Supplementary Table S4). The A allele of rs11648673 was associated with worse OS (P = 1.86 × 10−4; HR, 1.31; 95% CI, 1.14–1.52; Supplementary Fig. S5). Because of the lack of bioinformatics evidence for the three SNPs located between MGST1 and LMO3 (Supplementary Table S4), we performed functional assays on the AXIN1 SNP.

Luciferase assay of AXIN1 rs11648673

Because of the lack of functionality of rs11644916 through bioinformatics, further functional validation was performed using the SNP in high LD with it, rs11648673 (Supplementary Table S4). A luciferase reporter assay was used to determine differences in enhancer activity between the reference and variant allele of rs11648673 in AXIN1. The forward constructs did not show any enhancer activity. The reverse constructs showed enhancer activity in both alleles, with the G allele showing higher activity than the A allele in all three colorectal cancer cell lines tested (Fig. 4A).

Figure 4.

A, Luciferase assay of AXIN1 rs11648673. The G allele had statistically greater enhancer activity than the A allele of rs11648673 in SW480 (P = 0.0008), HCT-116 (P = 0.0042), and RKO (P = 0.0016) human colorectal cancer cell lines. RLU, relative luminometer units. B, TCF/LEF reporter assay. Cells were treated with 10 μmol/L CT99021 for WNT activation. Luciferase activity for each mutated HCT-116 cell line +/− CT99021 was normalized to WT HCT-116 (control) cell line at baseline (shown with an *). The log2 (fold change) is reported. N = 3 independent replicates.

Figure 4.

A, Luciferase assay of AXIN1 rs11648673. The G allele had statistically greater enhancer activity than the A allele of rs11648673 in SW480 (P = 0.0008), HCT-116 (P = 0.0042), and RKO (P = 0.0016) human colorectal cancer cell lines. RLU, relative luminometer units. B, TCF/LEF reporter assay. Cells were treated with 10 μmol/L CT99021 for WNT activation. Luciferase activity for each mutated HCT-116 cell line +/− CT99021 was normalized to WT HCT-116 (control) cell line at baseline (shown with an *). The log2 (fold change) is reported. N = 3 independent replicates.

Close modal

TCF/LEF reporter assay

We generated two AXIN1-knockout HCT-116 cell lines through the creation of frameshift mutations in exon 2, as well as a mutant HCT-116 cell line with a 32-bp deletion surrounding rs11648673 in intron 2. A TCF/LEF reporter assay measured β-catenin binding to its transcriptional domain, allowing for visualization of the effect of this deletion on WNT/β-catenin signaling. There was no change in TCF/LEF reporter activity in the AXIN1-knockout cell lines compared with WT. In the mutant cell line containing the 32-bp deletion spanning rs11648673, there was a significant increase in luciferase activity from the WNT reporter compared with the WT cell line, with and without treatment of a WNT activator (CT99021; Fig. 4B).

A genome-wide association study (GWAS) of OS in patients with mCRC treated with standard-of-care therapies has selected a common SNP in AXIN1 as a novel prognostic factor in mCRC. This result has been supported by a concordant effect on OS in patients with mCRC from TCGA, as well as bioinformatics and experimental evidence of functionality of the SNP in colorectal cancer cell lines. AXIN1 is a known tumor suppressor gene in the WNT pathway, but this is the first germline variation in this gene to be reported as associated with reduced survival in mCRC.

AXIN1 is negative regulator of the WNT/β-catenin pathway. The WNT/β-catenin pathway is central to colorectal cancer tumorigenesis and progression (13–15) and is highly dysregulated in more than 90% of colorectal cancers (16). In the WNT pathway, AXIN1 acts with APC, GSK3β, and CK1 to form a β-catenin destruction complex, which prevents transcriptional activation of targets of β-catenin, the ultimate effector of the pathway (13). In cell lines and transgenic mice, overexpression of AXIN1 leads to increased β‐catenin degradation and attenuation of WNT signaling, supporting its tumor suppressor activity (17, 18).

The luciferase assays in three human colorectal cancer cell lines clearly indicate a reduction in the enhancer activity in the variant construct versus WT (Fig. 4A). This result is consistent with the biological function of AXIN1 because a potential reduction in the ability of AXIN1 to efficiently activate a destruction complex would result in a more aggressive tumor phenotype, and eventually, worse survival. Moreover, despite the fact that WNT is constitutively activated in the HCT-116 colorectal cancer cell line due to a deletion in the β-catenin gene (19), the deletion of a 32-bp region surrounding the AXIN1 SNP in HCT-116 colorectal cancer cells increases WNT activation (measured by the TCF/LEF reporter activity assay), both at baseline and after stimulation with a WNT activator (Fig. 4B). These results indicate that the AXIN1 variant might interact with a functional element in intron 2 of AXIN1. The bioinformatic analysis (Supplementary Table S4) already indicates the presence of a few regulatory transcription factors, and we hypothesize that their activity is modulated by the presence of the AXIN1 variant.

While the experimental evidence provided above is in-line with the clinical effect on OS, it does not exactly recapitulate the biological mechanism by which a SNP in AXIN1 alters OS. Our study provides the rationale for developing animal models to test the underlying mechanism. At the same time, it should be considered that alterations of the WNT pathway in cancer cells might alter the function of other cell systems. In addition to its central role in colorectal cancer oncogenesis and progression, WNT/β-catenin was the first pathway to be causally linked to T-cell inflammation (20). WNT/β-catenin signaling promotes cancer progression by regulating the tumor immune cycle in most of its nodes, including dendritic, T, and tumor cells (21). Recent publications in colorectal cancer have shown that the activation of tumor intrinsic WNT/β-catenin signaling was correlated with the absence of T-cell infiltration in colorectal cancer (22). This mechanism is a pan-tumor effect observed in 90% of tumor types from TCGA, including colorectal cancer (23). Tumor, immunocompetent xenografts should be used for testing the function of genetic manipulations of the AXIN1 gene.

The potential tumor suppressor function of AXIN1 has long been hypothesized based largely on its role in the β-catenin destruction complex, but the in vivo importance of AXIN1 protein for the inhibition of tumor development or progression is not well established, and the functional consequences of AXIN1 mutations in cancer remain largely undefined (24).

This study provides further evidence of the tumor suppressor role of AXIN1, indicating that germline alterations in the AXIN1 gene affect the survival of patients with mCRC. With all the caveats of using survival data from TCGA, the effect of the AXIN1 SNP was concordant between TCGA and CALGB/SWOG 80405 (Fig. 1A and B), even after a resampling, permutation test in TCGA (due its smaller sample size). The results of this study suggest that this could be a prognostic effect, as it seems independent from the biologic used (bevacizumab or cetuximab). In this regard, an uncertainty remains in relation to establishing differences between FOLFIRI and FOLFOX, as only one-third of patients were treated with FOLFIRI. Cytotoxicity experiments with 5-FU, irinotecan, SN-38, or oxaliplatin in a colorectal cancer cell line rule out that this AXIN1 variant can change the efficacy of these agents in vitro (results not shown).

None of the SNP–OS associations had a P value less than 10−8, the conventional genome-wide statistical significance. Although the lack of adjustment for genome-wide testing can increase the likelihood of false positivity, as stated in the Materials and Methods section of this article, the P values were used for feature selection to direct further validation and experimental studies. When a GWAS is not preplanned and is not the primary endpoint of the clinical study, as in this case, losses in power and/or statistical significance should be compensated by triangulation of evidence from other results. This limitation is, therefore, in part, mitigated by the well-established role of AXIN1 in colorectal cancer, the replication of this signal in TCGA, and the functional assessment of the effect of the AXIN1 variants in colorectal cancer cell lines.

Interestingly, the AXIN1 variant seems to colocalize, in part, with the right/transverse primary tumors, which are characterized by worse OS as compared with left tumors (7). However, both the AXIN1 variant and right/transverse primary tumors contribute, independently, to worse OS (Fig. 3C and D). If further confirmed in additional studies, AXIN1 germline variants could represent a modifier predisposing patients with colorectal cancer to have right/transverse-sided and more aggressive disease at presentation.

A.S. Etheridge reports grants from Alliance for Clinical Trials during the conduct of the study. F.-S. Ou reports grants from NCI during the conduct of the study as well as outside the submitted work. S. J. Plummer reports grants from NIH/NCI during the conduct of the study and NIH outside the submitted work. G. Casey reports grants from NCI during the conduct of the study as well as outside the submitted work. H.L. McLeod reports other from Clariifi (cofounder) and Interpares Biomedicine (cofounder) during the conduct of the study, personal fees from Pharmazam outside the submitted work, and is a visiting professor at the Clinical Pharmacology Institute at Xiangya Hospital (Changsha China), which supports travel and accommodations for in-person sessions with students and faculty. C.D. Blanke reports grants from NCI during the conduct of the study. A.P. Venook reports grants from NIH during the conduct of the study as well as other from Genentech/Roche (research support) and Bristol-Myers Squibb (research support) outside the submitted work. H.-J. Lenz reports personal fess from Bristol-Myers Squibb (ad board), Genentech Roche (ad board), Bayer (ad board lectures), and Merck KG (ad board lectures) during the conduct of the study. M.J. Ratain reports personal fees from multiple generic companies (patent litigation, including joint defense groups), Pneuma Respiratory, Aptevo, and Cyclacel, grants and personal fees from Genentech, grants from Bristol-Myers Squibb and Xencor, and other from BeiGene (meeting support) outside the submitted work, as well as listed as a coinventor on a patent regarding portfolio regarding UGT1A1 genotyping and irinotecan treatment owned by the University of Chicago and licensed to Mayo Medical, a pending patent regarding genomic prescribing system owned by the University of Chicago, and a pending patent for a low-dose tocilizumab for COVID owned by the University of Chicago, and is director and treasurer of the Value in Cancer Care Consortium. No disclosures were reported by the other authors.

The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institutes of Health.

F. Innocenti: Conceptualization, resources, formal analysis, supervision, funding acquisition, methodology, writing-original draft, project administration, writing-review and editing. A.B. Sibley: Data curation, formal analysis, visualization, methodology, writing-review and editing. S.A. Patil: Formal analysis, writing-review and editing. A.S. Etheridge: Visualization, methodology, writing-review and editing. C. Jiang: Data curation, formal analysis, visualization, writing-review and editing. F.-S. Ou: Resources, data curation, writing-review and editing. S.D. Howell: Formal analysis, investigation, methodology, writing-review and editing. S.J. Plummer: Investigation, visualization, methodology, writing-review and editing. G. Casey: Resources, supervision, investigation, writing-review and editing. M.M. Bertagnolli: Resources, funding acquisition, project administration, writing-review and editing. H.L. McLeod: Resources, writing-review and editing. J.T. Auman: Investigation, writing-review and editing. C.D. Blanke: Resources, project administration, writing-review and editing. Y. Furukawa: Resources, writing-review and editing. A.P. Venook: Resources, writing-review and editing. M. Kubo: Resources, writing-review and editing. H.-J. Lenz: Resources, writing-review and editing. J.S. Parker: Supervision, writing-review and editing. M.J. Ratain: Resources, project administration, writing-review and editing. K. Owzar: Data curation, formal analysis, supervision, investigation, writing-review and editing.

Research reported in this article was supported by the NCI of the NIH under award numbers U10CA180821, U10CA180882, and U24CA196171 (to the Alliance for Clinical Trials in Oncology), UG1CA233253, U10CA180820 (ECOG-ACRIN), and U10CA180888 (SWOG) and NIH R01CA143237. https://acknowledgments.alliancefound.org.

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.

1.
Venook
AP
,
Niedzwiecki
D
,
Lenz
HJ
,
Innocenti
F
,
Fruth
B
,
Meyerhardt
JA
, et al
Effect of first-line chemotherapy combined with cetuximab or bevacizumab on overall survival in patients with KRAS wild-type advanced or metastatic colorectal cancer: a randomized clinical trial
.
JAMA
2017
;
317
:
2392
401
.
2.
Innocenti
F
,
Owzar
K
,
Cox
NL
,
Evans
P
,
Kubo
M
,
Zembutsu
H
, et al
A genome-wide association study of overall survival in pancreatic cancer patients treated with gemcitabine in CALGB 80303
.
Clin Cancer Res
2012
;
18
:
577
84
.
3.
Hanahan
D
,
Weinberg
RA
. 
Hallmarks of cancer: the next generation
.
Cell
2011
;
144
:
646
74
.
4.
Price
AL
,
Patterson
NJ
,
Plenge
RM
,
Weinblatt
ME
,
Shadick
NA
,
Reich
D
. 
Principal components analysis corrects for stratification in genome-wide association studies
.
Nat Genet
2006
;
38
:
904
9
.
5.
Ullah
AZD
,
Oscanoa
J
,
Wang
J
,
Nagano
A
,
Lemoine
NR
,
Chelala
C
. 
SNPnexus: assessing the functional relevance of genetic variation to facilitate the promise of precision medicine
.
Nucleic Acids Res
2018
;
46
:
W109
13
.
6.
Therneau
TM
,
Grambsch
PM
.
Modeling survival data: extending the Cox model
.
New York, NY
:
Springer
; 
2000
.
7.
Innocenti
F
,
Ou
FS
,
Qu
X
,
Zemla
TJ
,
Niedzwiecki
D
,
Tam
R
, et al
Mutational analysis of patients with colorectal cancer in CALGB/SWOG 80405 identifies new roles of microsatellite instability and tumor mutational burden for patient outcome
.
J Clin Oncol
2019
;
37
:
1217
27
.
8.
Pruim
RJ
,
Welch
RP
,
Sanna
S
,
Teslovich
TM
,
Chines
PS
,
Gliedt
TP
, et al
LocusZoom: regional visualization of genome-wide association scan results
.
Bioinformatics
2010
;
26
:
2336
7
.
9.
R Core Team
.
R: a language and environment for statistical computing
.
Vienna, Austria
:
R Foundation for Statistical Computing
; 
2020.
10.
Therneau
T
. 
A package for survival analysis in S. version 2.38
.
Available from:
https://CRAN.R-project.org/package=survival.
11.
Xie
Y
.
Dynamic Documents with R and knitr
. 2nd ed.
Boca Raton, FL
:
Chapman and Hall/CRC
; 
2015
.
12.
The Cancer Genome Atlas Network
. 
Comprehensive molecular characterization of human colon and rectal cancer
.
Nature
2012
;
487
:
330
7
.
13.
Li
VS
,
Ng
SS
,
Boersema
PJ
,
Low
TY
,
Karthaus
WR
,
Gerlach
JP
, et al
Wnt signaling through inhibition of beta-catenin degradation in an intact Axin1 complex
.
Cell
2012
;
149
:
1245
56
.
14.
Flanagan
DJ
,
Vincan
E
,
Phesse
TJ
. 
Wnt signaling in cancer: not a binary ON:OFF switch
.
Cancer Res
2019
;
79
:
5901
6
.
15.
Zhan
T
,
Rindtorff
N
,
Boutros
M
. 
Wnt signaling in cancer
.
Oncogene
2017
;
36
:
1461
73
.
16.
Silva
AL
,
Dawson
SN
,
Arends
MJ
,
Guttula
K
,
Hall
N
,
Cameron
EA
, et al
Boosting Wnt activity during colorectal cancer progression through selective hypermethylation of Wnt signaling antagonists
.
BMC Cancer
2014
;
14
:
891
.
17.
Picco
G
,
Petti
C
,
Centonze
A
,
Torchiaro
E
,
Crisafulli
G
,
Novara
L
, et al
Loss of AXIN1 drives acquired resistance to WNT pathway blockade in colorectal cancer cells carrying RSPO3 fusions
.
EMBO Mol Med
2017
;
9
:
293
303
.
18.
Hsu
W
,
Shakya
R
,
Costantini
F
. 
Impaired mammary gland and lymphoid development caused by inducible expression of Axin in transgenic mice
.
J Cell Biol
2001
;
155
:
1055
64
.
19.
Miyamoto
M
,
Hayashi
T
,
Kawasaki
Y
,
Akiyama
T
. 
Sp5 negatively regulates the proliferation of HCT116 cells by upregulating the transcription of p27
.
Oncol Lett 2018
:
15
:
4005
9
.
20.
Spranger
S
,
Bao
R
,
Gajewski
TF
. 
Melanoma-intrinsic beta-catenin signalling prevents anti-tumour immunity
.
Nature
2015
;
523
:
231
5
.
21.
Wang
B
,
Tian
T
,
Kalland
KH
,
Ke
X
,
Qu
Y
. 
Targeting Wnt/beta-catenin signaling for cancer immunotherapy
.
Trends Pharmacol Sci
2018
;
39
:
648
58
.
22.
Grasso
CS
,
Giannakis
M
,
Wells
DK
,
Hamada
T
,
Mu
XJ
,
Quist
M
, et al
Genetic mechanisms of immune evasion in colorectal cancer
.
Cancer Discov
2018
;
8
:
730
49
.
23.
Luke
JJ
,
Bao
R
,
Sweis
RF
,
Spranger
S
,
Gajewski
TF
. 
WNT/beta-catenin pathway activation correlates with immune exclusion across human cancers
.
Clin Cancer Res
2019
;
25
:
3074
83
.
24.
Mazzoni
SM
,
Fearon
ER
. 
AXIN1 and AXIN2 variants in gastrointestinal cancers
.
Cancer Lett
2014
;
355
:
1
8
.