Identification of core SASP-related genes TGFBI and ANXA6 for diagnosi

Introduction

Ulcerative colitis (UC) is a chronic, recurrent inflammatory bowel disease with an unclear etiology, that may be closely related to multiple factors such as genetic, environmental, microbial and immune factors.1–4 UC has become a global disease, and in the early 21st century, the prevalence of UC is increasing in some newly industrialized countries, which leads to a large increase in the burden of disease.5,6 UC patients have a high risk of developing ulcerative colitis-related colorectal cancer (UCRCC), and a meta-analysis showed that the cumulative risk of UCRCC in UC patients with 30 years of onset was as high as 13.9%.7,8 Although recent research reports indicate that the incidence of colorectal cancer (CRC) has been declining for a long time, the risk for patients with chronic ulcerative colitis remains higher compared to the general population.9 Furthermore, there are still no effective indicators for diagnosing UCRCC at present.

Cellular senescence is a stable terminal state of cells characterized by the loss of proliferative and differentiation capabilities.10 Cells can transition to a senescent state due to various stressors, such as oncogene-induced senescence, therapy-induced senescence, mitochondrial dysfunction-induced senescence, and immune-induced senescence.11–15 The senescence-associated secretory phenotype (SASP) is the primary mediator secreted by senescent cells in the tissue microenvironment, including chemokines, growth factors, cytokines, and stromal metalloproteinases.16 Increasing evidence suggests that senescent cells mainly promote disease progression by secreting SASP.16–18 In UC and CRC, various methods can induce cellular senescence, including therapy-induced senescence, oncogene-induced senescence, and immune-induced senescence. Among these, the senescence induced by changes in the immune microenvironment mediated by chronic inflammation has received considerable attention from researchers in recent years.19–22 It is reported that the SASP-related genes (SASPRGs) can promote the malignant progression of CRC.23,24 However, due to the lack of research on the SASPRGs in UC and CRC, it is currently unclear which SASPRGs contribute to their progression, so identifying core SASPRGs may provide a theoretical basis for the diagnosis, prevention, and formulation of personalized treatment strategies for these diseases.

In this study, we evaluated SASP-related differentially expressed genes in UC and CRC patients using bulk RNA transcriptomics. We determined the core SASPRGs using 113 combinations of 10 machine learning algorithms and identified their diagnostic potential and clinical value. Furthermore, we also screened potential therapeutic drugs for patients with high expression of TGFBI and ANXA6, and constructed a molecular regulatory network, providing new insights for more effective diagnostic and therapeutic strategies for UCRCC.

Materials and Methods

Datasets and Processing

Gene expression matrices and related clinical information about UC and CRC patients were obtained from the TCGA (https://portal.gdc.cancer.gov/) and GEO (https://www.ncbi.nlm.nih.gov/geo/) databases. Supplementary Table 1 presents the full information about the dataset used in this study.

Differential Expression and Functional Enrichment Analysis

The R package “limma” was used to identify differential expression genes (DEGs) in the training datasets. Any gene with adjusted values of P < 0.05 and |log2 (FC)| of > 0.8 was regarded as DGEs. Gene Ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) and Gene Set Enrichment Analysis (GSEA) analyses were performed using the R package “clusterProfiler”. The visualization of this result was plotted by R package “ggplot2”.

Weight Gene Co-Expression Network Analysis (WGCNA)

The weight gene co-expression network analysis (WGCNA) was performed using the R package “WGCNA”. Firstly, the hierarchical clustering method is used to filter the dataset, and the “pickSoftThreshould” function is used to determine the appropriate soft threshold β. Finally, the corresponding modules are generated. Among these modules, the one most relevant to the clinical phenotype was selected.

Identification of Core SASPRGs Associated with UC and CRC

Human putative SASP factor genes were extracted from the Gene Ontology Consortium and QuickGO database.25 DEGs and key module genes identified by WGCNA were intersected with SASP to obtain SASPRGs associated with UC and CRC. To further identify core SASPRGs, we calculated the AUC values of 113 combinations of 10 machine learning algorithms based on these results, respectively. Ultimately, the genes included in the combination with the highest AUC value were considered core SASPRGs. The median was used to categorize the core immune-related genes into low and high expression groups, and then Kaplan–Meier analysis and Cox regression analysis were performed with the TCGA-COAD dataset downloaded from TCGA databases.

Immune Cell Infiltrates Analysis

We evaluated the immune cell infiltration in tissue from UC and CRC patients by performing ssGSEA, CIBERSORT and MCPcounter analyses using the R package “GSVA”, “CIBERSORT”, and “MCPcounter”. Finally, we explored the relationship between the 28 immune cells and core SASPRGs using Spearman correlation analysis.

Single-Cell Analysis

The GSE162335 and GSE231559 datasets were used for single-cell analysis. The R package “Seurat” was used for subsequent data processing. TSNE analysis was used to identify different clusters. Using the “BlueprintEncodeData” dataset from the R package “celldex” as a reference, each cluster was initially annotated with “SingleR”. The R package “UCell” is used to evaluate the expression of core SASPRGs in different cell types, the UCell score of each cell was mapped to the TSNE embedding with “ggplot2” for visualization. Bubble plots were used to represent the expression levels of core SASPRGs in different cell types. In addition, expression patterns of core SASPRGs in different patient tissues were validated at the single-cell level.

Prediction of the Immunotherapy Response and Drug Sensitivity Analysis

CRC patients were divided into two groups according to the average expression level of SASPRGs, respectively. The immunotherapy response scores of different groups were calculated by TIDE (http://tide.dfci.harvard.edu/login/), and the differences in TIDE scores between different groups were compared by Wilcoxon test. A higher TIDE score indicates a greater likelihood of tumor immune escape.

Drug sensitivity data were obtained from GDSC (https://www.cancerrxgene.org/) and CTRP (https://portals.broadinstitute.org/ctrp/) databases, and then the R package “oncoPredict” was used to calculate the IC50 value of drugs in different subgroups of CRC patients. The Wilcoxon test was used to compare the difference of IC50 values between different groups.

To further analyze the correlation between the identified drugs and SASP-related proteins, we first obtained the 3D structure of the drugs from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/) and the 3D structure of SASP-related proteins from the PDB database (https://www.ebi.ac.uk/pdbe/) and then used the Autodock tool to predict the binding ability of the drugs to the proteins.

Construction of a Molecular Regulator Network

To predict the regulatory mechanisms of ANXA6 and TGFBI, we constructed a molecular regulatory network using the NetworkAnalyst (https://www.networkanalyst.ca/NetworkAnalyst/) platform. We firstly obtained gene-miRNA interaction information from miRTarBase and TF-gene interaction information from JASPAR. Finally, Cytoscape was used to visualize the molecular regulatory network.

Construction of the Animal Model

The UC animal model was established according to the previous experimental methods. First, eight-week-old C57BL/6 male mice were acclimatized in a specific pathogen-free animal facility for 7 days. Normal drinking water was subsequently replaced with a 3% dextran sulfate sodium (DSS, meilunbio, MB5535) solution for 7 days to induce UC. The mice in the control group were provided with equal volumes of water throughout the course of the experiment. Mice were euthanized on day 8, and colon tissues were collected for further studies. Body weight, fecal occult blood, and fecal viscosity were recorded daily to assess disease activity index (DAI). The CRC animal model previously constructed by our team was used in this study. The specific methods were as follows: Step 1: Eight-week-old C57BL/6 male mice were placed in a specific pathogen-free animal facility. The mice in the experimental group were intraperitoneally injected with azomethane (AOM) at a dose of 10mg/kg and fed with normal drinking water for 1 week. Step 2: Normal drinking water feeding was replaced with 2.5% DSS for one week. Step 3: Feeding with normal drinking water for 2 weeks. Step 4: Repeat steps 2 and 3 three times. All experiments followed the guidelines approved by the Ethical Review Committee for Animal Experiments of Nanjing Medical University.

RNA Isolation and qRT-PCR

Total RNA was extracted from colon tissues by Trizol reagent (Invitrogen) and used to synthesize cDNA with the PrimeScriptTM RT reagent Kit (TaKaRa, Dalian, China) according to the manufacturer’s instructions. qRT-PCR was performed on an ABI 7900HT PCR sequencer (Applied Biosystems, Massachusetts, USA) using TB GreenR Premix Ex Taq™ II (TaKaRa). PCR sequencer (Applied Biosystems, Massachusetts, USA) using TB GreenR Premix Ex Taq™ II (TaKaRa). Fold changes in mRNA expression were calculated using 2−ΔΔCt and standardized based on GAPDH. Primer sequences used for qRT-PCR were as follows: ANXA6, forward: 5′-AGAGCTACAAGTCCCTCTACG-3′, reverse: 5′-CCCACAATCAACCGTTCAAAC-3′. TGFBI, forward: 5′-CACTCTCAAACCTTTACGAGACC-3′, reverse: 5′-CGTTGCTAGGGGCGAAGATG-3′. GAPDH, forward: 5′-TGCACCACCAACTGCTTAGC-3′, reverse: 5′-GGCATGGACTGTGGTCATGAG-3′.

Histopathological Staining

Mouse colon tissues were fixed overnight in universal tissue fixative solution (Servicebio, G1101). The fixed tissues were embedded in paraffin and sliced in a gradient concentration of alcohol after dehydration. After staining with hematoxylin and eosin, the slices were dehydrated and sealed; observations were made using an upright light microscope (NIKON ECLIPSE E100, Nikon, Japan) and imaging system (NIKON DS-U3, Nikon, Japan).

Immunostaining

Antigens of tissue slides were repaired using EDTA antigen repair solution (Servicebio, G1203). After completion of repair, the tissues were blocked with 10% donkey serum for 30 min. After the addition of the primary antibody, the slices were incubated in a wet box overnight at 4° C. The next day, tissue slices were washed three times using 1*PBS (PH7.4), then secondary antibodies were added and incubated for 50 min at room temperature in the dark. Finally, tissue slices were sealed using an anti-fluorescence quench sealant (Servicebio, G1221). Images were acquired by an Upright fluorescence microscope (NIKON ECLIPSE C1, Nikon, Japan) and scanner (Pannoramic MIDI, 3DHISTECH).

Data Statistics

All statistical analyses were performed with R version 4.2 software and its resource packages. The differences between the two groups were compared using t-test (normally distributed data) or Wilcoxon test (non-normally distributed data). p < 0.05 was considered statistically significant.

Results

Identification of Shared DEGs in UC and CRC

Firstly, we identified 105 shared genes in UC and CRC through DEGs and WGCNA analyses. Next, we identified 35 SASPRGs from the intersection of 105 genes with human putative SASP factor genes. Using machine learning algorithms and K-M survival analyses, we finally determined 2 core SASPRGs related to poor prognosis. The flowchart of identifying hub genes is shown in Figure 1. In GSE87466 dataset of UC patients, a total of 1642 DEGs were identified, including 1061 up-regulated genes and 581 down-regulated genes (Figure 2A). In GSE44076 dataset of CRC patients, 2950 DEGs were identified, including 1546 up-regulated genes and 1404 down-regulated genes (Figure 2B). In order to identify shared DEGs, we intersected 1642 DEGs in GSE87211 dataset with 2950 DEGs in GSE44076 dataset. Finally, 632 genes were identified.

Figure 1 Flowchart of the research.

Abbreviations: UC, ulcerative colitis; CRC, colorectal cancer; DEGs, differentially expressed genes; WGCNA, weighted gene co-expression network analysis; SASP, senescence-associated secretory phenotype; SASPRGs, SASP-related genes.

Figure 2 Identification of shared DEGs genes in UC and CRC. (A and B) Volcano plots of all DEGs in GSE87466 and GSE44076, green indicates downregulated DEGs, and red indicates upregulated DEGs. (C and D) The selection of a soft threshold in UC (GSE87466) and CRC (GSE44076). (E and F) Gene cluster dendrogram for UC and CRC by dynamic tree cut algorithm. (G and H) Heatmap of the association between modules and clinical traits in UC and CRC. (I) Venn diagram demonstrates the intersection of genes obtained by WGCNA and DEGs. (J) The enrichment analysis results of the GO and KEGG pathways.

Abbreviations: BP, Biological Process; CC, Cell Component; MF, Molecular Function.

To explore the biological functions and pathways of 632 shared DEGs. GO and KEGG analyses were performed based on these genes. As for GO, these genes were significantly enriched in the processes of leukocyte migration, cell chemotaxis, extracellular matrix organization, leukocyte chemotaxis, and humoral immune response. The KEGG pathway enrichment analysis indicated that these genes were significantly enriched in pathways, such as viral protein interaction with cytokine and cytokine receptor, protein digestion and absorption, cytokine-cytokine receptor interaction, and IL-17 signaling pathway (Figure 2J).

Identification and Validation of Core SASPRGs Associated with UC and CRC

WGCNA was conducted to identify biologically significant gene modules, enhancing our understanding of genes associated with the pathogenesis of UC and CRC in training datasets. We set the optimal soft threshold to 2 in GSE87466 dataset and 9 in GSE44076 dataset (Figure 2C and D). Subsequently, 23 modules in the GSE87466 dataset and 13 modules in the GSE44076 dataset were identified (Figure 2E and F). We selected the modules most relevant to clinical features as key modules. The result showed that the turquoise module was significantly associated with UC (correlation coefficient=0.72, p=3e-18) and the blue module was significantly associated with CRC (correlation coefficient=−0.93, p=1e-107) (Figure 2G and H). We selected the genes based on MM>0.6 and GS>0.4 as the candidate genes. Using Venn diagrams, we identified 105 co-pathogenic genes for UC and CRC in the training datasets (Figure 2I).

Next, we intersected 105 genes with human putative SASP factor genes, resulting in 35 SASPRGs (Figure 3A). To further identify core SASPRGs, we used 113 combinations of 10 machine learning algorithms and calculated their AUC values to evaluate the diagnostic efficacy of each combination (Figure 3B and C). We found that the Enet[alpha=0.2] and Lasso+GBM models had the highest average AUC value across three UC datasets (0.991) and three CRC datasets (0.995). Therefore, we selected the genes identified by the Enet[alpha=0.2] and Lasso+GBM algorithms as candidate SASPRGs. We obtained 12 core SASPRGs through a Venn diagram (MMP7, S100A11, TGFBI, EGFL6, CXCL3, NPY, TNFSF11, CTSH, ANXA6, STC2, LRRN2, SRPX2) (Figure 3D).

Figure 3 Identification of core SASPRGs through multiple machine learning algorithms. (A) Venn diagram demonstrates the intersection of genes obtained by shared genes and SASP. (B and C) A total of 113 combinations of machine learning algorithms for the core SASPRGs. The AUC of each model was calculated across 3 datasets. (D) Venn diagram demonstrates the intersection of core SASPRGs obtained by Enet[alpha=0.2] and Lasso+GBM algorithms.

To further evaluate the predictive power of core SASPRGs for CRC patient prognosis, Cox regression and K-M survival analyses were performed using TCGA-COAD data, which revealed that TGFBI, ANXA6, and NPY among the 12 core SASPRGs were potentially associated with poorer prognosis (All p values were less than 0.05) (Figure 4A and D). Interesting, the results of the Kaplan–Meier Plotter platform (http://kmplot.com/analysis/) showed that there was no significant difference in recurrence-free survival probability of CRC patients in the high NPY group compared with the low NPY group (p=0.26) (Figure 4E). Consequently, we excluded NPY from further validation. Subsequently, we constructed a clinical characteristic map for CRC patients. Patients with varying expression levels of ANXA6 and TGFBI exhibited distinct patterns of clinical characteristics (Supplementary Figure 1). Furthermore, the GSE92415 dataset was used to detect the diagnostic efficacy of TGFBI and ANXA6 in assessing disease risk for UC patients. The results showed that both TGFBI and ANXA6 were significantly higher in the high Mayo score group (>5 score) compared to the low Mayo score group (<=5 score), with AUC values of 0.772 (TGFBI) and 0.773 (ANXA6), respectively (Figure 4B and C). These results indicate that TGFBI and ANXA6 have clinical significance in stratifying disease risk for UC and CRC patients.

Figure 4 Clinical value of ANXA6 and TGFBI in UC and CRC patients. (A) Univariate Cox regression analysis for disease-free survival in TCGA-COAD cohorts. (B and C) Clinical value of ANXA6 and TGFBI in UC patients. (D) K-M survival analysis of high-expression versus low-expression groups of TGFBI, ANXA6, and NPY in TCGA-COAD dataset. (E) K-M survival analysis of high-expression versus low-expression groups of TGFBI, ANXA6, and NPY in Kaplan–Meier Plotter platform.

Furthermore, we evaluated expression level and diagnostic efficacy of ANXA6 and TGFBI in validation datasets of UC and CRC. The results showed that the expression levels of TGFBI and ANXA6 were significantly higher in the disease group compared to the normal group. Additionally, TGFBI and ANXA6 effectively distinguished between UC patients, CRC patients, and healthy controls (AUC>0.75) (Figure 5A and B). These findings indicate that TGFBI and ANXA6 have excellent diagnostic performance for UC and CRC.

Figure 5 Validation of the core SASPRGs in UC and CRC datasets. (A) Expression levels and diagnostic potential of ANXA6 and TGFBI in the validation datasets for UC. (B) Expression levels and diagnostic potential of ANXA6 and TGFBI in the validation datasets for CRC. *p < 0.05; ** p < 0.01; *** p < 0.001. **** p < 0.0001.

Biological Functions Associated with ANXA6 and TGFBI in UC and CRC

To explore the biological functions related to ANXA6 and TGFBI, we conducted a GSEA analysis using the genes that are most closely associated with them (p<0.05, |r|>0.5). The results show that ANXA6 and TGFBI are associated with primary immunodeficiency, NF-kappa B signaling pathway, Th1 and Th2 cell differentiation, and B cell receptor signaling pathway in UC (Figure 6A). In CRC, ANXA6 is associated with focal adhesions, protein digestion and absorption, ECM receptor interaction, cell adhesion molecules, and NF-kappa B signaling pathway, while TGFBI is related to the hedgehog signaling pathway, protein digestion and absorption, ECM receptor interaction, Wnt signaling pathway, and carbohydrate digestion and absorption (Figure 6B).

Figure 6 Biological functions associated with ANXA6 and TGFBI in UC and CRC. (A) GSEA analysis shown the biological processes related to ANXA6 and TGFBI in UC. (B) GSEA analysis shown the biological processes related to ANXA6 and TGFBI in CRC.

Expression Patterns of ANXA6 and TGFBI in scRNA-Seq

We analyzed single-cell transcriptome data from the colon tissues of UC and CRC patients to further validate the expression patterns of core SASPRGs. We identified 35241 and 45270 cells in UC and CRC colonic tissues after data quality control, respectively (Figure 7A and B). Based on the expression of marker genes, we annotated 6 cell types for UC and 10 cell types for CRC patients (Figure 7C and D). Subsequently, we evaluated the expression levels of core SASPRGs in different cell types, and the results showed that ANXA6 and TGFBI were mainly expressed in monocytes and fibroblasts (Figure 7E–H). Furthermore, we further validated the expression levels of ANXA6 and TGFBI in both the disease group and the control group at the single-cell level. As shown in Figure 7I–J, expression levels of ANXA6 and TGFBI were significantly higher in the colon tissues of UC and CRC patients (All p values were less than 0.05).

Figure 7 Single-cell transcriptomic characterization of colon tissues from UC and CRC patients. (A and B) The tSNE visualization of different cell groups in UC (GSE162335) and CRC (GSE231559) patients. (C and D) The tSNE visualization after annotation of cell groups in UC and CRC patients. (E and F) The tSNE visualization of expression levels of ANXA6 and TGFBI in UC and CRC patients. (G and H) The dot plot illustrates the expression levels of ANXA6 and TGFBI across different cell populations in UC and CRC patients. (I and J) Violin plot showing expression levels of ANXA6 and TGFBI in Normal vs UC and Normal vs CRC.

Immune Cell Infiltration and Prediction of the Immunotherapy Response

Previous studies have shown that the occurrence and development of UC and CRC are closely related to immune response pathways. Therefore, we performed immune cell infiltration analysis to understand the immune infiltration status of UC and CRC. We found a significant difference in the level of immune cell infiltration between the disease group and the control group (Figure 8A and B and Supplementary Figure 2AD). Compared with the control group, the UC group showed higher levels of immune cells, including Activated B cell, Natural killer cell, Macrophage, Monocyte, T cell CD4 memory activated, Fibroblasts, etc. In CRC group, Activated CD8 T cell, Gamma delta T cell, Activated B cell, Immature B cell, CD56dim natural killer cell, Plasma cells, Macrophages M2, and T cells regulatory showed higher levels, while Activated CD4 T cell, Central memory CD4 T cell, Immature dendritic cell, Macrophages M1, and Fibroblasts showed lower levels. Figure 8C and E show the correlation heatmaps among different immune cells in UC and CRC, further indicating the immune characteristics exhibited by the disease group, as well as highlighting the interactions among various types of immune cells. In addition, we investigated the relationship between core SASPRGs and immune cell components in UC and CRC. The results showed that TGFBI and ANXA6 were significantly positively correlated with Activated CD4 T cell, Central memory CD4 T cell, Type 2 T helper cell, Regulatory T cell, Natural killer cell, Natural killer T cell, and Immature dendritic cell, but negatively correlated with Activated CD8 T cell, Activated B cell, and CD56dim natural killer cell (Figure 8D and F). In order to further explore the relationship between ANXA6 and TGFBI and immune cells. We divided them into high-expression and low-expression groups based on the median. The results showed that compared with the low-expression group, the infiltration level of CD56 bright natural killer cells in the high-expression groups was significantly increased in UC (Supplementary Figure 2E and F), while in CRC, the infiltration levels of T follicular helper cells and Memory B cells were significantly increased (Supplementary Figure 2G and H). This suggests that TGFBI and ANXA6 may play a key role in the pathogenesis of UC and CRC by interfering with immune cell infiltration.

Figure 8 Analysis of immune cell infiltration and prediction of the immunotherapy response. (A and B) The ssGSEA analysis shown enrichment score of immune cells between the UC and CRC disease groups and the control group. (C and E) Heatmap of correlations between different immune cells in colon tissues of patients with UC and CRC. (D and F) Correlation analysis between infiltrating immune cells and core SASPRGs in colon tissues of patients with UC and CRC. (G and I) Expression of the immune checkpoint genes PDCD1, CD274, and CTLA4 in different SASPRGs subgroups. (H and J) TIDE and T-cell exclusion scores in different SASPRGs subgroups. *p < 0.05; ** p < 0.01; *** p < 0.001. **** p < 0.0001.

To evaluate the sensitivity of the core SASPRGs to predict immunotherapy response in CRC patients, we first evaluated the expression levels of PDCD1, CD274, and CTLA4 in different SASPRGs subgroups (Figure 8G and I). High expression of immune checkpoints often suggests that these patients may be suitable for immunotherapy. We then used TIDE to predict immunotherapy efficacy in different SASPRGs subgroups. As expected, the high TGFBI group had lower expression levels of immune checkpoint and higher TIDE scores and T-cell exclusion scores, suggesting that patients did not benefit from immunotherapy (All p values were less than 0.05). Notably, the high ANXA6 group had higher expression levels of immune checkpoint, while the TIDE scores and T-cell exclusion scores were also high (Figure 8H and J).

Screening of Candidate Drugs for High ANXA6 and TGFBI Patients

To screen candidate drugs for CRC patients with high expression of ANXA6 and TGFBI, we evaluated drug sensitivity in different subgroups of patients using drug response data provided by GDSC and CTRP databases. Subsequently, we selected the top three drugs or compounds with the highest sensitivity as candidates (Figure 9A and B). To further analyze the interactions between candidate drugs and SASP-associated proteins, we performed molecular docking. The results demonstrated that ANXA6 (binding energy: Staurosporine: −7.11 kcal/mol, AT7867: −6.31 kcal/mol, BRD_K44224150: −5.93 kcal/mol) and TGFBI (binding energy: SCH772984_1564: −7.47 kcal/mol, VX-11e_2096: −6.21 kcal/mol, ERK_6604_1714: −6.47 kcal/mol) exhibited a certain degree of binding affinity for the three drugs (Figure 9C and D). These candidate drugs have the potential to become therapeutic drugs for CRC.

Figure 9 Drug sensitivity analysis and Construction of a molecular regulator network. (A and B) Outcome of drug sensitivity analysis of three potential therapeutic drugs in different SASPRGs subgroups. (C and D) 3D docking model of three candidate drugs and SASP-associated proteins. (E) Molecular regulatory network based on ANXA6 and TGFBI. Green hexagons represent miRNAs, and blue arrows signify TFs. *p < 0.05; ** p < 0.01; *** p < 0.001. **** p < 0.0001.

Construction of a Molecular Regulator Network

A thorough understanding of the regulatory mechanisms of core SASPRGs is conducive to identifying key biomarkers and enhancing the understanding of diseases. Therefore, we constructed a molecular regulatory network based on ANXA6 and TGFBI to identify their transcriptional and post-transcriptional regulatory molecules. Ultimately, we identified 25 miRNAs targeting TGFBI, 23 miRNAs targeting ANXA6, 13 transcription factors regulating TGFBI, and 9 transcription factors regulating ANXA6. Notably, four transcription factors, USF2, E2F1, TFAP2A, and FOXC1, jointly regulate both ANXA6 and TGFBI (Figure 9E).

Validation of Core SASPRGs Associated with UC and CRC in vivo

To confirm the previous findings, UC model was constructed using DSS. CRC model previously established by our team was used in this study (Figure 10A and B). Colon tissues of mice were collected for histological staining. Compared with the control group, the UC group had disordered colonic mucosa structure, inflammatory cell infiltration, and reduced goblet cells, while the CRC group had disordered colonic gland arrangement and increased cell atypia (Figure 10C). Further evaluation showed that the colon length of the UC group mice was significantly shortened, and the intestinal mucosa showed erosion, bleeding and edema (Figure 10D). Compared with the control group, the body weight of UC group mice decreased significantly with the increase of treatment time, and the DAI score continued to increase (p<0.05) (Figure 10E and F). In addition, we also observed that UC group mice gradually developed symptoms such as depression, hair loss, and diarrhea during the treatment period. At the same time, we determined the mRNA expression levels of ANXA6 and TGFBI in the colon tissues of each group of mice by qRT-PCR. The results showed that the mRNA levels of ANXA6 and TGFBI were significantly increased in the UC and CRC groups compared with the control group (All p values were less than 0.05) (Figure 10H). Immunofluorescence staining further revealed the expression of ANXA6 and TGFBI in the colon tissues of each group of mice (Figure 10G and I).

Figure 10 Expression of core SASPRGs in mouse models of UC and CRC. (A and B) Flow chart of UC mice and CRC mice model establishment. (C) HE staining of colonic tissues from normal, UC, and CRC mice. (D) Colon length of normal and UC mice. (E and F) Body weight and DAI score of normal and UC mice. (G and I) Immunofluorescence staining and relative fluorescence intensity quantification of colonic tissues from normal, UC, and CRC mice. (H) The mRNA expression of ANXA6 and TGFBI in colon tissue detected by qRT-PCR. *p < 0.05; ** p < 0.01; *** p < 0.001. **** p < 0.0001.

Abbreviations: DSS, dextran sulfate sodium; AOM, azomethane; DAI, disease activity index.

Discussion

UC is a common chronic inflammatory bowel disease.26 There is evidence that active UC is prone to progressing to UCRCC, and the risk of UCRCC in UC patients is 2.4 times higher than that in the general population.27 SASP, which has been widely studied in the fields of inflammatory diseases and cancer in recent years, plays a key role in the development of a variety of diseases.28–31 According to the report, GDF15 promotes the malignant progression of CRC by activating the MAPK and PI3K signaling pathways as an essential SASP factor.23 Furthermore, research indicates that SPP1 may be associated with the inhibitory immune microenvironment mediated by senescent cells in CRC.24 However, there is still a lack of SASP-based studies on the molecular mechanisms of UC and CRC. Therefore, exploring the pathogenic genes based on SASP of UC and CRC patients not only enhances our understanding of the pathogenesis of these diseases, but also helps clinicians to make better decisions and improve patient prognosis.

In this study, we identified 35 overlapping SASPRGs in UC and CRC patients by WGCNA as well as differential expression analysis. To comprehensively evaluate the processing performance of different machine learning algorithms and address the specific characteristics of our datasets, we strategically selected 10 representative algorithms. These diverse sets of algorithms allow for a robust comparison across different models, including Lasso, support vector machine (SVM), elastic net (Enet), Ridge, generalized linear model (GLM), random forest (RF), generalized boosted regression modeling (GBM), linear discriminant analysis (LDA), extreme gradient boosting (XGBoost), Naïve Bayes. Subsequently, we identified 12 core SASPRGs using 113 combinations of 10 machine learning algorithms, among which ANXA6 and TGFBI expression were significantly increased in disease groups in multiple datasets and were significantly associated with poor prognosis in CRC patients. Similar results were obtained in mouse models of UC and CRC. Therefore, these two genes were selected for further analysis. By analyzing the scRNA-seq data, ANXA6 and TGFBI were found to be mainly expressed in fibroblasts and monocytes. In addition, we also found that CRC patients with high ANXA6 and TGFBI expression had higher TIDE scores and T cell exclusion scores, indicating that these patients seem to have a poor response to immunotherapy. Further investigation showed that candidate drugs from both pharmacogenomic databases exhibited better therapeutic sensitivity in the high ANXA6 group (Staurosporine, AT7867, BRD_K44224150) and the high TGFBI group (SCH772984_1564, VX-11e_2096, ERK_6604_1714).

Annexin-a6 (ANXA6), a highly conserved Ca2+ -dependent membrane binding protein, is the largest of all annexin families.32 It is located in the cytoplasm and attached to the phospholipid membrane of cells, and is involved in a variety of cellular functions, mainly vesicle trafficking and membrane repair.33 Studies have shown that ANXA6 is differentially expressed in gastric cancer, colorectal cancer, melanoma, cervical cancer, and breast cancer, and plays an important role in the formation, development, and drug resistance of tumors.34–37 ANXA6 shows dual effects in different tumors. Generally, ANXA6 is considered as a tumor suppressor gene mainly because ANXA6 negatively regulates the phosphorylation of EGFR and the downstream MAPK and PI3K-AKT signaling pathways.38–40 Some studies have found that ANXA6 is detected to be elevated in the feces of patients with colorectal cancer, suggesting that ANXA6 is a biomarker for early detection of colorectal cancer.37 This is similar to our conclusion. In addition, ANXA6 can promotes the secretion of pro-inflammatory cytokines that mediate the malignant progression of triple-negative breast cancer.41 Increasing evidence suggests that ANXA6 mediates tumor drug resistance and malignant progression by participating in cholesterol metabolism in LEs, autophagy regulation, interacting with TGFR, and immune-inflammatory responses.36,42–45

TGFBI (transforming growth factor beta induced) is a secreted protein stimulated by TGF-β (transforming growth factor β).46 As a component of extracellular matrix, TGFBI is involved in mediating tumor cell adhesion and motility, metabolic reprogramming, angiogenesis, and immune microenvironment homeostasis.47–49 In recent years, many studies have found that TGFBI may promote tumor cell proliferation by enhancing intercellular signaling, activating the PPARγ or PI3K/AKT/mTOR signaling pathway, and is significantly associated with poor prognosis of patients, including cholangiocarcinoma, urothelial carcinoma, Wilms tumor, and glioma.50–53 A study on metastatic colon cancer found that TGFBI promoted angiogenesis by activating spliceosome and lysosome-related pathways.54 In addition, there was a significant correlation between the expression level of TGFBI and the degree of immune cell infiltration.55 TGFBI secreted by immune cells promotes tumor progression by creating an inhibitory immune microenvironment.56 Interestingly, high TGFBI expression was significantly associated with natalizumab resistance in lung cancer patients.57

Currently, the most commonly used chemotherapeutic agents for CRC are fluorouracil, oxaliplatin, and irinotecan.58 Our results showed that both ANXA6 and TGFBI subgroups showed differences in drug sensitivity to fluorouracil, indicating heterogeneity among the subgroups. Staurosporine, AT7867 and BRD_K44224150 showed good treatment sensitivity in the high ANXA6 group. Staurosporine is an ATP-competitive, nonselective protein kinase inhibitor. Staurosporine has been shown to inhibit the progression of gastric cancer by inhibiting CTNNB1 transcription and Wnt/β-catenin signaling activation.59 AT7867 is a potent AKT inhibitor with demonstrated ability to inhibit tumor cell proliferation and induce tumor cell apoptosis, especially in colorectal cancer.60 Interestingly, all three therapeutically sensitive drugs screened in the high TGFBI group were the ERK inhibitors. Studies have shown that ERK signaling pathway can promote tumor cell migration and invasion by regulating the expression of TGFBI.61,62 These results provide a solid rationale for TGFBI as a novel biomarker for CRC.

In this study, we utilized multiple UC and CRC datasets from the GEO and the TCGA database for a comprehensive and integrated analysis of SASPRGs associated with UC and CRC. However, it is important to acknowledge that this study contains certain limitations. Firstly, our research indicates that ANXA6 and TGFBI, as the functional executors of senescent cells, have significant clinical value in patients with UC and CRC. However, further in vitro and in vivo experiments are needed for verification. Second, the limited sample size may affect the generalization ability of the conclusions. Thirdly, the effects of the selected drugs have not been verified in the CRC model.

Conclusion

In conclusion, NXA6 and TGFBI, as core SASPRGs associated with UC and CRC, showed good diagnostic efficacy. Notably, ANXA6 and TGFBI appear to be effective therapeutic targets in CRC patients. In the future, multi-center and large-sample data sets are needed to verify our conclusions, and subsequent in vivo and in vitro experiments are needed to further explore the possible molecular mechanisms.

Data Sharing Statement

The datasets supporting the conclusions of this article are available in the TCGA (https://portal.gdc.cancer.gov/) and GEO (https://www.ncbi.nlm.nih.gov/geo/) databases. [https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE87466/GSE107499/GSE59071/GSE92415/GSE44076/GSE17536/GSE87211]. The analysis methods and used packages are illustrated in the “Methods” section. All other R code and analyses are available from the author upon request.

Ethics Approval

In accordance with national legislative guidelines on human data, our research does not require approval, is derived from Article 32, items (1) and (2) of the “Measures for the Ethical Review of Life Science and Medical Research Involving Human Subjects” issued by China on February 18, 2023. The legislation in these two items is that research using legally obtained public data, data generated through observation without interfering with public behavior, or research using anonymized information data can be exempted from ethical approval. This study complied with the Declaration of Helsinki and was approved by the Ethical Review Committee for Animal Experiments of the Nanjing Medical University (2011031).

Author Contributions

All authors made a significant contribution to the work reported, whether that is in the conception, study design, execution, acquisition of data, analysis and interpretation, or in all these areas; took part in drafting, revising or critically reviewing the article; gave final approval of the version to be published; have agreed on the journal to which the article has been submitted; and agree to be accountable for all aspects of the work.

Funding

This study was supported by Scientific Research Fund of Henan Provincial Science and Technology Department (No. 182102310302).

Disclosure

The authors report no conflicts of interest in this work.

References

1. Delcò F, Sonnenberg A. Exposure to risk factors for ulcerative colitis occurs during an early period of life. Am J Gastroenterol. 1999;94(3):679–684. doi:10.1111/j.1572-0241.1999.00936.x

2. Du L, Ha C. Epidemiology and pathogenesis of ulcerative colitis. Gastroenterol Clin North Am. 2020;49(4):643–654. doi:10.1016/j.gtc.2020.07.005

3. Keshteli AH, Madsen KL, Dieleman LA. Diet in the pathogenesis and management of ulcerative colitis; a review of randomized controlled dietary interventions. Nutrients. 2019;11(7):1498. doi:10.3390/nu11071498

4. Ordás I, Eckmann L, Talamini M, et al. Ulcerative colitis. Lancet. 2012;380(9853):1606–1619. doi:10.1016/S0140-6736(12)60150-0

5. Kaplan GG, Ng SC. Globalisation of inflammatory bowel disease: perspectives from the evolution of inflammatory bowel disease in the UK and China. Lancet Gastroenterol Hepatol. 2016;1(4):307–316.

6. Ng SC, Shi HY, Hamidi N, et al. Worldwide incidence and prevalence of inflammatory bowel disease in the 21st century: a systematic review of population-based studies. Lancet. 2017;390(10114):2769–2778. doi:10.1016/S0140-6736(17)32448-0

7. Bopanna S, Ananthakrishnan AN, Kedia S, et al. Risk of colorectal cancer in Asian patients with ulcerative colitis: a systematic review and meta-analysis. Lancet Gastroenterol Hepatol. 2017;2(4):269–276. doi:10.1016/S2468-1253(17)30004-3

8. Jess T, Simonsen J, Jørgensen KT, et al. Decreasing risk of colorectal cancer in patients with inflammatory bowel disease over 30 years. Gastroenterology. 2012;143(2):375–381.e371. quiz e313-374. doi:10.1053/j.gastro.2012.04.016

9. Shah SC, Itzkowitz SH. Colorectal cancer in inflammatory bowel disease: mechanisms and management. Gastroenterology. 2022;162(3):715–730.e713. doi:10.1053/j.gastro.2021.10.035

10. Gorgoulis V, Adams PD, Alimonti A, et al. Cellular senescence: defining a path forward. Cell. 2019;179(4):813–827.

11. Serrano M, Lin AW, McCurrach ME, et al. Oncogenic ras provokes premature cell senescence associated with accumulation of p53 and p16INK4a. Cell. 1997;88(5):593–602.

12. Demaria M, O’Leary MN, Chang J, et al. Cellular senescence promotes adverse effects of chemotherapy and cancer relapse. Cancer Discovery. 2017;7(2):165–176. doi:10.1158/2159-8290.CD-16-0241

13. Schmitt CA, Wang B, Demaria M. Senescence and cancer – role and therapeutic opportunities. Nat Rev Clin Oncol. 2022;19(10):619–636.

14. Wiley CD, Velarde MC, Lecot P, et al. Mitochondrial dysfunction induces senescence with a distinct secretory phenotype. Cell Metab. 2016;23(2):303–314.

15. Faust HJ, Zhang H, Han J, et al. IL-17 and immunologically induced senescence regulate response to injury in osteoarthritis. J Clin Invest. 2020;130(10):5493–5507.

16. Wang B, Han J, Elisseeff JH, et al. The senescence-associated secretory phenotype and its physiological and pathological implications. Nat Rev Mol Cell Biol. 2024;25(12):958–978.

17. Rao SG, Jackson JG. SASP: tumor suppressor or promoter? Yes! Trends Cancer. 2016;2(11):676–687.

18. Dong Z, Luo Y, Yuan Z, et al. Cellular senescence and SASP in tumor progression and therapeutic opportunities. Mol Cancer. 2024;23(1):181. doi:10.1186/s12943-024-02096-7

19. Zhang F, Guo J, Yu S, et al. Cellular senescence and metabolic reprogramming: unraveling the intricate crosstalk in the immunosuppressive tumor microenvironment. Cancer Communications. 2024;44(9):929–966. doi:10.1002/cac2.12591

20. Zingoni A, Antonangeli F, Sozzani S, et al. The senescence journey in cancer immunoediting. Mol Cancer. 2024;23(1):68. doi:10.1186/s12943-024-01973-5

21. Hsu SK, Li CY, Lin IL, et al. Inflammation-related pyroptosis, a novel programmed cell death pathway, and its crosstalk with immune therapy in cancer treatment. Theranostics. 2021;11(18):8813–8835. doi:10.7150/thno.62521

22. Chibaya L, Snyder J, Ruscetti M. Senescence and the tumor-immune landscape: implications for cancer immunotherapy. Semi Cancer Biol. 2022;86(Pt 3):827–845. doi:10.1016/j.semcancer.2022.02.005

23. Guo Y, Ayers JL, Carter KT, et al. Senescence-associated tissue microenvironment promotes colon cancer formation through the secretory factor GDF15. Aging Cell. 2019;18(6):e13013. doi:10.1111/acel.13013

24. Yu S, Chen M, Xu L, et al. A senescence-based prognostic gene signature for colorectal cancer and identification of the role of SPP1-positive macrophages in tumor senescence. Front Immunol. 2023;14:1175490. doi:10.3389/fimmu.2023.1175490

25. Sturmlechner I, Zhang C, Sine CC, et al. p21 produces a bioactive secretome that places stressed cells under immunosurveillance. Science. 2021;374(6567):eabb3420. doi:10.1126/science.abb3420

26. Ungaro R, Mehandru S, Allen PB, et al. Ulcerative colitis. Lancet. 2017;389(10080):1756–1770. doi:10.1016/S0140-6736(16)32126-2

27. Jess T, Rungoe C, Peyrin-Biroulet L. Risk of colorectal cancer in patients with ulcerative colitis: a meta-analysis of population-based cohort studies. Clin Gastroenterol Hepatol. 2012;10(6):639–645. doi:10.1016/j.cgh.2012.01.010

28. Li X, Li C, Zhang W, et al. Inflammation and aging: signaling pathways and intervention therapies. Signal Transduction Targeted Ther. 2023;8(1):239. doi:10.1038/s41392-023-01502-8

29. Yasuda T, Koiwa M, Yonemura A, et al. Inflammation-driven senescence-associated secretory phenotype in cancer-associated fibroblasts enhances peritoneal dissemination. Cell Rep. 2021;34(8):108779. doi:10.1016/j.celrep.2021.108779

30. Faget DV, Ren Q, Stewart SA. Unmasking senescence: context-dependent effects of SASP in cancer. Nat Rev Cancer. 2019;19(8):439–453. doi:10.1038/s41568-019-0156-2

31. Wang B, Kohli J, Demaria M. Senescent cells in cancer therapy: friends or foes? Trends Cancer. 2020;6(10):838–857. doi:10.1016/j.trecan.2020.05.004

32. Moss SE, Morgan RO. The annexins. Genome Biol. 2004;5(4):219. doi:10.1186/gb-2004-5-4-219

33. Boye TL, Jeppesen JC, Maeda K, et al. Annexins induce curvature on free-edge membranes displaying distinct morphologies. Sci Rep. 2018;8(1):10309. doi:10.1038/s41598-018-28481-z

34. Mussunoor S, Murray GI. The role of annexins in tumour development and progression. J Pathol. 2008;216(2):131–140. doi:10.1002/path.2400

35. Grewal T, Wason SJ, Enrich C, et al. Annexins – insights from knockout mice. Biol Chem. 2016;397(10):1031–1053. doi:10.1515/hsz-2016-0168

36. Qi H, Liu S, Guo C, et al. Role of annexin A6 in cancer. Oncol Lett. 2015;10(4):1947–1952. doi:10.3892/ol.2015.3498

37. Komor MA, Bosch LJ, Coupé VM, et al. Proteins in stool as biomarkers for non-invasive detection of colorectal adenomas with high risk of progression. J Pathol. 2020;250(3):288–298. doi:10.1002/path.5369

38. Koese M, Rentero C, Kota BP, et al. Annexin A6 is a scaffold for PKCα to promote EGFR inactivation. Oncogene. 2013;32(23):2858–2872. doi:10.1038/onc.2012.303

39. de Muga S V, Timpson P, Cubells L, et al. Annexin A6 inhibits Ras signalling in breast cancer cells. Oncogene. 2009;28(3):363–377. doi:10.1038/onc.2008.386

40. Grewal T, Koese M, Rentero C, et al. Annexin A6-regulator of the EGFR/Ras signalling pathway and cholesterol homeostasis. Int J Biochem Cell Biol. 2010;42(5):580–584. doi:10.1016/j.biocel.2009.12.020

41. Sakwe NI, Vuong NB, Black PJ, et al. Annexin A6 modulates the secretion of pro-inflammatory cytokines and exosomes via interaction with SNAP23 in triple negative breast cancer cells. bioRxiv: Preprint Serv Biol. 2024. doi:10.1101/2024.10.22.619710

42. Korolkova OY, Widatalla SE, Williams SD, et al. Diverse roles of annexin A6 in triple-negative breast cancer diagnosis, prognosis and EGFR-targeted therapies. Cells. 2020;9(8):1855. doi:10.3390/cells9081855

43. Grewal T, Rentero C, Enrich C, et al. Annexin animal models-from fundamental principles to translational research. Int J Mol Sci. 2021;22(7). doi:10.3390/ijms22073439

44. García-Melero A, Reverter M, Hoque M, et al. Annexin A6 and late endosomal cholesterol modulate integrin recycling and cell migration. J Biol Chem. 2016;291(3):1320–1335. doi:10.1074/jbc.M115.683557

45. Cao J, Wan S, Chen S, et al. ANXA6: a key molecular player in cancer progression and drug resistance. Discover Oncol. 2023;14(1):53. doi:10.1007/s12672-023-00662-x

46. Peng D, Fu M, Wang M, et al. Targeting TGF-β signal transduction for fibrosis and cancer therapy. Mol Cancer. 2022;21(1):104. doi:10.1186/s12943-022-01569-x

47. Corona A, Blobe GC. The role of the extracellular matrix protein TGFBI in cancer. Cell Signalling. 2021;84:110028. doi:10.1016/j.cellsig.2021.110028

48. Nielsen NS, Poulsen ET, Lukassen MV, et al. Biochemical mechanisms of aggregation in TGFBI-linked corneal dystrophies. Prog Retinal Eye Res. 2020;77:100843. doi:10.1016/j.preteyeres.2020.100843

49. Liu M, Iosef C, Rao S, et al. Transforming growth factor-induced protein promotes NF-κB-mediated angiogenesis during postnatal lung development. Am J Respir Cell Mol Biol. 2021;64(3):318–330. doi:10.1165/rcmb.2020-0153OC

50. Wang BJ, Chi KP, Shen RL, et al. TGFBI promotes tumor growth and is associated with poor prognosis in oral squamous cell carcinoma. J Cancer. 2019;10(20):4902–4912. doi:10.7150/jca.29958

51. Lee J, Lee J, Sim W, et al. Soluble TGFBI aggravates the malignancy of cholangiocarcinoma through activation of the ITGB1 dependent PPARγ signalling pathway. Cellular Oncol. 2022;45(2):275–291. doi:10.1007/s13402-022-00668-7

52. Lang K, Kahveci S, Bonberg N, et al. TGFBI protein is increased in the urine of patients with high-grade urothelial carcinomas, and promotes cell proliferation and migration. Int J Mol Sci. 2019;20(18). doi:10.3390/ijms20184483

53. Shang D, Song B, Liu Y. Epirubicin suppresses proliferative and metastatic potential by downregulating transforming growth factor-β-induced expression in urothelial carcinoma. Cancer Science. 2018;109(4):980–987. doi:10.1111/cas.13403

54. Chiavarina B, Costanza B, Ronca R, et al. Metastatic colorectal cancer cells maintain the TGFβ program and use TGFBI to fuel angiogenesis. Theranostics. 2021;11(4):1626–1640. doi:10.7150/thno.51507

55. Wen L, Han Z, Du Y. Identification of gene biomarkers and immune cell infiltration characteristics in rectal cancer. J Gastrointestinal Oncol. 2021;12(3):964–980. doi:10.21037/jgo-21-255

56. Lecker LSM, Berlato C, Maniati E, et al. TGFBI production by macrophages contributes to an immunosuppressive microenvironment in ovarian cancer. Cancer Res. 2021;81(22):5706–5719. doi:10.1158/0008-5472.CAN-21-0536

57. Liu SY, Zhu RH, Wang ZT, et al. Landscape of immune microenvironment in epithelial ovarian cancer and establishing risk model by machine learning. J Oncol. 2021;2021:5523749. doi:10.1155/2021/5523749

58. Dekker E, Tanis PJ, Vleugels JLA, et al. Colorectal cancer. Lancet. 2019;394(10207):1467–1480. doi:10.1016/S0140-6736(19)32319-0

59. Ye G, Yang Q, Lei X, et al. Nuclear MYH9-induced CTNNB1 transcription, targeted by staurosporin, promotes gastric cancer cell anoikis resistance and metastasis. Theranostics. 2020;10(17):7545–7560. doi:10.7150/thno.46001

60. Grimshaw KM, Hunter LJ, Yap TA, et al. AT7867 is a potent and oral inhibitor of AKT and p70 S6 kinase that induces pharmacodynamic changes and inhibits human tumor xenograft growth. Mol Cancer Ther. 2010;9(5):1100–1110. doi:10.1158/1535-7163.MCT-09-0986

61. Yeh LT, Lin CW, Lu KH, et al. Niclosamide suppresses migration and invasion of human osteosarcoma cells by repressing TGFBI expression via the ERK signaling pathway. Int J Mol Sci. 2022;23(1):484. doi:10.3390/ijms23010484

62. Maeng YS, Aguilar B, Choi SI, et al. Inhibition of TGFBIp expression reduces lymphangiogenesis and tumor metastasis. Oncogene. 2016;35(2):196–205. doi:10.1038/onc.2015.73

Continue Reading