An anoikis-related signature predicts prognosis and immunotherapy response in gastrointestinal cancers

RNA-sequencing (RNA-seq) data and complete clinical information for esophageal carcinoma (ESCA), stomach adenocarcinoma (STAD), colon adenocarcinoma (COAD) and rectum adenocarcinoma (READ) samples from a total of 1255 patients were extracted from The Cancer Genome Atlas (TCGA) official website (https://portal.gdc.cancer.gov/↗). The RNA-seq data were downloaded as Fragments Per Kilobase of transcript per Million mapped reads (FPKM) values and log2-transformed for subsequent analysis. Additionally, we obtained the RNA-seq data and clinical data of 483 GC patients in the GSE84437↗ dataset and 585 CRC patients in the GSE39582↗ dataset from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/gds↗). Furthermore, 647 ARGs in total were acquired from GeneCards (https://www.genecards.org/↗) and Harmonizome (https://maayanlab.cloud/Harmonizome/↗).

The differential expression analysis between normal and tumor samples was performed using Wilcoxon rank-sum test. The criteria for identifying differentially expressed ARGs were set as |log2FC| > 1 and false discovery rate (FDR) < 0.05. Gene expression values were normalized prior to analysis. Univariate Cox proportional hazards regression analysis was conducted to identify ARGs significantly associated with overall survival (OS) (P < 0.05). In accordance with the expression profiles of these ARGs, the R package ‘ConsensusClusterPlus’ was utilized to perform clustering and identify anoikis-related patterns. The consensus matrix was used to acquire the optimal number of clusters. Principal component analysis (PCA), t-distributed stochastic neighbor embedding (t-SNE) and uniform manifold approximation and projection (UMAP) were adopted to evaluate the reliability of anoikis-related patterns. Kaplan-Meier (KM) curves were generated to evaluate the prognostic differences between clusters. Single-sample gene set enrichment analysis (ssGSEA) was applied to obtain the immune cell infiltration abundance. Gene set variation analysis (GSVA) and gene set enrichment analysis (GSEA) were used to conduct Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis, characterizing the functions and pathways of genes in different clusters.

All GI patients (n = 2220) from the TCGA and GEO databases were randomly divided into two equal-sized cohorts using the “caret” R package: a training cohort (n = 1110) and a test cohort (n = 1110). The training cohort was utilized to develop the risk score model, named Anoscore, while the test cohort was used for internal validation. To ensure the robustness of our validation process, we maintained similar distributions of key clinical characteristics between the two cohorts through stratified random sampling. The training cohort was utilized to develop Anoscore through a two-step selection process. First, least absolute shrinkage and selection operator (LASSO) regression was applied to reduce dimensionality and select the most significant candidate genes among the ARGs associated with OS. Multivariate Cox regression analysis was then performed on the LASSO-selected genes to further refine the model. The Anoscore was calculated using the following formula: Anoscore = Σ in(Coefi * Xi). X represents the expression level of each gene, and Coefi represents the weighted coefficient of the corresponding gene. The accuracy of the Anoscore in the training cohort was evaluated using KM and receiver operating characteristic (ROC) curves.

Subsequently, the performance of Anoscore was evaluated in the test cohort through multiple approaches for internal validation. KM and ROC curves were plotted to validate the accuracy of the prognostic signature in the test cohort. Moreover, the independent prognostic value of the Anoscore was determined via multivariate Cox regression analysis. Furthermore, the associations among cluster A, cluster B, the Anoscore and patient survival were analyzed and visually represented.

An Ano-nomogram was constructed using the ‘rms’ R package to predict OS by integrating Anoscore with key clinical features, including age, sex, and pathologic N stage, in GI tumor patients. The Cox proportional hazards regression model was used as the basis for the nomogram construction. The performance of the Ano-nomogram was evaluated through various validation methods. Calibration curves were generated to assess the agreement between predicted and observed survival probabilities at 1, 3, and 5 years. Additionally, decision curve analysis (DCA) was performed to evaluate the clinical utility of the Ano-nomogram compared with individual predictors, including Anoscore, age, gender, T stage, and N stage, at different time points. Cumulative hazard curves were also generated by stratifying patients into high- and low-risk groups based on Ano-nomogram scores to visualize the model’s discriminative ability.

The tumor immune single-cell hub 2 (TISCH2) database (http://tisch.comp-genomics.org/↗), which provides detailed cell type annotations at the single-cell level, was used to analyze the expression of S100A11 and TLR3 at the single-cell level. The correlations between the expression of genes and infiltrating immune cells were also assessed, enabling us to explore the TME in detail.

Drug response information and information on drug targeting pathways were obtained from the Genomics of Drug Sensitivity in Cancer (GDSC) database (https://www.cancerrxgene.org/↗). Drug sensitivity was predicted using the oncoPredict algorithm with empirical Bayes batch correction. Correlation analyses were performed between the Anoscore and the predicted drug sensitivity to identify potential therapeutic implications.

The TME of GI cancers was comprehensively characterized. The XCELL, TIMER, QUANTISEQ, MCPOUNTER, EPIC, CIBERSORT-ABS and CIBERSORT algorithms were used to estimate the correlation between immune cell infiltration and the Anoscore. The ESTIMATE algorithm was employed to calculate immune scores, stromal scores and ESTIMATE scores for both low- and high-Anoscore groups. Additionally, the expression of common immune checkpoints (ICPs) was assessed between the two subgroups based on the Anoscore. ICPs with significant differential expression (P < 0.05) were visualized using boxplots.

The human EC cell line KYSE-150, the GC cell line HGC-27 and the CRC cell line Caco-2 were acquired from the China Center for Type Culture Collection (CCTCC; Wuhan, China). KYSE-150 cells and HGC-27 cells were cultured in RPMI 1640 media (Gibco, New York, USA) supplemented with 10% fetal bovine serum (FBS, Gibco, New York, USA) and 1% antibiotics (Zqxzbio, Shanghai, China), while Caco-2 cells were cultured in DMEM (Gibco, New York, USA) supplemented with 15% FBS and 1% antibiotics. All the cells were incubated at 37°C with 5% CO₂.

Since S100A11 and TLR3 were identified as the most significant contributors to the Anoscore, we investigated their biological functions in detail. S100A11 was knocked down with a small interfering RNA (siRNA) in KYSE-150, HGC-27 and Caco-2 cells (the sequences were listed in,). The siRNAs and the negative control were transfected into cells using Lipofectamine 3000 (Invitrogen, Carlsbad, CA, USA). The pcDNA3.1-control and pcDNA3.1-TLR3 plasmids were designed and synthesized (Boshang Biotechnology, Jinan, China). The plasmids were transfected into the aforementioned cells using Lipofectamine 3000 (Invitrogen, Carlsbad, CA, USA). Quantitative real-time polymerase chain reaction (qRT-PCR) and western blotting (WB) were used to confirm the efficiency of knockdown and overexpression. 1 1

Cell proliferation was detected using the Cell Counting Kit-8 (CCK-8) (Bioss, Beijing, China). Transfected cells (3×10³ cells/well) and the corresponding controls were seeded into 96-well plates (Corning Incorporated, Corning, NY, USA). After culture for 24 h, 48 h, 72 h or 96 h, 10 microliters (μl) of the CCK-8 solution was added to each well, which contained 90 μl of serum-free medium. The cells were incubated for 2 h at 37 °C. Finally, we measured the absorbance with a spectrophotometer (Tecan, Männedorf, Switzerland) at 450 nm and generated growth curves.

Transfected cells (KYSE-150, HGC-27 and Caco-2 cells) (1×10⁵ cells/well) with 200 μl of serum-free RPMI 1640 media or DMEM were seeded into the upper chamber of a 24-well transwell system (Corning Incorporated, Corning, NY, USA) in triplicate, while 800 μl of RPMI 1640 media or DMEM containing 20% FBS was added to the lower chamber. After a certain number of hours, the migrated cells were fixed with methanol, and then a crystal violet solution (0.1%) was used to stain them. The migrated cells were observed and counted under a microscope (100×) (Olympus, Tokyo, Japan).

S100A11-knockdown and TLR3-overexpressing cells were plated into 6-well plates at a density of 1000 cells/well. After 10 days of incubation, the cells were fixed with 4% paraformaldehyde and stained with 0.5% crystal violet. Visible colonies (more than 50 cells) were counted.

Total RNA was extracted using the RNA Fast 2000 (Fastagen, Shanghai, China). The RNA concentration and purity were measured, with the quality criterion of A260/A280 ratio between 1.8 and 2.0. For cDNA synthesis, 2.0 micrograms (μg) of total RNA was reverse transcribed in a 20 μl reaction using Evo M-MLV RT Premix for qPCR (Accurate Biology, Hunan, China). The reverse transcription reaction was carried out at 37 °C for 15 minutes, followed by 85 °C for 5 seconds to inactivate the enzyme. qPCR was performed using the BlazeTaq SYBR Green qPCR mix 2.0 Kit (Accurate Biology, Hunan, China). The cycling conditions were as follows: initial denaturation at 95°C for 30 seconds, with 40 cycles consisting of denaturation at 95°C for 5 seconds and annealing and extension at 60°C for 30 seconds. All RNA extraction and qRT-PCR experiments were conducted according to the manufacturers’ instructions. The relative expression levels were calculated using the relative 2^-ΔΔCT method with β-actin as the internal control. The sequences of primers used were listed in Additional File 1, Supplementary Table S2.

The cells were lysed on ice with radioimmunoprecipitation assay (RIPA) lysis buffer (Solarbio, Beijing, China) containing phenylmethylsulfonyl fluoride (1 mM) and then homogenized via ultrasonication. The protein concentration was determined with Omni-Easy™ Instant BCA Protein Assay Kit (EpiZyme, Shanghai, China). The cell lysate was mixed with loading buffer (Beyotime, Shanghai, China) and denatured at 95°C for 10 minutes. The protein samples were electrophoresed on sodium dodecylsulfate−polyacrylamide gels and transferred to polyvinylidene fluoride membranes (Millipore, MA, USA). The membranes were blocked with 5% milk and then incubated with primary antibodies against S100A11 (1:1000, Proteintech, 10237-1-AP, China), TLR3 (1:1000, ABclonal, A11778, China), β-actin (1:2000, HUABIO, EM21002, China), BAX (1:20000, HUABIO, ET1702-53, China), Bcl-2 (1:2000, HUABIO, ET1702-53, China), Caspase-3 (1:2000, HUABIO, ET1602-39, China), active Caspase-3 (1:1000, HUABIO, ET1602-47, China), Caspase-7 (1:1000, HUABIO, JE59-36, China), active Caspase-7 (1:1000, HUABIO, ER60002, China), PARP1 (1:10000, Proteintech, 66520-1-Ig, China), and cleaved PARP1 (1:10000, Proteintech, 60555-1-Ig, China) overnight at 4°C. After an incubation with a secondary antibody (1:5000) for 1 hour at room temperature, exposure and visualization were conducted.

Statistical analysis was conducted using R (version 4.2.2), GraphPad Prism (version 8.0) and SPSS (version 23.0) software. The Wilcoxon test was used for the comparison of two independent samples. T-tests and one-way ANOVA were used to compare parametric data. The Kruskal-Wallis test was used to compare nonparametric data across multiple samples. A p-value < 0.05 was considered statistically significant.

We obtained RNA transcriptome data and corresponding clinical data from ESCA (N = 185), STAD (N = 443), COAD (N = 459) and READ (N = 168) patients in TCGA database. Among the 647 ARGs extracted from the GeneCards and Harmonizome databases, 142 differentially expressed ARGs were recognized between tumor and normal tissues, including 109 upregulated and 33 downregulated genes in tumors (Figures 1A, B). After combining the expression data and clinical information of TCGA, GSE39582↗ and GSE84437↗ cohorts, we performed a univariate Cox analysis to identify differentially expressed ARGs associated with the prognosis. The results suggested that 87 ARGs were associated with the survival status and survival events of patients. Among them, the expression of 49 genes, such as CHEK2, BID, CDC25C and PTRH2, was negatively correlated, and the expression of 38 genes, such as MAPK10, NTRK3, PDGFB and BDNF, was positively correlated with survival (Figure 1C). A network plot was constructed to visualize the intricate relationship between ARGs and their prognostic value for GI tumors (Figure 1D). Copy number variations (CNVs) in ARGs were visualized on chromosomes (Figure 1E). Our investigation of 87 ARGs revealed that CNV-related mutations were prevalent, with 49 genes exhibiting widespread CNV amplification and 38 genes exhibiting CNV deletions (Figure 1F).

We performed an unsupervised clustering analysis to examine whether anoikis-related patterns could be used to classify patients with GI cancers based on the expression of the 87 prognostic ARGs. The results showed that k = 2 was the optimal parameter, which meant that patients were classified into two anoikis-related patterns, termed ARG clusters A and B (Figure 2A). ARG clusters A and B were distinguished dramatically via PCA (Figure 2B), t-SNE (Figure 2C), and UMAP (Figure 2D) analyses, which validated our anoikis-related patterns. A boxplot was generated to visualize the differentially expressed ARGs in the two clusters, in which 37 ARGs, such as CHEK2, CDC25C and TRAF2, were upregulated in ARG cluster A, and 48 ARGs, such as PDK4, CXCL12, MAPK10 and NOTCH3, were upregulated in ARG cluster B (Figure 2E). A heatmap was generated to visualize the detailed expression of ARGs and the clinical characteristics of clusters A and B (Figure 2F). The KM analysis revealed a difference in patient survival between ARG clusters A and B. The prognosis of patients in ARG cluster B was dramatically worse than that of patients in ARG cluster A (p < 0.001) (Figure 2G). By analyzing immune cell infiltration via ssGSEA, we observed that immune cells, such as CD8-positive T lymphocytes (CD8⁺ T cells), activated natural killer (NK) cells, dendritic cells (DCs), macrophages, and activated B cells, were more abundant in ARG cluster B (Figure 2H). These results indicated that anoikis-related patterns were correlated with the immune landscape of GI cancers. ARG cluster A tended to represent cold tumors, whereas ARG cluster B tended to represent hot tumors. Thus, patients in ARG cluster B might be more sensitive to immunotherapy. GSVA revealed that multiple cancer-related pathways, such as ECM-receptor interaction, focal adhesion, regulation of the actin cytoskeleton, the mitogen-activated protein kinase (MAPK) signaling pathway, gap junctions and cell adhesion molecules, were significantly activated in ARG cluster B (Figure 2I). GSEA revealed that DNA replication and peroxisomes were more strongly activated in ARG cluster A (Figure 2J). In addition, ECM-receptor interaction, focal adhesion and cell adhesion molecules were more enriched in ARG cluster B, which was consistent with the results of GSVA (Figure 2J).

We systematically divided patients into a training cohort (n =1110) and a test cohort (n = 1110). Then we used the LASSO algorithm (Figure 3A) and multivariate Cox analysis (Figure 3B) to further identify the significant prognostic ARGs for patients in the training cohort. Finally, we obtained 12 ARGs with p < 0.05 that were utilized to construct the prognostic signature. The risk score of each patient was calculated as follows: Anoscore = 0.179 × ExpPDK4 + 0.078 × ExpCDKN2A − 0.180 × ExpMYC + 0.130 × ExpFGF2 – 0.291 × ExpTLR3 + 0.124 × ExpIL1RAP + 0.120 × ExpCD36 − 0.222 × ExpTRAF2 + 0.145 × ExpTNFRSF12A + 0.134 × ExpNOTCH3 + 0.260 × ExpS100A11 − 0.109 × ExpWNT2 (Additional File 1, Supplementary Table S3). The heatmap depicted the expression levels of the ARGs that comprised the Anoscore in different risk groups (Figure 3C). FGF2, PDK4, CD36, CDKN2A, S100A11, IL1RAP and NOTCH3 were significantly upregulated in patients with high Anoscores (Figure 3C).

The KM analysis revealed a significant relationship between the Anoscore and the prognosis of patients with GI cancers. ROC curves were constructed to evaluate the sensitivity and specificity of the risk model. The outcomes were assessed based on the area under the ROC curve (AUC). The model demonstrated robust performance in both cohorts. In the training cohort, significant prognostic stratification was observed, with the high-Anoscore group showing markedly shorter survival times compared to the low-Anoscore group (p < 0.001, Figure 3D). Furthermore, the model exhibited strong discriminative power, achieving AUCs of 0.689, 0.707, and 0.695 for 1-, 3-, and 5-year survival predictions, respectively (Figure 3E). Validation in the test cohort confirmed the model’s reliability, with significant prognostic stratification maintained (p < 0.001) and the low-Anoscore group consistently demonstrating better clinical outcomes (Figure 3F). The discriminative performance of Anoscore remained stable, with AUCs of 0.628, 0.660, and 0.637 for 1-, 3-, and 5-year survival predictions (Figure 3G). The similar performance metrics between the training and test cohorts indicated good model generalizability.

The KM analysis further revealed that all patients with high Anoscores had generally worse prognoses (p < 0.001) (Figure 3H), with 1-, 3-, and 5-year AUCs of 0.656, 0.683, and 0.666, respectively (Figure 3I). The robustness of Anoscore was further confirmed through cancer-type specific analyses, showing consistent performance. We evaluated the relationships between OS and the Anoscore in patients with ESCA, STAD, COAD and READ using KM analyses. Patients with high Anoscores had shorter OS times, consistent with the analysis of all patients with GI cancers (Additional File 2, Supplementary Figure S1). Similarly, the Anoscore proved highly reliable in patients with STAD, COAD and READ (Additional File 2, Supplementary Figure S1). However, perhaps due to the small sample size in the ESCA cohort, the survival analysis was not statistically significant. Moreover, the multivariate Cox analysis confirmed the Anoscore as an independent prognostic factor (hazard ratio [HR] 1.45, 95% confidence interval [CI] 1.35–1.56) (Figure 3J). Considering these results, a high Anoscore can be recognized as an effective indicator of unfavorable OS among patients with GI cancers.

Consistent with our previous findings, we observed the biological consistency between Anoscore and the ARG clusters (Figure 3K). We also illustrated the interrelationship among the two ARG clusters, Anoscore typing, and the survival status of patients in a Sankey diagram. Patients in ARG cluster A predominantly corresponded to the low-Anoscore group, of which surviving patients constituted the majority (Figure 3L).

We constructed a novel nomogram capable of predicting survival probabilities at 1, 3, and 5 years that comprised the clinical factors and the Anoscore (Figure 3M). The hybrid nomogram was stable and accurate, potentially offering valuable clinical utility in patient management. Strong agreement between the actual survival probabilities and the Ano-nomogram-predicted survival probabilities was exhibited in the calibration curve (Figure 3N). The cumulative hazards of different groups were calculated by analyzing the Nelson-Aalen cumulative risk curve. The cumulative risk of the low Ano-nomogram group increased at a slower rate over time, indicating that patients with low Ano-nomogram scores might have better survival outcomes, as shown in Figure 3O. The 1-year DCA curve revealed that most of the Ano-monogram curve was in the area above the other lines, which demonstrated that using the Ano-nomogram to reach a decision on the prognosis was more reliable and could result in a greater net benefit (Figure 3P). The 3-year DCA curve (Figure 3Q) and 5-year DCA curve (Figure 3R) also revealed that the Ano-nomogram had outstanding performance.

S100A11 is a calcium-binding protein that regulates cellular processes including cell cycle and differentiation (8), while TLR3 is a pattern recognition receptor crucial for immune response and cell death signaling (9). According to our previous LASSO algorithm and multivariate Cox analyses, S100A11 had the largest positive coefficient, and TLR3 had the largest negative coefficient. Therefore, we investigated the expression of S100A11 and TLR3 in single-cell sequencing data from TISCH2. In EC patients (GSE160269↗), the expression of these two genes could be clustered into 31 groups (Figure 4A), and the clustering scatter plot showed the clustering results (Figure 4B). The pie chart showed the number of different cell types, and the stacked bar chart showed the proportions of different cells in each patient (Figure 4C). S100A11 was highly distributed and expressed in malignant cells, while TLR3 expression was lower (Figures 4D, E). The distribution and expression of S100A11 and TLR3 in different cell types of the dataset were clearly revealed in the violin plot (Figure 4F). The results revealed that S100A11 was expressed mainly in tumor cells, monocytes/macrophages and fibroblasts in GSE160269↗ (Figure 4F). In contrast, TLR3 was expressed at very low levels in all cell types, especially in tumor cells (Figure 4F). In GSE134520↗ (GC), the clustering analysis revealed 22 clusters (Figures 4G, H), of which pit mucous cells were the majority (Figure 4I). The distribution and expression of S100A11 and TLR3 expression in GSE134520↗ were shown in the map (Figures 4J, K). S100A11 was enriched in mast cells, DCs, malignant cells and myofibroblasts (Figure 4L). Meanwhile, TLR3 was enriched in DCs and myofibroblasts (Figure 4L). With respect to CRC patients, the clustering analysis clustered the GSE166555↗ dataset into 32 groups (Figure 4M). The distribution of each cell type was shown in Figure 4N. The proportion of each cell type was shown in Figure 4O. S100A11 was expressed at high levels in monocytes/macrocytes, DCs and malignant cells, whereas TLR3 was expressed at slightly higher levels in fibroblasts and endothelial cells and at low levels in other cells (Figures 4P–R).

Immunotherapy plays an important role in tumor treatment. Considering the relatively low response rate to immunotherapy and the lack of promising predictors, we conducted an immune-related analysis to evaluate whether the Anoscore was associated with the immune status and could predict the immunotherapy response in GI tumors. We used seven different software packages to estimate the correlation between the Anoscore and tumor-infiltrating immune cells. We found that the high Anoscore was significantly related to increased immune cell infiltration, as shown in the immune cell bubble chart (Figure 5A). In general, B cells were positively correlated with the Anoscore across the above cancer types according to CIBERSORT, CIBERSORT-ABS, EPIC, MCPCOUNTER, QUANTISEQ, and XCELL (Figure 5B). DCs were positively correlated with the Anoscore in CIBERSORT, CIBERSORT-ABS, MCPCOUNTER, QUANTISEQ, TIMER and XCELL (Figure 5C). Macrophages and neutrophils were positively correlated with the Anoscore in 6 and 4 datasets, respectively (CIBERSORT-ABS, EPIC, MCPCOUNTER, QUANTISEQ, TIMER, and XCELL for macrophages and CIBERSORT, CIBERSORT-ABS, MCPCOUNTER, and TIMER for neutrophils; Additional File 2, Supplementary Figures S2A, B). In addition, with CIBERSORT-ABS, EPIC and XCELL, the number of CD4-positive T lymphocytes (CD4⁺ T cells) was positively correlated with the Anoscore (Additional File 2, Supplementary Figure S2C). Moreover, CD8⁺ T cells were positively correlated with the Anoscore according to the QUANTISEQ, TIMER and XCELL tools (Additional File 2, Supplementary Figure S2D). NK cells (Additional File 2, Supplementary Figure S2E) and monocytes (Additional File 2, Supplementary Figure S2F) were also positively correlated with the Anoscore. Overall, tumors in the high-Anoscore group exhibited much greater immune cell infiltration, indicating that these tumors were hot and could respond to immunotherapy more sensitively. According to TME scores, the infiltration of stromal cells and immune cells was greater in the high-Anoscore group than in the low-Anoscore group, indicating that the TME differed between the low- and high-Anoscore groups (Figure 5D). We also analyzed the relationships between the Anoscore and ICPs. Cooccurrences of ligands and receptor expression were identified. Nine pairs, including PDCD1LG2/CD274-PDCD1, CD86-CD28/CTLA4, TNFSF18-TNFRSF18, ICOSLG-ICOS, TNFRSF14-CD160/TNFSF14/BTLA, TNFSF9-TNFRSF9, CD40-CD40LG, LGALS9-HAVCR2 and CD200-CD200R1, were highly expressed in the high-Anoscore group (Figure 5E), suggesting that appropriate immune checkpoint inhibitors could be effectively applied, based on the Anoscore for GI cancers. Next, we evaluated the associations between the Anoscore and the infiltration of immune cells in ESCA (Figure 5F), STAD (Figure 5G), COAD (Figure 5H) and READ (Figure 5I). Among these four tumors, the high Anoscore was associated with increased immune cell infiltration, which is consistent with the previous overall data.

Drug treatment, including chemotherapy and targeted drugs, plays a crucial role in cancer management. However, the therapeutic outcomes vary significantly due to differences in drug sensitivity. We evaluated the value of the Anoscore in predicting drug sensitivity in GI cancers. The low- and high-Anoscore groups were differentially sensitive to 165 types of drugs in the GDSC2 database. Among them, the high-Anoscore group was more sensitive to 141 types of drugs, whereas the low-Anoscore group was more sensitive to 24 types of drugs (Additional File 1, Supplementary Table S4). The high-Anoscore group was more sensitive to multiple chemotherapy drugs, including cisplatin, paclitaxel, oxaliplatin, irinotecan, gemcitabine, cyclophosphamide, epirubicin and 5-fluorouracil (Figures 6A–H, Additional File 2, Supplementary Figure S3). Similar results were obtained with targeted medicine. The patients in the high-Anoscore group exhibited increased sensitivity to various targeted medicines, including dabrafenib, erlotinib, foretinib, gefitinib, nilotinib, sorafenib, osimertinib and palbociclib (Figures 6I–P, Additional File 2, Supplementary Figure S3). The low-Anoscore group was more sensitive to drugs that included dasatinib, doramapimod, IGF1R_3801, AMG-319, AZD1332, AZD2014 and AZD5582 (Additional File 2, Supplementary Figure S3, Additional File 1, Supplementary Table S4). Our analysis suggested that the two groups had different sensitivities to medicine. Overall, patients with high Anoscores were more sensitive to drug treatment, especially chemotherapy and targeted therapy.

We conducted experiments to explore the biological functions of S100A11 and TLR3 in vitro. KYSE-150, HGC-27 and Caco-2 cells with S100A11 knockdown and TLR3 overexpression were constructed (Figures 7A, B). Afterward, we evaluated the expression of apoptosis-related proteins via WB analysis. Ultralow-attachment 6-well plates were used to suspend EC, GC and CRC cells, analogous to the conditions of anoikis. After cells were cultured in suspension for 24 h, S100A11 knockdown and TLR3 overexpression increased the expression of active caspase 3, active caspase 7, cleaved PARP and Bax in KYSE-150, HGC-27 and Caco-2 cells but reduced the expression of Bcl-2 (Figures 7C, D, Additional File 2, Supplementary Figure S4). CCK-8 assays revealed that S100A11 knockdown and TLR3 overexpression inhibited cell proliferation (Figures 7E, F). Furthermore, we conducted colony formation assays, which revealed fewer colonies in the S100A11-knockdown and TLR3-overexpressing groups than in the negative control groups (Figures 7G, H). The results of transwell migration experiments revealed that S100A11 knockdown and TLR3 overexpression decelerated the migration of cells (Figures 7I, J). These results indicated that S100A11 might promote the proliferation and migration of tumor cells, whereas TLR3 might have the opposite effects.

In summary, S100A11 knockdown and TLR3 overexpression promoted anoikis and diminished the malignancy of tumor cells, indicating that S100A11 expression might be positively correlated with anoikis resistance and a poor prognosis, whereas TLR3 expression was negatively correlated. The in vitro validation of S100A11 and TLR3 functions was consistent with the results of our bioinformatics analysis described above, enhancing the reliability of our research.

The original contributions presented in the study are included in the article/. Further inquiries can be directed to the corresponding authors. 1