Skip to main content
  • Primary research
  • Open access
  • Published:

Identification of key genes involved in tumor immune cell infiltration and cetuximab resistance in colorectal cancer

Abstract

Background

The anti-epidermal growth factor receptor (EGFR) antibody introduces adaptable variations to the transcriptome and triggers tumor immune infiltration, resulting in colorectal cancer (CRC) treatment resistance. We intended to identify genes that play essential roles in cetuximab resistance and tumor immune cell infiltration.

Methods

A cetuximab-resistant CACO2 cellular model was established, and its transcriptome variations were detected by microarray. Meanwhile, public data from the Gene Expression Omnibus and The Cancer Genome Atlas (TCGA) database were downloaded. Integrated bioinformatics analysis was applied to detect differentially expressed genes (DEGs) between the cetuximab-resistant and the cetuximab-sensitive groups. Then, we investigated correlations between DEGs and immune cell infiltration. The DEGs from bioinformatics analysis were further validated in vitro and in clinical samples.

Results

We identified 732 upregulated and 1259 downregulated DEGs in the induced cellular model. Gene Ontology and Kyoto Encyclopedia of Genes and Genomes pathway enrichment analyses, along with Gene Set Enrichment Analysis and Gene Set Variation Analysis, indicated the functions of the DEGs. Together with GSE59857 and GSE5841, 12 common DEGs (SATB-2, AKR1B10, ADH1A, ADH1C, MYB, ATP10B, CDX-2, FAR2, EPHB2, SLC26A3, ORP-1, VAV3) were identified and their predictive values of cetuximab treatment were validated in GSE56386. In online Genomics of Drug Sensitivity in Cancer (GDSC) database, nine of twelve DEGs were recognized in the protein-protein (PPI) network. Based on the transcriptome profiles of CRC samples in TCGA and using Tumor Immune Estimation Resource Version 2.0, we bioinformatically determined that SATB-2, ORP-1, MYB, and CDX-2 expressions were associated with intensive infiltration of B cell, CD4+ T cell, CD8+ T cell and macrophage, which was then validated the correlation in clinical samples by immunohistochemistry. We found that SATB-2, ORP-1, MYB, and CDX-2 were downregulated in vitro with cetuximab treatment. Clinically, patients with advanced CRC and high ORP-1 expression exhibited a longer progression-free survival time when they were treated with anti-EGFR therapy than those with low ORP-1 expression.

Conclusions

SATB-2, ORP-1, MYB, and CDX-2 were related to cetuximab sensitivity as well as enhanced tumor immune cell infiltration in patients with CRC.

Background

In the United States, colorectal cancer (CRC) is estimated to be the third-most frequently occurring malignant tumor and third leading cause of cancer-related mortality, leading to the deaths of about 53,200 people in 2020 [1]. The incidence and mortality rates in 2015 were 37,000 cases and 19,000 deaths, with a 25- to 28-months overall survival time for Chinese patients suffering advanced CRC [2]. Targeting epidermal growth factor receptor (EGFR), monoclonal antibodies cetuximab and panitumumab are effective therapeutics [3, 4] that shrink tumors for resectability [5,6,7] or relieve symptoms from unresectable masses by repressing tumor growth [8,9,10].

Mutations in the RAS family suggest continuous activation of EGFR downstream signaling [11,12,13,14,15] and RAS genes (KRAS, NRAS) are used as biomarkers to predict poor response to anti-EGFR therapy [16, 17]. However, tumor microenvironment re-plasticity, characterized by tumor immune cell infiltration and upregulation of immune checkpoint-related proteins, has recently been recognized as a novel explanation for drug resistance and treatment failure [18,19,20]. In particular, the interaction between tumor and immune microenvironment triggered by cetuximab or panitumumab in CRC accounts for immunogenic cell death and immune treatment resistance [21, 22].

In this study, we firstly performed high-throughput screening of transcriptomic alterations before and after induction cetuximab administration in a colon cancer cell model (CACO2). Then, we conducted integrated bioinformatics analysis of the transcriptome variations using public datasets from the Gene Expression Omnibus (GEO) database to identify several core differentially expressed genes (DEGs) and investigated whether they were correlated with tumor immune cell infiltration. Finally, we validated the identified core DEGs in vitro and in clinical samples.

Methods

Cellular model and culture

Cell lines of human CRC (CACO2, HCT116, HT29, NCIH508, and RKO) were obtained from the Cell Bank of the Chinese Academy of Sciences (Shanghai, China). All cell lines were authenticated by short tandem repeat (STR) profiling (Genetic Testing Biotechnology Corporation, Suzhou, China) and routinely tested for mycoplasma using MycoAlertā„¢ Mycoplasma Detection Kit (Lonza; LT07-218, Rockland, ME, USA). We established an EGFR antagonist-resistant cellular model using cetuximab sensitive cell line (CACO2-CS) and validated the transcriptional changes induced by cetuximab in the other cell lines. Using stepwise induction, we started with a low dose causing 50% cell growth inhibition (IC50). The dose was increased to 10Ā Āµg/mL after about two months, 50Ā Āµg/mL after another two months, and lastly 300Ā Āµg/mL. Cell Counting Kit 8 (CCK8) (Dojindo, Shanghai, China) was used to test the cell viability. For colonies formation assay, cells were cultured in six-well plates with 1000 cells and 2Ā ml media including cetuximab (300Ā Āµg/ml) per well. After 14 days, colonies were fixed with 4% paraformaldehyde and stained with 0.1% crystal violet. Finally, the established cetuximab resistant cell line (CACO2-CR) was cultured with the maximal dose of cetuximab. CACO2-CS, HCT116, HT29, NCIH508, and RKO in the exponential growth phase were precultured in 12-well tissue culture plates for 24Ā h. Different concentrations of cetuximab (S20130004, Merck KGaA, Darmstadt, Germany) were added to cells which were incubated for 72Ā h (0Ā Āµg/mL and 25Ā Āµg/mL for CACO2, 0Ā Āµg/mL and 1Ā Āµg/mL for NCIH508, and 0Ā Āµg/mL and 50Ā Āµg/mL for HT29, HCT116, and RKO). CACO2-CS and CACO2-CR cells were cultured in Dulbeccoā€™s Modified Eagleā€™s Medium (DMEM) (HyClone, Logan, UT, USA) containing 20% fetal bovine serum (FBS) (Gibco, Paisley, UK) 100 U/mL penicillin, and 100 U/mL streptomycin (Gibco, Paisley, UK). HT-29 cells were cultured in McCoyā€™s 5A (Gibco, Grand Island, NY, USA) supplemented with 10% FBS (Gibco, Paisley, UK), 100 U/mL penicillin, and 100 U/mL streptomycin (Gibco, Paisley, UK). HCT116 and NCIH508 cells were cultured in Roswell Park Memorial Institute Modified Medium (Hyclone) supplemented with 10% FBS (Gibco, Paisley, UK), 100 U/mL penicillin, and 100 U/mL streptomycin (Gibco, Paisley, USA). RKO cells were cultured in high glucose DMEM (Hyclone) supplemented with 10ā€‰% FBS (Gibco, Paisley, UK), 100 U/mL penicillin, and 100 U/mL streptomycin (Gibco, Paisley, UK). All cells were cultured in a humidified 5% CO2 incubator at 37Ā Ā°C. All cell lines experiments were repeated at least three times with three to six replicates.

Microarray screening for transcriptional variation after cetuximab induction

Total RNA was extracted and purified with the miRNeasy Mini Kit (Cat. # 217,004 QIAGEN GmBH, Hilden, Germany) according to the manufacturerā€™s instructions. The RNA integrity number was determined with the Agilent Bioanalyzer 2100 (Agilent Technologies, Santa Clara, CA, USA), and the RNA was amplified and labeled with the Low Input Quick Amp WT Labeling Kit (Cat. # 5190āˆ’2943, Agilent Technologies, Santa Clara, CA, USA). Labeled cRNA was purified with the RNeasy mini kit (Cat. # 74,106). Samples were hybridized with about 1.65Ā Āµg Cy3-labeled cRNA by the Gene Expression Hybridization Kit (Cat. # 5188ā€“5242, Agilent Technologies) in a hybridization oven (Cat. # G2545A, Agilent Technologies). After about 17 h, the slides were washed in staining dishes (Cat. # 121, Thermo Shandon, Waltham, MA, USA) with the Gene Expression Wash Buffer Kit (Cat. # 5188ā€“5327, Agilent Technologies) and scanned by the Agilent Microarray Scanner (Cat. # G2565CA, Agilent Technologies). Data were retrieved with Feature Extraction software 10.7 (Agilent Technologies). Quantile normalization of the raw data was achieved by the limma package in R (version 3.4.1) [23].

Real-time quantitative PCR after cetuximab treatment in vitro

Total cellular RNA of cell lines treated with cetuximab was extracted with TRIzol reagent (Invitrogen, Carlsbad, CA, USA) under the manufacturerā€™s instructions. Recombinant DNase I (Takara Bio, Beijing, China) was used to remove potential genomic DNA contamination. cDNA was generated with the PrimerScriptā„¢ RT master mix kit (Takara Bio, Dalian China). Real-time quantitative PCR analysis of the identified DEGs was performed using the TB Green Premix Ex Taqā„¢ kit (Takara Bio, Dalian China). All results were normalized to human GAPDH mRNA expression. The primers were listed in Additional file 1: Table S1. The relative threshold cycle (Ct) method was used to display the results.

Identification of differentially expressed genes (DEGs) in the cellular model

For the microarray screening of CACO2 cells treated with cetuximab, we visualized the principal component analysis (PCA) results using the R package ggbiplot (https://github.com/vqv/ggbiplot) after converting the raw signals into gene expression levels. We set thresholds of pā€‰<ā€‰0.05 and | log2 (fold change) | > 1 to identify significant DEGs using the limma package and created heat maps using the heatmap package.

Gene expression omnibus (GEO) microarray data and DEGs analysis

To further narrow down candidate DEGs in CRC under EGFR antagonist pressure, we downloaded three gene expression datasets including cell lines and clinical tissue sample data from GEO (http://www.ncbi.nlm.nih.gov/geo). GSE59857 included the transcriptional and pharmacological profiles of 155 CRC cell lines. GSE5851 contained data from 80 clinical advanced CRC samples obtained before cetuximab monotherapy. GSE56386 contained data from eight primary CRC tumor tissue samples, comprising four from responders to cetuximab therapy and four from non-responders. From GSE59857, we downloaded the data for 20 cell lines, including 10 cetuximab-sensitive cases (OXCO2: GSM1448146, NCIH508: GSM1448142, DIFI: GSM1448175, COCM1: GSM1448167, CCK81: GSM1448097, C75: GSM1448201, HCA46: GSM1448177, C99: GSM1448204, HDC82: GSM1448128, and COGA1: GSM1448099) and 10 cetuximab-insensitive cases (SNU1047: GSM1448085, COLO320DM: GSM1448173, HUTU80: GSM1448180, KM12: GSM1448073, HDC143: GSM1448185, KM12SM: GSM1448188, KM12C: GSM1448186, COLO320: GSM1448152, KM12L4: GSM1448187, and C10: GSM1448196), among which HDC143 was deemed invalid and removed after analysis. We then converted the probe IDs into gene symbols using illuminaHumanv4.db. We retrieved the data for eight patients with the wild-type KRAS gene, including four responders to cetuximab treatment (GSM136609, GSM136593, GSM136654, and GSM136626) and four non-responders (GSM136635, GSM136646, GSM136607, and GSM136640) from GSE5851 and converted the probe IDs into gene symbols using the annotation library hgu133a2.db. The procedure for analyzing DEGs was the same as that for analyzing the CACO2 microarray data.

Gene Ontology (GO) term and Kyoto Encyclopedia of genes and genomes (KEGG) pathway enrichment analyses

For the enrichment analyses of the DEGs, a Metascape (https://metascape.org/) online analysis tool was used [24]. The parameters were set as the following: min overlapā€‰=ā€‰3, p-value cutoffā€‰=ā€‰0.01, and min enrichmentā€‰=ā€‰1.5. Related terms were selected from the top 20, according to their p-values.

Gene sets Enrichment analysis (GSEA) and Gene set variation analysis (GSVA)

GSEA was used to screen the biological states and processes associated with significantly upregulated or downregulated genes in the resistance group, with H.all.v7.1.symbols.gmt as the reference set [25]. The number of permutations was set to 1000. A p-valueā€‰<ā€‰0.05 was considered indicative of significant enrichment. A normalized enrichment score was established to evaluate the degree of enrichment. GSVA was used to describe the enrichment degree of a biological state or process in a sample [26]. Enrichment scores were compared, and those meeting the standards of |log2 (fold change) | > 1 and pā€‰<ā€‰0.05 were considered significant results.

Construction of proteinā€protein interaction (PPI) network

The PPIs of the central DEGs were generated using the STRING online database (https://string-db.org). To construct a PPI network, we retrieved cetuximab resistance candidates from the Genomics of Drug Sensitivity in Cancer (GDSC) database (https://www.cancerrxgene.org) [27], confirmed their interactions using the STRING database [28], and finally visualized the result using Cystoscape [29].

Prediction of tumor immune cell infiltration

To investigate the relationship between hub DEG expression and tumor immune cell infiltration, the transcriptome landscapes of colonic or rectal adenocarcinoma samples (nā€‰=ā€‰177) in The Cancer Genome Atlas (TCGA) database were analyzed using the online tool Tumor Immune Estimation Resource Version 2.0 (https://cistrome.shinyapps.io/timer/) [30, 31].

Validation in clinical samples by immunohistochemistry staining

Patients with advanced CRC (nā€‰=ā€‰102) were diagnosed by pathological examination and administered systemic chemotherapeutics at Zhongshan Hospital of Fudan University. Fifty-two patients received therapy containing cetuximab. The progression-free survival time (PFS) was recorded to evaluate the efficacy of cetuximab. Patients signed informed consent documents, and the ethics committee of Zhongshan Hospital approved the study.

Tumor samples were fixed with 4% paraformaldehyde, embedded in paraffin, cut into sections of about 5Ā Āµm, and placed onto glass slides. After the samples were deparaffinized with xylene, hydrophilized, and unmasked, they were blocked with bovine serum albumin, immunostained with primary antibodies against SATB-2 (21307-1-AP, 1:100, Proteintech Group, Inc., Wuhan, China), ORP-1 (bs-17514R, 1:200, Bioss Biological Technology Co., Ltd., Beijing, China), MYB (bs-5978R, 1:200, Bioss Biological Technology Co., Ltd.), CDX-2 (bsm-33063m, 1:200, Bioss Biological Technology Co., Ltd.), CD8 (GB13068, 1:100, Servicebio Technology, Wuhan, China), CD19 (GB11061, 1:500, Servicebio Technology), CD4 (GB13064, 1:100, Servicebio Technology), and CD68 (GB13067-M-2, 1:100, Servicebio Technology) in a humidified environment at about 4Ā Ā°C overnight and incubated with goat anti-rabbit or anti-mouse secondary antibody (1:200) for about 30Ā min at about 20Ā Ā°C. Subsequently, the slides were stained with 3,3ā€™-diaminobenzidine and counterstained with hematoxylin. Antigenā€“antibody complexes in the whole sample were detected using a panoramic slice scanner (3DHISTECH, Budapest, Hungary), recorded in a file, and viewed using CaseViewer 2.2 (3DHISTECH). To evaluate gene expression in tissues, the following formula was used to calculate the H-score using Quant Center 2.1 (3DHISTECH): H-SCORE = āˆ‘ (PI Ɨ I) = (percentage of cells of weak intensity Ɨ 1) + (percentage of cells of moderate intensity Ɨ 2)ā€‰+ā€‰percentage of cells of strong intensity Ɨ 3), where PI is the proportion of the positive signal pixel area and I is the coloring intensity.

Statistical analysis

Most statistical analyses were completed using bioinformatic tools mentioned above in R (version 3.4.1). The Benjamini and Hochberg False Discovery Rate method was utilized to adjust the p- values in screening DEGs from GEO profiles. Fisherā€™s exact test was employed to identify the significant GO terms and KEGG pathways. The correlation significance was examined by Spearman and Pearson correlation analyses. Differential expression levels of identified DEGs were assessed by a two-tailed Studentā€™s t-test. Kaplan - Meier survival curve and Log - rank test analysis was applied to investigated the predictive valued of identified DEGs for patients with CRC and receiving cetuximab contained therapy. The value of pā€‰<ā€‰0.050 was considered statistically significant.

Results

Model establishment and screening of DEGs

The study design is shown in Fig.Ā 1. After the cetuximab resistant cellular model was established successfully (Fig.Ā 2a), preprocessing of the raw data revealed a uniform distribution of DEGs between CACO2-CS and CACO2-CR (Fig.Ā 2b). The PCA results had acceptable reproducibility (Fig.Ā 2c). A total of 1991 DEGs were identified, with 732 upregulated genes and 1259 downregulated genes in CACO2-CR. These DEGs are shown in a volcano plot (Fig.Ā 2d) and a heat map (Fig.Ā 2e) and are listed in TableĀ 1. In GO analysis, the most enriched GO terms determined by the cellular component function (CC) were ā€œextracellular matrix,ā€ ā€œapical plasma membraneā€ and ā€œextracellular matrix componentā€ (Fig.Ā 2h). KEGG signaling analysis suggested that the DEGs were considerably enriched in ā€œprotein digestion and absorption,ā€ ā€œcytokine-cytokine receptor interaction,ā€ ā€œretinol metabolismā€ shown in Fig.Ā 2i and TableĀ 2.

Fig. 1
figure 1

Flow chart of the study

Fig. 2
figure 2

Identification and enrichment analysis of DEGs related to anti-EGFR antibody resistance in CACO2. a Establishment of a cetuximab-resistant CACO2 cell line (CACO2-CR) from its parental sensitive cell line (CACO2-CS) and validation of the resistance using cell viability tests. b Comparison of total gene expression between duplicate samples from the innate sensitive group (blue) and adaptive resistant (red) group. c Principal component analysis results suggested favorable reproducibility between the sensitive (blue) and resistant (red) groups. d The volcanic map reveals the distributions of downregulated (blue) and upregulated (orange) genes in the resistant group versus the sensitive group. Top 10 of DEGs between two groups were shown. e The heat map indicates the upregulated (red) and downregulated (blue) DEGs between the sensitive (blue) and resistant (red) groups. Each column is a sample and each row is a gene. Top 10 of DEGs between two groups were shown. f GO biological process term enrichment analysis of DEGs (left) and their interrelationships (right). g GO molecular function term enrichment analysis of DEGs (left) and their interrelationships (right). h GO cellular component term enrichment analysis of DEGs (left) and their interrelationships (right). i KEGG pathway enrichment analysis of DEGs (left) and their interrelationships (right). DEGs: Differentially expressed genes. GO: Gene Ontology. KEGG: Kyoto Encyclopedia of Genes and Genomes

Table 1 The expression of the top 20 upregulated and downregulated genes in the cetuximab resistant CACO2
Table 2 GO term and KEGG pathway enrichment analyses of differentially expressed genes

Gene set enrichment

GSEA and GSVA of the DEGs (TablesĀ 3 and 4, respectively) indicated that the downregulated genes were enriched in ā€œHALLMARK_COMPLEMENTā€ (Fig.Ā 3a), ā€œHALLMARK_BILE_ACID_METABOLISMā€ (Fig.Ā 3b), ā€œHALLMARK_IL2_STATS_SIGNALINGā€ (Fig.Ā 3c), ā€œHALLMARK_UV_RESPONSE_DNā€ (Fig.Ā 3d), ā€œHALLMARK_ESTROGEN_RESPONSE_EARLYā€ (Fig.Ā 3e), ā€œHALLMARK_APOPTOSISā€ (Fig.Ā 3f), ā€œHALLMARK_INFLAMMATORY_RESPONSEā€ (Fig.Ā 3g), and ā€œHALLMARK_ESTROGEN_RESPONSE_LATEā€ (Fig.Ā 3h). The upregulated genes were enriched in ā€œHALLMARK_PEROXISOMEā€ (Fig.Ā 3i). The GSVA results suggested that the downregulated genes were enriched in ā€œHALLMARK_MYC_TARGETS_V1 and V2ā€ and upregulated genes were enriched in ā€œHALLMARK_MYOGENESIS,ā€ ā€œHALLMARK_ANGIOGE-NESIS,ā€ ā€œHALLMARK_INTERFERON_ALPHA_RESPONSE,ā€ and ā€œHALLMARK_TGF_ BETA_SIGNALINGā€ (Fig.Ā 3j).

Table 3 Gene set enrichment analysis of upregulated and downregulated gene sets in the cetuximab resistant group
Table 4 Gene set variation analysis of upregulated and downregulated gene sets in the cetuximab resistant group
Fig. 3
figure 3

GSEA and GSVA of DEGs. a-h GSEA of downregulated and (i) upregulated DEGs involved in biological states and processes in the resistant group, compared with those in the sensitive group. j Volcanic diagram of the GSVA results for the DEGs, which were enriched in the resistant group or the sensitive group (highlighted orange points and blue points, respectively). DEGs: Differentially expressed genes. GSEA: Gene Set Enrichment Analysis. GSVA: Gene Set Variation Analysis

Common DEGs and PPI network construction

We identified 708 (Fig.Ā 4a) and 298 (Fig.Ā 4b) DEGs in GSE5851 and GSE59857, respectively. After integrating the DEGs from our own screenings with those from the two GEO datasets, 12 common DEGs were identified (Fig.Ā 4c). These DEGs did not represent any protein-level interactions (Fig.Ā 4d). Then, we combined recognized cetuximab resistance-related genes from the GDSC database (TableĀ 5) with the 12 common DEGs and constructed a PPI network (Fig.Ā 4e). The following six core DEGs were identified: SATB-2, ORP-1, MYB, CDX-2, SLC26A3, and EPHB2. Consequently, we selected the GSE56386 dataset to validate the variations in the expression levels of the 12 DEGs between cetuximab responders and non-responders in clinical terms. SATB-2, MYB, CDX-2, SLC26A3, and FAR2 were downregulated (Fig.Ā 5a, cā€“e, g) and AKR1B10 was upregulated in the non-responder group (Fig.Ā 5l). Although a trend of downregulation was observed in the other six genes (Fig.Ā 5b, f, hā€“k), they exhibited no significant differences between two groups, possibly owing to the small sample size. These 12 genes had an optimal prediction accuracy in GSE56386. The receiver operating characteristic curves for the 12 genes are presented in Additional file 2: Fig. S1.

Fig. 4
figure 4

Identification of common DEGs and construction of PPI network. a, b Identification of DEGs between the resistant (red) and sensitive (blue) groups from the GSE5851 and GSE59857 datasets, illustrated by a heat map. Each column is a sample and each row is a gene. c The venn diagram shows the intersection of DEGs among GSE5851, GSE59857 and CACO2-CR cellular model. d The Network diagram illustrates the interactions of common DEGs. e The PPI network including 12 DEGs and recognized cetuximab resistance-related genes from the Genomics of Drug Sensitivity in Cancer database. DEGs: Differentially expressed genes. PPI: protein-protein interaction network

Table 5 Cetuximab resistance-related genes from the Genomics of Drug Sensitivity in Cancer database
Fig. 5
figure 5

Comparison of the expression levels of core DEGs in GSE56386. aĀ SATB-2, cĀ MYB, dĀ CDX-2, eĀ SLC26A3, and gĀ FAR2 were significantly downregulated and (l)AKR1B10 was apparently upregulated among cetuximab non-responders compared with cetuximab responders. No significant differences were seen in (b)ORP-1, fĀ EPHB2, hĀ VAV3, iĀ ATP10B, jĀ ADH1A, kĀ ADH1C. DEGs: Differentially expressed genes

Core DEGs and tumor immune cell infiltration

We observed significant correlations between the expression levels of the core DEGs SATB-2 (Fig.Ā 6a), ORP-1 (Fig.Ā 6b), MYB (Fig.Ā 6c), and CDX-2 (Fig.Ā 6d) and tumor immune cell infiltration represented by the expression of B cell, CD4+ T cell, CD8+ T cell, and macrophage markers in TCGA. Immunohistochemical staining of a tissue microarray containing 102 CRC clinical samples (Fig.Ā 7a) demonstrated that the expression of these four genes was positively associated with tumor-infiltrating immune cell markers such as CD4, CD8, CD19, and CD68 (Fig.Ā 7b; TableĀ 6).

Fig. 6
figure 6

Association between the expression level of common DEGs and immune cell infiltration. The relationship between expression levels of (a) SATB-2, (b) ORP-1, (c) MYB, and (d) CDX-2 and tumor purity, infiltrating B cells, CD4+ T cells, CD8+ T cells, and macrophages was investigated by the online tool Tumor Immune Estimation Resource Version. DEGs: Differentially expressed genes. TIMER: Tumor Immune Estimation Resource Version

Fig. 7
figure 7

Relationship between core DEGs and tumor-infiltrating immune cell markers in CRC. a The representative immunohistochemical staining images of SATB-2, ORP-1, MYB, CDX-2, CD19, CD4, CD8, and CD68. b Spearman analysis results indicated that core DEGs were positively associated with tumor-infiltrating immune cell markers. DEGs: Differentially expressed genes. CRC: Colorectal cancer

Table 6 Spearman correlation analysis of differentially expressed genes and tumor-infiltrating immune cell markers in clinical samples

Exploring the expression of identified core DEGs in vitro

We conducted preliminary experiments in several CRC cell lines and detected that ORP-1, MYB, and CDX-2 were downregulated in NCIH508wtRAS/RAF (a cell line sensitive to cetuximab), while SATB-2 expression decreased in CACO2 wtRAS/RAF (a cell line partially sensitive to cetuximab) and HT29 wtRAS/mtRAF (a cell line resistant to cetuximab) at the half-maximal inhibition concentration of cetuximab. A reduction in ORP-1 and CDX-2 expression was observed in most of these cell lines regardless of RAS and RAF status (Fig.Ā 8a).

Fig. 8
figure 8

Exploring the expression of core DEGs in vitro and their prognostic roles in CRC. a SATB-2, ORP-1, MYB, and CDX-2 were downregulated after cetuximab treatment. b Patients with high expression levels of core DEGs in primary tumors exhibited a tendency to experience a longer progression-free survival time. **, pā€‰<ā€‰0.01, ***, pā€‰<ā€‰0.001, ns, pā€‰>ā€‰0.05. DEGs: Differentially expressed genes. CRC: Colorectal cancer

The roles of identified core DEGs in clinical samples with anti-EGFR therapy

The clinicopathological characteristics of the patients were recorded and are shown in Additional file 3: Table S2. Median follow-up time was 33.1 months (interquartile range, 17.2 to 52.5 months). Compared with patients with low ORP-1 expression, those with high ORP-1 expression in CRC experienced a significantly longer PFS time (median times were 9.0 months and 11.0 months, respectively, hazard rate [HR]ā€‰=ā€‰1.901, pā€‰=ā€‰0.047) after administration of chemotherapeutics containing cetuximab (Fig.Ā 8b).

Discussion

Alterations in genes that function in the epidermal growth factor (EGF) signaling pathway, such as increases in EGFR copy number, amplification of ERBB family, overexpression of IGF1 or VEGF, or novel mutations such as point mutations in RAS, BRAF, PI3KCA, or MEK in the EGFR extracellular domain or in the downstream pathway, result in EGFR antagonist resistance. Recently, disturbances in the tumor microenvironment caused by EGFR antibody have also been recognized as factors in treatment failure. Garvey et al. reported that cetuximab causes cancer-associated fibroblasts to secrete more EGF and reactivate mitogen-activated protein kinase signaling in para-CRC cells [32]. Critically, it is important to understand the interaction between tumor cells and the tumor microenvironment in anti-EGFR therapy resistance.

We conducted an integrated bioinformatic analysis to identify nine common DEGs between the cetuximab sensitive and resistant groups combining high-throughput data of the cetuximab resistant cell line model and data of cell lines and clinical samples from GEO profiles. The relationship between DEGs identified in our study and tumor immune cell infiltration was evaluated in TCGA and validated in clinical samples from our hospital. We found that four (SATB-2, ORP-1, MYB, and CDX-2) of nine DEGs were associated with infiltrated T cells, B cells, and macrophages in CRC. The decreasing trend of expression levels of ORP-1, MYB and CDX-2 under cetuximab pressure was observed in vitro (Additional file 4: Fig. S2), which was consistent with the expression changes in the established cetuximab resistant cell line (CACO2-CR). Moreover, patients with high expression levels of these genes, especially ORP-1, exhibited prolonged survival receiving anti-EGFR therapy.

SATB-2 encodes a DNA binding protein and mediates transcription regulation as well as chromatin remodeling. SATB-2 attenuates the activity of MEK5/ERK5 and suppresses tumorigenesis and metastasis [33]. The upregulation of SATB-2 via DNA demethylation of the promoter region and H3K4me3 increases TH1-type chemokine expression and immune cell density in CRC [34]. ORP-1, as a member of the oxysterol-binding protein family involve in human innate immune system, binds to phosphatidylinositol 3-phosphate by interacting with RAB7A and stabilizing GTP-RAB7A and regulates the MHC class I -mediated antigen processing and presentation pathway [35, 36]. In addition, ORP-1 suppresses tumorigenesis via metabolism-associated pathway [37]. MYB is a protein-encoding oncogene that functions as a transcriptional activator. Paradoxically, patients with CRC and high MYB expression exhibit low incidence of distant metastases [38] and favorable clinical prognosis [39]. Millen et al. demonstrated that a high level of CD8+ tumor immune infiltrating cells and a clinical history of longer relapse-free survival were related to high expression of MYB. An immunomodulatory effect conferred by MYB was also observed in CD8+ TILs in the murine CRC model [40]. CDX-2 encodes a regulator of intestine-specific genes involved in cell growth and differentiation. CDX-2Ā is downregulated in the invasive part of tumor tissues and is associated with tumor-stroma protein expression as well as inflammatory cytokine release in CRC [41]. Low levels of CDX-2 expression indicate a particularly poor survival prognosis, especially in patients with tumors that have a high stromal content [42].

Zanella et al. reported that EGFR antagonists inhibit colorectal tumor growth and simultaneously protect the tumor from inhibition by transcriptional regulation [43]. Woolston et al. supported this finding and revealed that cetuximab prompts a transformation from a mutational variation to a mesenchymal transition representative of tumor-associated fibroblast enrichment [20]. It seems that cetuximab evokes immune inflammation via antibody-dependent cell-mediated cytotoxicity but paradoxically weakens this effect via a form of IgG1 antibody-mediated immunogenic cell death [21] while upregulating immunosuppressive TGF-Ī² expression in CRC [44]. In fact, in a stage Ib/II trial of combined cetuximab and pembrolizumab treatment for patients with advanced CRC, the combination treatment of cetuximab with pembrolizumab significantly increased the density of CD3+, CD8+, and CTLA-4+ lymphocytes and natural killer cells in tumors. In the peripheral blood, the overall density of CD4+ and CD8+ lymphocytes decrease, especially that of the PD1+ memory T cells [45, 46]. In the process of acquiring secondary resistance, the tumor might experience the transition from the immune-inflamed phenotype to the immune-desert phenotype, which features by key target genes involved in tumor immune cell infiltration [47,48,49]. We believed that the core DEGs identified in this study, SATB-2, ORP-1, MYB, and CDX-2, might play critical roles in the transition.

Nevertheless, there were some limitations lying in this study. As models of anti-EGFR antibody resistance are not easy to established, only one cellular model was established and screened in the study. More cetuximab resistant cellular models with wild-type RAS and RAF genotypes are required to be established for validation. The DEGs identified using in silico methods should be further validated in vitro and in vivo. The results drawn from the retrospective clinical cohort with a small sample size also need further verification in a larger population. The specific immune cell types involved in secondary anti-EGFR antibody resistance should be more clearly elucidated.

Conclusions

In summary, we distinguished cetuximab-induced DEGs associated with variations of tumor immune cell infiltration in CRC by establishing a cetuximab resistant cellular model and integrated bioinformatics analysis. The results suggested that transcriptomic alterations and immune landscape remodeling should receive additional scrutiny during anti-EGFR antibody treatment. Furthermore, immunotherapy could be considered in the early stages of cetuximab treatment rather than after resistance has already occurred.

Availability of data and materials

Public datasets of colon adenocarcinoma and rectum adenocarcinoma analyzed during the current study can be retrieved from The Cancer Genome Atlas (TCGA) at https://portal.gdc.cancer.gov/. GSE5851, GSE56368, and GSE59857 datasets can be retrieved from the Gene Expression Omnibus (GEO) at https://www.ncbi.nlm.nih.gov/geo/. Other data that support the findings of this study are available from the corresponding author on reasonable request.

Abbreviations

CRC:

Colorectal cancer

EGFR:

Epidermal growth factor receptor

GEO:

Gene Expression Omnibus

DEGs:

Differentially expressed genes

PCA:

Principal component analysis

GO:

Gene Ontology

KEGG:

Kyoto Encyclopedia of Genes and Genomes

GSEA:

Gene Sets Enrichment Analysis

GSVA:

Gene Set Variation Analysis

PPI:

Protein-protein interaction

GDSC:

Genomics of Drug Sensitivity in Cancer

TCGA:

The Cancer Genome Atlas

OSBP:

Oxysterol-binding protein

TILs:

Tumor immune infiltrating cells

COAD:

Colon adenocarcinoma

READ:

Rectum adenocarcinoma

References

  1. Siegel RL, Miller KD, Jemal A. Cancer statistics, 2020. CA Cancer J Clin. 2020;70(1):7ā€“30.

    ArticleĀ  PubMedĀ  Google ScholarĀ 

  2. Chen W, Zheng R, Baade PD, Zhang S, Zeng H, Bray F, Jemal A, Yu XQ, He J. Cancer statistics in China, 2015. CA Cancer J Clin. 2016;66(2):115ā€“32.

    ArticleĀ  PubMedĀ  Google ScholarĀ 

  3. Heinemann V, von Weikersthal LF, Decker T, Kiani A, Vehling-Kaiser U, Al-Batran SE, Heintges T, Lerchenmuller C, Kahl C, Seipelt G, et al. FOLFIRI plus cetuximab versus FOLFIRI plus bevacizumab as first-line treatment for patients with metastatic colorectal cancer (FIRE-3): a randomised, open-label, phase 3 trial. Lancet Oncol. 2014;15(10):1065ā€“75.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  4. Avallone A, Budillon A. Impact of subsequent therapies on outcome of the FIRE-3/AIO KRK0306 trial. J Clin Oncol. 2016;34(13):1564.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  5. Primrose J, Falk S, Finch-Jones M, Valle J, Oā€™Reilly D, Siriwardena A, Hornbuckle J, Peterson M, Rees M, Iveson T, et al. Systemic chemotherapy with or without cetuximab in patients with resectable colorectal liver metastasis: the New EPOC randomised controlled trial. Lancet Oncol. 2014;15(6):601ā€“11.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  6. Levi FA, Boige V, Hebbar M, Smith D, Lepere C, Focan C, Karaboue A, Guimbaud R, Carvalho C, Tumolo S, et al. Conversion to resection of liver metastases from colorectal cancer with hepatic artery infusion of combined chemotherapy and systemic cetuximab in multicenter trial OPTILIV. Ann Oncol. 2016;27(2):267ā€“74.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  7. Bridgewater JA, Pugh SA, Maishman T, Eminton Z, Mellor J, Whitehead A, Stanton L, Radford M, Corkhill A, Griffiths GO, et al. Systemic chemotherapy with or without cetuximab in patients with resectable colorectal liver metastasis (New EPOC): long-term results of a multicentre, randomised, controlled, phase 3 trial. Lancet Oncol. 2020;21(3):398ā€“411.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  8. Garattini S, Torri V, Floriani I. Cetuximab for metastatic colorectal cancer. N Engl J Med. 2009;361(1):96. author reply 96ā€“97.

    CASĀ  PubMedĀ  Google ScholarĀ 

  9. Venook AP, Niedzwiecki D, Lenz HJ, Innocenti F, Fruth B, Meyerhardt JA, Schrag D, Greene C, Oā€™Neil BH, Atkins JN, 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(23):2392ā€“401.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  10. Adam R, Aloia T, Levi F, Wicherts DA, de Haas RJ, Paule B, Bralet MP, Bouchahda M, Machover D, Ducreux M, et al. Hepatic resection after rescue cetuximab treatment for colorectal liver metastases previously refractory to conventional systemic therapy. J Clin Oncol. 2007;25(29):4593ā€“602.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  11. Bregni G, Sciallero S, Sobrero A. HER2 amplification and anti-EGFR sensitivity in advanced colorectal cancer. JAMA Oncol. 2019;5(5):605ā€“6.

    ArticleĀ  PubMedĀ  Google ScholarĀ 

  12. Waddell T. Targeting EGFR in colorectal cancer: beyond KRAS exon 2. Lancet Oncol. 2014;15(6):540ā€“1.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  13. Diaz LA Jr, Williams RT, Wu J, Kinde I, Hecht JR, Berlin J, Allen B, Bozic I, Reiter JG, Nowak MA, et al. The molecular evolution of acquired resistance to targeted EGFR blockade in colorectal cancers. Nature. 2012;486(7404):537ā€“40.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  14. Sahin IH, Oktay AA. EGFR signalling and oncogenic pathway signature in colorectal cancer. Gut. 2013;62(11):1670.

    ArticleĀ  PubMedĀ  Google ScholarĀ 

  15. Rachiglio AM, Lambiase M, Fenizia F, Roma C, Cardone C, Iannaccone A, De Luca A, Carotenuto M, Frezzetti D, Martinelli E, et al. Genomic profiling of KRAS/NRAS/BRAF/PIK3CA wild-type metastatic colorectal cancer patients reveals novel mutations in genes potentially associated with resistance to anti-EGFR agents. Cancers (Basel). 2019;11(6):859.

    ArticleĀ  CASĀ  Google ScholarĀ 

  16. Douillard JY, Oliner KS, Siena S, Tabernero J, Burkes R, Barugel M, Humblet Y, Bodoky G, Cunningham D, Jassem J, et al. Panitumumab-FOLFOX4 treatment and RAS mutations in colorectal cancer. N Engl J Med. 2013;369(11):1023ā€“34.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  17. Sorich MJ, Wiese MD, Rowland A, Kichenadasse G, McKinnon RA, Karapetis CS. Extended RAS mutations and anti-EGFR monoclonal antibody survival benefit in metastatic colorectal cancer: a meta-analysis of randomized, controlled trials. Ann Oncol. 2015;26(1):13ā€“21.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  18. Bertotti A, Papp E, Jones S, Adleff V, Anagnostou V, Lupo B, Sausen M, Phallen J, Hruban CA, Tokheim C, et al. The genomic landscape of response to EGFR blockade in colorectal cancer. Nature. 2015;526(7572):263ā€“7.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  19. Misale S, Di Nicolantonio F, Sartore-Bianchi A, Siena S, Bardelli A. Resistance to anti-EGFR therapy in colorectal cancer: from heterogeneity to convergent evolution. Cancer Discov. 2014;4(11):1269ā€“80.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  20. Woolston A, Khan K, Cunningham D, Gerlinger M. Genomic and transcriptomic determinants of therapy resistance and immune landscape evolution during anti-EGFR treatment in colorectal cancer. Cancer Cell. 2019;36(1):35-50 e39.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  21. Pozzi C, Cuomo A, Spadoni I, Magni E, Silvola A, Conte A, Sigismund S, Ravenda PS, Bonaldi T, Zampino MG, et al. The EGFR-specific antibody cetuximab combined with chemotherapy triggers immunogenic cell death. Nat Med. 2016;22(6):624ā€“31.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  22. Giordano GR, Porras A, Pancione AM. Immune resistance and EGFR antagonists in colorectal cancer. Cancers (Basel). 2019;11(8):1089.

    ArticleĀ  CASĀ  Google ScholarĀ 

  23. Ritchie ME, Phipson B, Wu D, Hu Y, Law CW, Shi W, Smyth GK. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 2015;43(7):e47.

    ArticleĀ  PubMedĀ  PubMed CentralĀ  CASĀ  Google ScholarĀ 

  24. Zhou Y, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, Benner C, Chanda SK. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun. 2019;10(1):1523.

    ArticleĀ  PubMedĀ  PubMed CentralĀ  CASĀ  Google ScholarĀ 

  25. Subramanian A, Kuehn H, Gould J, Tamayo P, Mesirov JP. GSEA-P: a desktop application for Gene Set Enrichment Analysis. Bioinformatics. 2007;23(23):3251ā€“3.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  26. Ferreira MR, Santos GA, Biagi CA, Silva Junior WA, Zambuzzi WF. GSVA score reveals molecular signatures from transcriptomes for biomaterials comparison. J Biomed Mater Res A 2020.

  27. Yang W, Soares J, Greninger P, Edelman EJ, Lightfoot H, Forbes S, Bindal N, Beare D, Smith JA, Thompson IR, et al. Genomics of Drug Sensitivity in Cancer (GDSC): a resource for therapeutic biomarker discovery in cancer cells. Nucleic Acids Res. 2013;41:D955-61 Database issue.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  28. Szklarczyk D, Gable AL, Nastou KC, Lyon D, Kirsch R, Pyysalo S, Doncheva NT, Legeay M, Fang T, Bork P, et al. The STRING database in 2021: customizable protein-protein networks, and functional characterization of user-uploaded gene/measurement sets. Nucleic Acids Res. 2021;49(D1):D605ā€“12.

    ArticleĀ  PubMedĀ  Google ScholarĀ 

  29. Otasek D, Morris JH, Boucas J, Pico AR, Demchak B. Cytoscape automation: empowering workflow-based network analysis. Genome Biol. 2019;20(1):185.

    ArticleĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  30. Li T, Fu J, Zeng Z, Cohen D, Li J, Chen Q, Li B, Liu XS. TIMER2.0 for analysis of tumor-infiltrating immune cells. Nucleic Acids Res. 2020;48(W1):W509ā€“14.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  31. Li T, Fan J, Wang B, Traugh N, Chen Q, Liu JS, Li B, Liu XS. TIMER: a web server for comprehensive analysis of tumor-infiltrating immune cells. Cancer Res. 2017;77(21):e108ā€“10.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  32. Garvey CM, Lau R, Sanchez A, Sun RX, Fong EJ, Doche ME, Chen O, Jusuf A, Lenz HJ, Larson B, et al. Anti-EGFR therapy induces EGF secretion by cancer-associated fibroblasts to confer colorectal cancer chemoresistance. Cancers (Basel). 2020;12(6):1393.

    ArticleĀ  CASĀ  Google ScholarĀ 

  33. Mansour MA, Hyodo T, Ito S, Kurita K, Kokuryo T, Uehara K, Nagino M, Takahashi M, Hamaguchi M, Senga T. SATB2 suppresses the progression of colorectal cancer cells via inactivation of MEK5/ERK5 signaling. FEBS J. 2015;282(8):1394ā€“405.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  34. Xu M, Xu X, Pan B, Chen X, Lin K, Zeng K, Liu X, Xu T, Sun L, Qin J, et al. LncRNA SATB2-AS1 inhibits tumor metastasis and affects the tumor immune cell microenvironment in colorectal cancer by regulating SATB2. Mol Cancer. 2019;18(1):135.

    ArticleĀ  PubMedĀ  PubMed CentralĀ  CASĀ  Google ScholarĀ 

  35. Johansson M, Lehto M, Tanhuanpaa K, Cover TL, Olkkonen VM. The oxysterol-binding protein homologue ORP1L interacts with Rab7 and alters functional properties of late endocytic compartments. Mol Biol Cell. 2005;16(12):5480ā€“92.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  36. Suchanek M, Hynynen R, Wohlfahrt G, Lehto M, Johansson M, Saarinen H, Radzikowska A, Thiele C, Olkkonen VM. The mammalian oxysterol-binding protein-related proteins (ORPs) bind 25-hydroxycholesterol in an evolutionarily conserved pocket. Biochem J. 2007;405(3):473ā€“80.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  37. Motazacker MM, Pirhonen J, van Capelleveen JC, Weber-Boyvat M, Kuivenhoven JA, Shah S, Hovingh GK, Metso J, Li S, Ikonen E, et al. A loss-of-function variant in OSBPL1A predisposes to low plasma HDL cholesterol levels and impaired cholesterol efflux capacity. Atherosclerosis. 2016;249:140ā€“7.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  38. Tichy M, Knopfova L, Jarkovsky J, Pekarcikova L, Veverkova L, Vlcek P, Katolicka J, Capov I, Hermanova M, Smarda J, et al. Overexpression of c-Myb is associated with suppression of distant metastases in colorectal carcinoma. Tumour Biol. 2016;37(8):10723ā€“9.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  39. Tichy M, Knopfova L, Jarkovsky J, Vlcek P, Katolicka J, Capov I, Hermanova M, Smarda J, Benes P. High c-Myb expression associates with good prognosis in colorectal carcinoma. J Cancer. 2019;10(6):1393ā€“7.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  40. Millen R, Malaterre J, Cross RS, Carpinteri S, Desai J, Tran B, Darcy P, Gibbs P, Sieber O, Zeps N, et al. Immunomodulation by MYB is associated with tumor relapse in patients with early stage colorectal cancer. Oncoimmunology. 2016;5(7):e1149667.

    ArticleĀ  PubMedĀ  PubMed CentralĀ  CASĀ  Google ScholarĀ 

  41. Olsen J, Espersen ML, Jess P, Kirkeby LT, Troelsen JT. The clinical perspectives of CDX2 expression in colorectal cancer: a qualitative systematic review. Surg Oncol. 2014;23(3):167ā€“76.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  42. Sandberg TP, Sweere I, van Pelt GW, Putter H, Vermeulen L, Kuppen PJ, Tollenaar R, Mesker WE. Prognostic value of low CDX2 expression in colorectal cancers with a high stromal contentā€”a short report. Cell Oncol (Dordr). 2019;42(3):397ā€“403.

    ArticleĀ  CASĀ  Google ScholarĀ 

  43. Zanella ER, Galimi F, Sassi F, Migliardi G, Cottino F, Leto SM, Lupo B, Erriquez J, Isella C, Comoglio PM, et al. IGF2 is an actionable target that identifies a distinct subpopulation of colorectal cancer patients with marginal response to anti-EGFR therapies. Sci Transl Med. 2015;7(272):272ra212.

    ArticleĀ  CASĀ  Google ScholarĀ 

  44. Tauriello DVF, Palomo-Ponce S, Stork D, Berenguer-Llergo A, Badia-Ramentol J, Iglesias M, Sevillano M, Ibiza S, Canellas A, Hernando-Momblona X, et al. TGFbeta drives immune evasion in genetically reconstituted colon cancer metastasis. Nature. 2018;554(7693):538ā€“43.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  45. Boland PM, Ma WW. Immunotherapy for colorectal cancer. Cancers (Basel). 2017;9(5):50.

    ArticleĀ  CASĀ  Google ScholarĀ 

  46. Le DT, Kim TW, Van Cutsem E, Geva R, Jager D, Hara H, Burge M, Oā€™Neil B, Kavan P, Yoshino T, et al. Phase II open-label study of pembrolizumab in treatment-refractory, microsatellite instability-high/mismatch repair-deficient metastatic colorectal cancer: KEYNOTE-164. J Clin Oncol. 2020;38(1):11ā€“9.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  47. Inoue Y, Hazama S, Suzuki N, Tokumitsu Y, Kanekiyo S, Tomochika S, Tsunedomi R, Tokuhisa Y, Iida M, Sakamoto K, et al. Cetuximab strongly enhances immune cell infiltration into liver metastatic sites in colorectal cancer. Cancer Sci. 2017;108(3):455ā€“60.

    ArticleĀ  CASĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

  48. Wang L, Wei Y, Fang W, Lu C, Chen J, Cui G, Diao H. Cetuximab enhanced the cytotoxic activity of immune cells during treatment of colorectal cancer. Cell Physiol Biochem. 2017;44(3):1038ā€“50.

    ArticleĀ  CASĀ  PubMedĀ  Google ScholarĀ 

  49. Sunakawa Y, Yang D, Cao S, Zhang W, Moran M, Astrow SH, Hsiang J, Stephens C, Tsuji A, Takahashi T, et al. Immune-related genes to dominate Neutrophil-lymphocyte Ratio (NLR) associated with survival of cetuximab treatment in metastatic colorectal cancer. Clin Colorectal Cancer. 2018;17(4):e741-9.

    ArticleĀ  PubMedĀ  PubMed CentralĀ  Google ScholarĀ 

Download references

Acknowledgements

Not applicable.

Funding

This study received funding from the Outstanding Youth Foundation of Zhongshan Hospital (No. 2019ZSYQ21), CSCO-Youth Innovation Research Fund (Y-young2019-057), and National Natural Science Foundation of China (81602038, 81772511 and 81900482). The funding sources had no role in the design of the study; collection, analysis, and interpretation of data; or writing the manuscript.

Author information

Authors and Affiliations

Authors

Contributions

TL: conceptualization, methodology, supervision, resources. YC: methodology, software, resources. LL: data curation, formal analysis, writingā€”original draft preparation. ML: investigation, validation, visualization, writingā€” review & editing. XS, YY, KP and KR: data curation, writingā€”review & editing. YY and YC: project administration, writingā€”reviewing and editing. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Yanjie Chen or Tianshu Liu.

Ethics declarations

Ethics approval and consent to participate

All procedures performed in this study involving human participants were in accordance with the ethical standards of the institutional and/or national research committee and with the 1964 Helsinki declaration and its later amendments or comparable ethical standards. The Ethics Committee of Zhongshan Hospital Fudan University approved this study.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisherā€™s note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1: Table S1.

Sequences of primers used for real-time quantitative PCR.

Additional file 2: Figure S1.

ROC curves for the cetuximab sensitivity of the 12 DEGs in GSE56386. DEGs: Differentially expressed genes. ROC: Receiver operating characteristic.

Additional file 3: Table S2.

Clinicopathological data of patients with colorectal cancer in the tissue microarray.

Additional file 4: Figure S2.

SATB-2, ORP-1, MYB, and CDX-2 were downregulated in the resistant CACO2. ##, pā€‰<ā€‰0.01, ###, pā€‰<ā€‰0.001.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Liang, L., Liu, M., Sun, X. et al. Identification of key genes involved in tumor immune cell infiltration and cetuximab resistance in colorectal cancer. Cancer Cell Int 21, 135 (2021). https://doi.org/10.1186/s12935-021-01829-8

Download citation

  • Revised:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12935-021-01829-8

Keywords