In this study, four independent cohort samples were included to demonstrate the role of TNRFSF11B in the progression of colon cancer. The transcriptome and clinical characteristics of 514 cases in a colon adenocarcinoma (COAD) dataset were downloaded from the TCGA database (11). Single-cell sequencing data and corresponding single-cell functional status from GSE81861 (12), which included the RNA expression profiles of 290 single CRC cells, were downloaded from CancerSEA (13). The transcriptional levels of TNFRSF11B in different normal tissues were acquired from the Harmonizome database (https://maayanlab.cloud/Harmonizome/) (14). The mRNA expression profiles of thirty-one pairs of colon cancer and normal colon epithelium as counterparts were obtained from samples collected during emergency surgery at the Fujian Medical University Union Hospital from June 2018 to January 2019. The presence of TNFRSF11B in paraffin-embedded colorectal cancer tissues sampled during general surgery was assessed at the Fujian Medical University Union Hospital from January 2013 to December 2017. The study protocol was approved by the Institutional Review Board of the Fujian Medical University Union Hospital, and all patients provided written informed consent. All experiments complied with standard biosecurity and institutional safety procedures. The workflow of this study has been presented in Supplementary Figure 2 .

To determine the immune-related genes (IRGs) and transcription factors (TFs) involved in the progression of colon cancer, analysis of differentially expressed genes between cancer and paired noncarcinoma samples was performed via the R software edgeR package (15). A false discovery rate (FDR) <0.01 and a log₂|fold change| >2 as the criteria values were used to select differential genes. Differentially expressed IRGs and TFs, downloaded from the Immunology Database and Analysis Portal (ImmPort) database (16), were then extracted from all differentially expressed genes. To explore the interactions between IRGs and TFs, a PPI network was constructed based on the STRING database (https://cytoscape.org/). Hub gene analysis of the PPI network was further performed using Cytoscape software version 3.8.2. Functional enrichment analyses were performed via the GO and KEGG pathways.

Total RNA was extracted from 62 colon cancer samples using TRIzol (Invitrogen, Paisley, UK), and mRNA libraries were prepared and sequenced at BGI (HiSequation 2000) (17). Reads were aligned to the human reference genome (GRCh38/hg38) by SOAP2.32. The number of transcripts per million was used to detect gene expression levels. The values of transcripts per million were transformed by log2 and normalized by the Limma R package to remove the batch effect (18). To further elucidate the function of TNFRSF11B in colon cancer, GSEA was performed on a transcriptional dataset in our center for the top 25% and bottom 25% of the expressed genes. C2.cp.kegg.v5.2.symbols.gmt from the KEGG pathway database was selected as the reference gene set. Using a Venn diagram, we determined the common altered signaling pathways between cancer paired noncarcinoma analysis and TNFRSF11B expression analysis strategies (19).

Eighty-six patients diagnosed with colon cancer who underwent radical colectomy from January 2013 to January 2017 during emergency surgery at the Fujian Medical University Union Hospital were included in this study. Differential analysis of the protein expression levels of TNFRSF11B in these colon cancer tissues revealed the correlation between TNFRSF11B and colon cancer tumorigenesis. Anti-osteoprotegerin (TNFRSF11B) monoclonal antibody (ab9986, Abcam, UK) was used at a working concentration of 1:200. We built IHC systematic scores, integrating the intensity and percentage of positive cells, to evaluate the expression level of TNFRSF11B. The staining intensity was separated as follows: 0, no staining; 1, light yellow staining; 2, yellow-brown staining; and 3, deep brown staining. The percentage of positive cells was scored as follows: 0, 0–5%; 1, 6–25%; 2, 26–50%; 3, 51–75%; and 4, >75%. The final score was calculated as follows: positive cell score × staining intensity score. Finally, patients were sorted into strong- and weak-TNFRSF11B expression groups, and differential analysis was performed between TNFRSF11B expression and clinicopathological features. To further detect the expression levels of TNFRSF11B in normal lymph nodes, TNFRSF11B expression in normal lymph nodes was reviewed by using the immunohistochemical (IHC) staining data provided in the Human Protein Atlas (http://www.proteinatlas.org/) (20).

Peripheral blood mononuclear cells (PBMCs) were isolated from healthy volunteers using lymphocyte separation medium (Lonza) according to the manufacturer’s instructions. After centrifugation and purification, lymphocytes were obtained. Cells with 1 µg/ml mouse anti-human CD28 PerCP-Cyanine5.5 mAb (eBioscience #45-0088-42) were washed with PBS and 1640 culture medium twice and then incubated in a 96-well plate precoated with 0.2 µg/ml mouse anti-human CD3 mAb Alexa Fluor^® 488 (eBioscience #53-0037-41) for one day. Soluble TNFRSF11B was added to the culture medium of the test group, while mIgG was added to the control group. Three days later, the cells were analyzed by flow cytometry. The percentages of subgroups including naïve memory CD4⁺ T cells (CCR7⁺CD45RA⁺; CCR7 mAb (3D12), PE-Cyanine7, eBioscience #25-1979-42; mouse anti-human CD45RA BV421, eBioscience #25-1979-42), central memory CD4+ T cells (CCR7⁻CD45RA⁻) and effector memory CD4+ T cells (CCR7⁺CD45RA⁻) were further assessed by flow cytometry.

Our transcriptome data were evaluated by comparing CIBERSORT (21) with the LM22 gene matrix as a reference to identify 22 immune cell subtypes. Microarray probes were replaced with HUGO gene symbols. Genes with multiple probes collapsed to the gene with the highest mean expression. Expression profiles were normalized using fRMA. Using TIMER2.0, we further detected the abundance of immune cells in the microenvironment of TCGA-COAD patient samples and the correlation with TNFRSF11B expression (22).

All graphs and statistical visualization were performed by R software (4.0.1) and Hiplot (https://hiplot.com.cn). Qualitative variables were compared by the χ² test or Fisher’s exact test, and quantitative variables were compared via t-tests. Through Kaplan–Meier analysis, the long-term survival outcome was calculated, and the independent risk factors were identified by the Cox proportional hazards regression model.

A total of 208 immune-related genes (IRGs) ( Figures 1B, E ), including 111 with downregulated expression and 97 with upregulated expression, were identified from the TCGA-COAD database via comparison between tumor and normal tissue (n = 472 and 42, respectively), following the analysis of 2,471 differentially expressed genes ( Figures 1A, D ) (|log₂FC| >2 and P value <0.01). We also identified six downregulated and 16 upregulated differentially expressed genes by the analysis of 318 classical transcription factors (TFs) ( Figures 1C, F ). Functional enrichment analysis demonstrated that these differentially activated IRGs were mostly involved in the cytokine–cytokine receptor interaction KEGG pathway ( Figure 1G ). Furthermore, we revealed twenty-one IRGs closely related to the survival outcome via univariate Cox analysis ( Table 1 ). These survival-related IRGs were enriched in complement activation, the Fc-gamma receptor signaling pathway involved in phagocytosis, the immune response, the B cell receptor signaling pathway, antigen binding, serine-type endopeptidase activity, the extracellular region and the external side of the plasma membrane based on gene ontology (GO) analysis ( Figure 1H ). We also constructed an interaction network of the survival-related IRGs and TFs ( Figure 1I ). MYC, LEF1, KLF4, SALL4, and IRF4 were the five important TFs involved in the regulation of the immune response in colon cancer. TNFRSF11B may correlate with IL17A in the regulation of the colon cancer immune response.

A risk scoring system consisting of eight IRGs (FGFR2, ZC3HAV1L, TNFRSF11B, CD79A, IGHV3-11, IGHV3-21, IGKV2D-30, and IGKV6D-21) was constructed to evaluate the prognosis of colon cancer patients based on multivariate Cox analysis. Patients were separated into two groups, high risk and low risk, according to the risk scoring system. We demonstrated that the high-risk group had worse overall survival outcomes (p = 5.123∗e⁻⁷, Figure 2A ) than the low-risk group, and the area under the curve (AUC) of the receiver operating characteristic (ROC) curve was 0.619 ( Figure 2B ), which supported the potential utility of the risk scoring system. As the risk scores increased, the survival time decreased ( Figure 2C ). We integrated the risk scores with patients’ clinical characteristics, and the risk scoring system was identified as an independent factor via univariate and multivariate Cox analyses ( Figures 2D, E ). We analyzed the relationship between the survival-related IRGs and the characteristics of colon cancer patients and found that only TNFRSF11B was strongly associated with TNM stage and lymph node status. Overexpression of TNFRSF11B obviously correlated with late-stage lymph node metastasis and worse survival outcomes (p = 0.010, p = 0.014, and p = 0.0061, respectively) ( Figures 2F–H ).

To validate the role of TNFRSF11B in the progression of colon cancer, we used immunohistochemistry (IHC) to assess the expression of TNFRSF11B and its correlation with pathological features. In our cohort, 72.09% (62/86) of the colon cancer patients who received radical surgery had TNFRSF11B overexpression, which was associated with significantly shorter overall survival times (p = 0.072, Figures 3A, B ). Patients with high TNFRSF11B expression typically had a late TNM stage (p = 0.067), a high frequency of lymph node (p = 0.029) and lymphovascular invasion (p = 0.007), and a high incidence of pneumonia (p = 0.056) ( Table 2 and Figures 3C, D ). Conversely, the transcriptional expression of TNFRSF11B in normal lymph nodes remained relatively low levels, via data mining in the Harmonizome database ( Supplementary Figure 1A ). Moreover, it’s hard to detect the protein expression of TNFRSF11B in germinal center cells and non-germinal center cells from the Human Protein Atlas database ( Supplementary Figures 1B–D ).

TNFRSF11B, as a self-regulator of M-cell differentiation in colon mucosa, is also involved in the differentiation of osteoclasts in the bone marrow and secondary lymphoid tissues (7, 9, 23). We analyzed the potential role of TNFRSF11B in the differentiation of colon cancer via a GSE81861 array. Increased expression of TNFRSF11B significantly enhanced three cancer-related functional states, namely, metastasis, angiogenesis and inflammation, in colon cancer cells (cor = 0.360, cor = 0.224, and cor = 0283, respectively, all p values <0.001) ( Figures 4A–C ). TNFRSF11B expression was also enriched in the terminal differentiation of colon cancer cells through pseudotime analysis ( Figure 4 ). Thus, TNFRSF11B may drive the progression of early-stage colon cells to cells with full tumorigenicity.

To determine which genes and pathways are responsible for TNFRSF11B expression in colon cancer, we performed transcriptome-wide RNA-sequencing analysis of colon cancer tissues and paired colon epithelium. Two independent gene set enrichment analyses (GSEAs), namely, tumor-only and comparisons between tumors and paired colon epithelium, were applied to identify the involved signaling pathways and genes. Six pathways, including RNA polymerase, one carbon pool by folate, insulin signaling pathway, beta alanine metabolism, mTOR signaling pathway and lysine degradation, were inhibited by the upregulation of TNFRSF11B expression. Pathogenic Escherichia coli (E. coli) infection, drug metabolism, other enzymes and Vibrio cholerae infection were the other three pathways activated by the increased expression of TNFRSF11B. Consistently, only the pathogenic E. coli infection pathway was significantly elevated with the increasing level of TNFRSF11B among tumors Normalized Enrichment Score (NES) = 1.814, p = 0.004) and in the comparison of tumors with paired colon epithelium (NES = 1.523, p = 0.027) ( Figures 5A, C, D ). TNFRSF11B also presented stronger transcriptional activity in tumor tissues than in paired colon epithelium (p <0.001, Figure 5B ). Previous studies have revealed that pathogenic E. coli infection induces the overexpression of CEACAM6 in colon cancer, which, in turn, facilitates pathogenic E. coli adhesion (24). We demonstrated that CEACAM6 expression was significantly upregulated in cancer tissues that had high TNFRSF11B expression ( Figure 5E ). In addition, only six genes (KRT18, ARPC5L, ACTG1, ARPC2, EZR, and YWHAZ) related to pathogenic E. coli infection were simultaneously increased in both GSEA studies; these genes are mostly involved in focal adhesion, as determined via GO analysis ( Figures 5F, G ). Several KEGG pathways, such as regulation of actin cytoskeleton, shigellosis, bacterial invasion of epithelial cells, and Salmonella infection, were also affected by these six genes ( Figure 5H ). Next, we attempted to determine the interaction between TNFRSF11B expression and these six genes via hub gene analysis. TNFRSF11B may affect the transcriptional activity of YWHAZ to regulate pathogenic E. coli infection ( Figure 5I ), which was further validated in the TCGA-COAD database, in which increased TNFRSF11B expression correlated with increased TWHAZ expression (cor = 0.185, p <0.001) ( Figure 5J ).

To analyze the potential role of TNFRSF11B in the regulation of the immune response, we determined the different subsets of immunocytes between different TNFRSF11B-expressing colon cancer tissues in our cohort via CIBERSORT algorithms. Only activated memory CD4⁺ T cells (p = 0.017) were significantly decreased in the high TNFRSF11B expression group ( Figure 6A ), which was supported by the findings in the TCGA-COAD dataset. Using the TIMER 2.0 database, we also demonstrated that TNFRSF11B negatively affected the infiltration of activated memory CD4⁺ T cells into the colon cancer microenvironment (cor = −0.158, p = 0.008) ( Figure 6B ). We also confirmed the suppressive role of TNFRSF11B in the regulation of central memory CD4⁺ T cells and effector memory CD4⁺ T cells via the addition of soluble TNFRSF11B to PBMC culture medium. TNFRSF11B dramatically downregulated the immune activity of central memory CD4⁺ T cells and effector memory CD4⁺ T cells (all p <0.001) ( Figures 6C, D ).