Constructing a chemokine-based model and identifying CCL17 as a core biomarker associated with immune infiltrates in thyroid cancer

Introduction

Thyroid cancer (THCA) is a prevalent malignancy with an escalating global incidence (1). THCA is categorized as papillary, follicular, poorly differentiated, and anaplastic tumors. Conventional therapeutic approaches for THCA encompass surgical resection, radiation therapy, targeted therapy, and radioactive iodine therapy, tailored to the specific tumor stage and individual patient characteristics (2,3). Nevertheless, despite remarkable progress in THCA treatment, prognostic challenges persist, as some patients may encounter relapse or metastasis (4). Moreover, prognosis assessment in THCA conventionally relies on clinical features such as age, gender, disease duration, tumor size, and lymph node metastasis. However, traditional clinical models have limitations in precisely predicting the THCA prognosis. Therefore, there has been considerable interest in constructing a prognostic model for THCA based on risk prediction as well as identifying a potent biomarker as a therapeutic target.

Chemokines, a class of protein or peptide molecules, play a pivotal role in the progression and metastasis of cancer (5). They facilitate the directed migration and infiltration of cells by binding to chemokine receptors on the cell surface (6,7). The interplay between chemokines and tumors manifests in two main aspects. Firstly, chemokines can stimulate the proliferation and growth of THCA cells. Studies have observed significant upregulation of specific chemokines, such as tumor necrosis factor-α and interleukin-6, in THCA tissues (8-10). These chemokines activate multiple signaling pathways, such as JAK/STAT and MAPK, thereby promoting cell cycle progression, amplifying the expression of proliferation genes, and facilitating the growth and proliferation of THCA cells (11-13). Secondly, chemokines are intricately linked to tumor cell infiltration and metastasis. Research has demonstrated that chemokine-induced cell migration and infiltration are pivotal steps in the dissemination of THCA cells from the primary site to distant organs (14-18). Certain chemokines, including tumor-associated factor β (TGF-β) and vascular endothelial growth factor, enhance the migratory capacity of tumor cells and promote the formation of new blood vessels, which provide nourishment and pathways for tumor metastasis (19-21). Furthermore, chemokine-related genes (CRGs) have been identified as valuable prognostic indicators in numerous cancers (22-24). However, the role of CRGs in the prognosis of THCA remains unclear.

Thus, the primary aim of this study is to develop a precise risk prediction model and identify a core biomarker for predicting prognosis in THCA, leveraging the distinct characteristics of THCA. By integrating molecular biology insights, our objective is to enhance the accuracy of THCA prognosis prediction and equip physicians with robust scientific support for treatment decision-making, improving the survival rates and quality of life for individuals affected by THCA. We present this article in accordance with the TRIPOD reporting checklist (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2024-2654/rc).

Methods

Data collection

The messenger ribonucleic acid (mRNA) expression and clinical data of THCA were downloaded from The Cancer Genome Atlas (TCGA) database. Inclusion criteria: with matched expression and clinical data; with complete survival information. Finally, there were 512 tumor samples: 504 papillary thyroid cancer (PTC) samples, 8 samples of other THCA types, and 97 normal samples. The chemokine-related genes were determined in GeneCard (https://www.genecards.org/) (25). This study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments.

Identification of differentially expressed chemokine-related genes (DECRGs)

The limma package in R was utilized to identify differentially expressed genes (DEGs) between THCA and normal samples, with the threshold values as fold change (FC) >1.5 and P<0.05. Then, the DEGs intersected with chemokine-related genes were considered as DECRGs.

Construction and evaluation of a risk model

Following the identification of the DECRGs, univariate Cox regression analysis was employed to identify genes closely associated with the prognosis of THCA. Subsequently, a multivariate Cox regression analysis was conducted to determine the hub genes and their corresponding coefficients. Based on the multivariate Cox regression results, a risk model was constructed and a risk score was calculated. The patients were subsequently stratified into low and high-risk groups using the optimal cutoff value of risk score determined by the maxstat package. Finally, the stability and suitability of the risk model were validated through Kaplan-Meier (K-M) and receiver operating characteristic (ROC) curve analyses, utilizing the survival and qROC packages, respectively.

Gene set variation analysis (GSVA)

GSVA package was used for GSVA analysis of the differential marker genes between low- and high-risk groups. All reference gene sets were obtained in the Molecular Signatures Database (MSigDB, https://www.gseamsigdb.org/gsea/downloads.jsp).

Immune characteristics in the two risk groups

To quantify immune cell abundance, we employed a CIBERSORT algorithm, a powerful deconvolution tool for calculating the proportions of 22 immune cell types in THCA samples with the CIBERSORT package. To assess the potential effectiveness of immune checkpoint blockade therapy, we utilized the Tumor Immune Dysfunction and Exclusion (TIDE) computational approach (26). Unlike traditional biomarkers such as tumor mutation burden or PD-L1 expression, TIDE is a more reliable predictor of response to immunotherapy (27). We first normalized the raw expression profiles using z-scores. Subsequently, we calculated TIDE scores, microsatellite instability (MSI), and Dysfunction and Exclusion scores through the TIDE website (http://tide.dfci.harvard.edu/). A lower TIDE score indicates a decreased likelihood of immune escape, suggesting a more favorable response to immunotherapy.

Cell culture

Healthy human thyroid cell (Htori-3, CL-0817) and THCA cell (TPC-1, CL-0643) were acquired from Procell (Wuhan, China). The cell lines were inoculated in the Rosewell Park Memorial Institute-1640 medium (PM150110; Procell) supplemented with 40% fetal bovine serum (FBS, 164210-50; Procell) and 1% penicillin-streptomycin solution (PB180120; Procell) at 37 ℃ in 5% CO₂. The media were changed every 2 or 3 days.

Small interfering RNA (siRNA) transfection

The siRNA targeting CCL17 and the scramble siRNA were purchased from GenePharma Co., Ltd. (Wuhan, China). Transfection was performed using Lipofectamine™ RNAiMAX reagent (Thermo Fisher Scientific, San Jose, CA, USA) according to the manufacturer’s protocol. The siRNA sequences were as follows: scramble siRNA: GAGUGAAGAAUGCAGUUAA; CCL17 siRNA (antisense): UAACUGCAUUCUUCACUC.

Vector-CCL17 transfection

A CCL17 overexpression plasmid was constructed using the pcDNA3.1 vector (GeneCreate Bioengineering Co., Ltd., Wuhan, China) (Figure S1). TPC-1 cells at 75–85% confluence were collected and counted, and 1×10⁶ cells were transfected with 10 nM plasmid using Lipofectamine 2000 (Thermo Fisher Scientific) according to the manufacturer’s instructions. Cells transfected with the empty pcDNA3.1 vector served as the negative control. All subsequent experiments were performed using cells harvested 48 hours post-transfection.

Quantitative real-time polymerase chain reaction (qRT-PCR)

The Trizol (B511311-0500; Sangon, Shanghai, China) was utilized for extracting total RNA from cells. After that, Cary 3500 UV-Vis Spectrophotometer (Agilent Technologies, Inc., Santa Clara, CA, USA) was employed to assess the A260/A230 and A260/A280 ratios. Following this, reverse transcription of RNA samples and qRT-PCR were conducted using the One Step RT-qPCR Kit (B639277-0100; Sangon) within the Mx3000P qPCR system (Agilent Technologies, Inc.). The expression levels were determined utilizing the 2^−ΔΔCt method, with GAPDH serving as the internal reference.

PCR reaction conditions were: 5 minutes at 50 ℃, 3 minutes at 95 ℃, then 10 seconds at 95 ℃ and 30 seconds at 60 ℃ for 40 cycles. The primer sequences are shown in Table 1.

Western blot (WB)

WB was performed to assess the expression levels of target proteins. Briefly, treated cells were harvested and lysed on ice for 30 minutes using radioimmunoprecipitation assay (RIPA) lysis buffer supplemented with protease and phosphatase inhibitors. The lysates were then centrifuged at 12,000 ×g for 15 minutes at 4 ℃, and the supernatants were collected. Protein concentrations were determined using a BCA protein assay kit. sodium dodecyl sulfate–polyacrylamide gel electrophoresis (SDS-PAGE) was used to separate equal amounts of protein (40 µg per lane), and the proteins were subsequently transferred onto polyvinylidene difluoride (PVDF) membranes (0.45 µm pore size). Following transfer, the membranes were blocked with 5% non-fat milk in tris-buffered saline with Tween 20 (TBST) buffer for 1 hour at room temperature. The membranes were then incubated overnight at 4 ℃ with primary antibodies against the target proteins (1:1,000 dilution). On the following day, membranes were washed three times with TBST (10 minutes each) and incubated with horseradish peroxidase (HRP)-conjugated secondary antibodies (1:10,000 dilution) for 1 hour at room temperature. After thorough washing, signal detection was carried out using an enhanced chemiluminescence (ECL) substrate, and the signals were visualized and captured using a gel imaging system. Band intensities were quantified using Quantity One software, and normalized to GAPDH as a loading control. The antibodies used in this study included: (I) anti-TARC/CCL17 (ab182793, Abcam, Cambridge, UK); (II) hospho-JAK2 (Y1007 + Y1008); recombinant rabbit monoclonal antibody (SY24-03) (ET1607-34, Huabio, Hangzhou, China); (3) JAK2 rabbit polyclonal antibody (RT1343, Huabio); (IV) phospho-STAT3 (S727) recombinant rabbit monoclonal antibody (SY24-09) (ET1607-39, Huabio); (V) STAT3 rabbit polyclonal antibody (RT1591, Huabio); (VI) GAPDH polyclonal antibody (10494-1-AP, Proteintech, Chicago, USA).

Colony formation assay

Colony formation assay was performed to evaluate the long-term proliferative capacity and clonogenic potential of the cells. Log-phase single-cell suspensions were seeded at a low density (1×10³ cells/well) into 6-well plates. After cell attachment, appropriate treatments were applied [transfection with siRNA, siRNA negative control (si-NC), overexpression RNA (oeRNA), or overexpression negative control (oe-NC)], and the cells were cultured under standard conditions (37 ℃, 5% CO₂) for 14 days, with fresh medium replaced every other day. At the end of the incubation period, the culture medium was removed, and cells were gently washed twice with phosphate-buffered saline (PBS). Cells were then fixed with methanol for 15 minutes. After discarding the fixative and allowing the plates to air-dry, the cells were stained with 0.1% crystal violet solution (dissolved in PBS or water) at room temperature for 30 minutes. The stain was carefully aspirated, and excess dye was gently rinsed away with running water. Plates were air-dried before imaging. Colonies with a diameter greater than 100 µm (typically defined as ≥50 cells per colony) were counted using ImageJ software.

Transwell cell invasion assay

Cell migration ability was assessed using Transwell chambers with 8 µm pore-size polycarbonate membranes. The Transwell inserts were placed in 24-well plates, with the lower chambers filled with 600 µL of complete medium containing 10% FBS to serve as a chemoattractant. Cells in the logarithmic growth phase were serum-starved for 4 hours, and then resuspended as single-cell suspensions in serum-free medium. A total of 2×10⁵ cells in 200 µL of serum-free medium were added to the upper chamber of each insert. The plates were incubated at 37 ℃ in a humidified atmosphere with 5% CO₂ for 24 hours. After incubation, the inserts were carefully removed, and non-migrated cells on the upper surface of the membrane were gently wiped off with a cotton swab. The cells that had migrated to the lower surface of the membrane were fixed in methanol for 15 minutes and then stained with 0.1% crystal violet for 20 minutes. After washing with PBS, five random fields per insert were photographed and counted under a light microscope at 200× magnification. Each treatment group was performed in triplicate.

Wound healing assay

A wound healing assay was conducted to evaluate the horizontal migratory capacity of cells. Cells in the logarithmic growth phase were seeded at high density into 6-well plates and cultured for 24–48 hours until a confluent monolayer (~100%) was formed. A sterile 200 µL pipette tip was used to create three parallel scratches of uniform width across the cell monolayer in each well. Detached cells and debris were gently removed by washing the wells three times with PBS.

Subsequently, low-serum medium containing 1% FBS was added to minimize the influence of cell proliferation on migration. The plates were incubated at 37 ℃ in a humidified atmosphere with 5% CO₂. Images of the wound areas were captured at 0 hours and 24 hours post-scratch using an inverted microscope equipped with a live-cell imaging system. The same fields were recorded at each time point for comparison.

The Human Protein Atlas (HPA) database analysis

The HPA, a Swedish-based program, serves as a valuable resource for accessing public and freely available protein information across approximately 36 tissues. One of the key features of the HPA database is the ability to delineate protein distribution in tissues through immunohistochemistry (IHC) staining images. We utilized IHC staining images from the HPA database to compare the protein levels of ACKR3, CCL2, and CCL17 between samples of THCA and normal thyroid tissues.

The prognostic value of the hub genes

K-M curves were plotted to assess the survival differences between the two gene groups. ROC curves were plotted with the area under the curve (AUC) calculated to analyze the value of the hub genes and risk score in predicting 1-, 3-, and 5-year survival of patients with THCA. The Delong test was employed to examine the statistical differences between the AUCs.

Clinical relevance of CCL17 and mechanism exploration

Due to the superiority of CCL17 among the three hub genes, we further explored the clinical value of CCL17 by ROC and decision curve analysis (DCA) using the qROC and ggDCA packages, respectively.

To further reveal the biological function of CCL17, we divided the samples into low- and high-CCL17 groups based on the median value to identify significant DEGs with the following threshold: |log₂ FC| ≥1, P<0.05. Then, the Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis on the identified significant DEGs was conducted using the clusterProfiler package. Besides, we conducted the gene set enrichment analysis (GSEA) to investigate the potential pathways in which CCL17 might be involved in THCA.

Connectivity map (CMap) analysis

The CMap (https://clue.io/) provides a wealth of information on small molecule drugs in the “query” module by selecting the query parameters “gene expression (L1000)”. The top 50 up- and down-regulated DEGs between the low- and high-CCL17 groups were submitted for query. The molecular compounds, mechanism of activities, and connectivity score were obtained. A connectivity score from −1 to 1 reflects the proximity between expression profiles. A drug with a negative score indicates a potential therapeutic molecule. We selected five drugs with the top five connectivity scores as potential drugs.

Molecular docking analysis

The RCSB Protein Data Bank (PDB) and traditional Chinese medicine systems pharmacology database and analysis platform (TCMSP) databases were used to obtain the 3D structure files of key genes in pdb format and molecule compounds in mol2 format, respectively. Then, molecular docking of molecule compounds and hub genes was performed using AutoDock software (version 4.2.6). Before docking, water molecules were removed from the receptor file, and hydrogen was added. The molecular docking results were visualized using PyMOL software (version 2.1.0).

Statistical analysis

All statistical analyses were conducted using R software version 4.3.1 and SPSS software version 23.0. Categorical data were expressed as numbers (percentages) and group comparisons were conducted by Chi-squared test. Continuous variables with non-normal distribution were presented as median (interquartile range); the differences were compared using the Mann-Whitney test or Kruskal-Wallis test for two groups and at least three groups, respectively. The Spearman correlation test was carried out for correlation analysis. Generalized additive models (GAMs) and restricted cubic spline (RCS) analyses were used to evaluate the dose-response relationship between dependent and independent variables. Statistical significance was defined as a P value less than 0.05.

Results

Construction and evaluation of the risk model

As shown in Figure 1A, there are 3,234 upregulated genes and 4,245 downregulated genes in THCA samples compared with normal samples. Besides, the Venn plot exhibited that 27 chemokines had significantly differential expression between THCA and normal samples (Figure 1B). Upon the univariate and multivariate Cox analysis, the genes with independent prediction values were screened out and the risk model was constructed using the expression of hub genes and corresponding coefficients as follows: Risk score = 0.023*Exp (ACKR3) + 0.028*Exp (CCL2) + 0.019*Exp (CCL17) (Figure 2A,2B). The PCA results showed that the high-risk group and low-risk group could be separated based on the risk score (Figure 2C). The K-M curve displayed that the high-risk group was closely related to the poor prognosis (P<0.05) (Figure 2D). The ROC curve showed a good prediction ability of the risk score on patient prognosis with AUCs of 0.86, 0.76 and 0.68 at 1, 3, and 5 years, respectively (Figure 2E).

Figure 1 The identification of DECRGs. (A) The upregulated and downregulated genes in THCA compared with normal samples. (B) The intersection of DEGs and chemokine genes. DECRGs, differentially expressed chemokine-related genes; DEGs, differentiated expressed genes; FC, fold change; FDR, false discovery rate; THCA, thyroid cancer.

Figure 2 The construction and evaluation of a risk model. The univariate (A) and multivariate (B) Cox analysis of DECRGs. (C) The PCA curve exhibited the characteristics of low- and high-risk groups. (D) The survival analysis of low- and high-risk groups. (E) The ROC curve exhibited a great predictive value of the risk model for THCA prognosis. AUC, the area under the curve; CI, confidence interval; DECRGs, differentially expressed chemokine-related genes; PC, principal component; PCA, principal component analysis; ROC, receiver operating characteristic; THCA, thyroid cancer.

The functional enrichment analysis between low- and high-risk groups

To explore the characteristics of low- and high-risk groups, the GSVA was used to analyze the differentially enriched functions. The results exhibited that the high-risk group was highly enriched in angiogenesis, epithelial-mesenchymal transition (EMT), TNF signaling via NFkB and immune-related signal pathways, including IL2-STAT5 signal pathway and IL6-JAK-STAT3 signal pathways (Figure 3).

Figure 3 The differentially enriched signal pathways between high and low-risk groups. GSVA, gene set variation analysis; KEGG, Kyoto Encyclopedia of Genes and Genomes.

The immunological analysis

Then, the immune infiltration analysis was first used to assess the immune microenvironment. From Figure 4A, the infiltration levels of CD8 T cells and M2 macrophages were significantly decreased in the high-risk group, but the infiltration levels of follicular helper T (Tfh) cells, M1 macrophages, and resting dendritic cells (DCs) were significantly increased. As shown in Figure 4B, ACKR3 was negatively related to the DCs resting. CCL2 was positively related to the follicular helper T cells, M1 macrophages, and resting DCs but negatively related to M2 macrophages. CCL17 was negatively related to the CD8 T cells and M2 macrophages but positively related to the follicular helper T cells and resting DCs. Besides, the TIDE score in the high-risk group was higher than that in the low-risk group (Figure 4C). Although there was no significant difference in MSI between the two groups (Figure 4D), the Dysfunction and Exclusion scores in the high-risk group were higher than the low-risk group (Figure 4E,4F). These suggested that the high-risk group had a higher potential for immune escape and worse response to immunotherapy.

Figure 4 The immune characteristics in high- and low-risk groups. (A) The immune infiltration levels in two risk groups. (B) The correlation between the hub genes and immune cells. The TIDE score (C), MSI score (D), dysfunction score (E), and exclusion score (F) in two risk groups. *, P<0.05; **, P<0.01; ***, P<0.001. MSI, microsatellite instability; TIDE, tumor immune dysfunction and exclusion.

The correlation between hub genes and immune checkpoint genes (ICGs)

We further assess the function of hub genes in immune escape. The results exhibited that ACKR3 was not related to ICGs (Figure 5A), and CCL2 was positively related to CD209 (Figure 5B). Additionally, CCL17 was positively related to the CD274, HLA-DRB1, CTLA4, HLA-DMA, HLA-C, SIRPA, HLA-DQB1, HLA-DOB, HLA-DRA, TNFSF18, HLA-F, HLA-B but negatively related to BTNL9 (Figure 5C).

Figure 5 The correlation between hub genes and immune checkpoints. (A) ACKR3, (B) CCL2, (C) CCL17.

The mRNA and protein expression of hub genes

Next, to explore the potential mechanism of hub genes, we first calculated the difference in hub gene expression in different groups using TCGA data. ACKR3 and CCL2 expressions were downregulated, while CCL17 was upregulated in the tumor compared with normal samples (Figure 6A). In addition, three hub genes were all upregulated in the high-risk group in comparison with the low-risk group (Figure 6B). The qRT-PCR results showed that CCL17 was upregulated, but ACKR3 and CCL2 were downregulated in TPC-1, which was consistent with bioinformatics results (Figure 6C). The IHC results exhibited that the protein levels of ACKR3 and CCL2 were decreased in tumors compared with normal samples, while CCL17 was increased (Figure 6D-6F).

Figure 6 The mRNA and protein expression of hub genes. (A) The mRNA expression of hub genes in normal and THCA samples. (B) The mRNA expression of hub genes in low- and high-risk groups. (C) The expression of hub genes in TPC-1 and Htori-3 cell lines. Representative immunohistochemistry images of (D) ACKR3, (E) CCL2, and (F) CCL17 in both normal and THCA tissues sourced from the Human Protein Atlas database (https://www.proteinatlas.org/). Image credit goes to the Human Protein Atlas. The links to the individual normal and tumor tissues of each protein are provided for ACKR3 (https://www.proteinatlas.org/ENSG00000144476-ACKR3/tissue/thyroid+gland#img, https://www.proteinatlas.org/ENSG00000144476-ACKR3/cancer/thyroid+cancer#img), CCL2 (https://www.proteinatlas.org/ENSG00000108691-CCL2/tissue/thyroid+gland#img, https://www.proteinatlas.org/ENSG00000108691-CCL2/cancer/thyroid+cancer#img), and CCL17 (https://www.proteinatlas.org/ENSG00000102970-CCL17/tissue/thyroid+gland#img, https://www.proteinatlas.org/ENSG00000102970-CCL17/cancer/thyroid+cancer#img), respectively. Scale bar, 100 µm. **, P<0.01; ***, P<0.001. mRNA, messenger ribonucleic acid; THCA, thyroid cancer.

The prognosis value of hub genes

Then, we assessed the prognostic value of hub genes in THCA. The K-M curve showed that the high expression of CCL17 (P<0.001) and ACKR3 (P=0.005) was closely related to the poor prognosis of THCA patients (Figure 7A,7B), but CCL2 was related to the poor prognosis with P=0.01 (Figure 7C). The ROC curve exhibited that ACKR3 has a higher AUC value of 0.706 than CCL17 (AUC =0.580) and CCL2 (AUC =0.552), which was slightly lower than that of the risk score (AUC =0.740) (Figure 7D). However, the Delong test result showed that there was no significant difference in AUCs of ACKR3, CCL17 and the risk score (Table 2). These suggested the ACKR3 and CCL17 had similar clinical value to the risk score. Of note, ACKR3 was not related to ICGs, while CCL17 had an intimate relationship with multiple ICGs and immune cells. In addition, the sensitivity, specificity, Youden index, optimal threshold, and accuracy of CCL17 for predicting poor prognosis were 0.500, 0.790, 0.290, 3.112, and 0.968, respectively, which helps reduce false positives and avoid unnecessary treatment for low-risk patients (Table S1). This triggers us to further analyze the significance of CCL17 in THCA.

Figure 7 The prognosis value of hub genes. The survival analysis of CCL17 (A), ACKR3 (B), and CCL2 (C). The red line meant the high expression of hub genes, and the blue line meant the low expression of hub genes. The ROC curve of the hub genes and risk score (D). AUC, the area under the curve; ROC, receiver operating characteristic.

The relationship between CCL17 and clinical parameters

Participants were divided into low-CCL17 and high-CCL17 groups according to the median value of CCL17. As shown in Table 3, pathological N and stage were differentially distributed in the two CCL17 groups, with more N1 cases (63.374% vs. 35.616%) and stage IV cases (16.016% vs. 6.299%) in the high-CCL17 group (P<0.001). Besides, PTC-columnar cell type accounted for 12.891% of the cases in the high-CCL17 group and 1.563% in the low-CCL17 group (P<0.05). Neoplasm depth seems different between the two groups. More positive lymph nodes were detected in patients in the high-CCL17 group (P<0.001). Other parameters were not significantly different between the two groups.

Further, we analyzed the association of CCL17 as a continuous variable with pathological N, stage, and histological type. A higher CCL17 expression was observed in the N1 group (P<0.001) (Figure 8A). Compared with the stage I group, CCL17 expression was slightly downregulated in stage II and tended to increase in stage III and stage IV (P<0.001) (Figure 8B). Similarly, CCL17 expression was significantly higher in the PTC-columnar cell group and others than the PTC-classical and PTC-follicular group (Figure 8C). In addition, we found that CCL17 expression was positively correlated with positive lymph nodes, but had no significant relationship with neoplasm depth (Figure 8D,8E). Since the close relationship between CCL17 and histologic type, which was associated with a differentiated state (28), we further analyzed the association of CCL17 with thyroid differentiated score (TDS). TDS was decreased with the increase of CCL17 expression (r=−0.31) (Figure 8F).

Figure 8 The association of CCL17 expression with clinical parameters. (A) Pathologic N. (B) Stage. (C) Histologic type. The association of CCL17 expression with positive lymph nodes (D), neoplasm depth (E), and TDS (F). N, node; PTC-1, papillary thyroid cancer-classical; PTC-2, PTC-follicular; PTC-3, PTC-columnar cell; TDS, thyroid differentiated score.

The value of CCL17 in predicting lymph node metastasis and TDS

Moreover, we investigated the value of CCL17 in predicting lymph node metastasis and TDS. First, we had an initial overlook of the dose-response relationship between CCL17 expression and positive lymph nodes. Due to the skewed distribution of CCL17, we did a log transformation of the CCL17 expression. The GAM result showed that log CCL17 was positively related to positive lymph nodes, but negatively related to TDS (Figure 9A,9B). When taking pathologic N and TDS group (taking TDS median value as the cutoff) as outcomes, similar trends were observed in the RCS results (Figure 9C,9D).

Figure 9 The dose-response relationship between CCL17 expression and positive lymph nodes, TDS. The GAM revealed CCL17 positively related to positive lymph nodes (A) but negatively related to TDS (B). RCS curve showed the association of CCL17 expression with pathologic N (C) and TDS group (D). (A,B) The two dotted lines in the graphs represent the upper and lower confidence bounds of 95% CI, and the blue solid line represents the relationship between log CCL17 and positive lymph nodes; (C,D) the red solid line represents the relationship between log CCL17 and the odds ratio of TDS, and the red shaded area indicates the 95% CI. CI, confidence interval; GAM, generalized additive models; N, node; RCS, restricted cubic spline; TDS, thyroid differentiated score.

Subsequently, we evaluated the predictive value of CCL17 for pathological N and TDS group. The AUCs of CCL17 for predicting pathological N and TDS groups were 0.696 and 0.732, respectively (Figure 10A,10B). As shown in Figure 10C, when the threshold was at about 0.4–0.9, CCL17 had a clinical net benefit for predicting the pathological N. When the threshold was at about 0.4–0.7, CCL17 had a clinical net benefit for predicting the TDS group (Figure 10D).

Figure 10 The clinical value of CCL17 in predicting pathologic N and TDS group. ROC curve of CCL17 expression for pathologic N (A) and TDS group (B). The DCA of CCL17 expression for pathologic N (C) and TDS group (D). AUC, the area under the curve; DCA, decision curve analysis; N, node; ROC, receiver operating characteristic; TDS, thyroid differentiated score.

CCL17 was related to the JAK-STAT pathway

We have preliminarily explored the possible function of CCL17 in regulating immune cells and demonstrated the clinical relevance of CCL17 in tumorigenesis and development. We next identified the DEGs between low- and high-CCL17 groups for KEGG pathway enrichment analysis. There were 1,304 upregulated genes and 436 downregulated genes meeting the selection criterion (Figure 11A). The identified DEGs were mainly involved in cytokine-cytokine receptor interaction, cell adhesion molecules, chemokine signaling pathway, phagosome, and JAK-STAT signaling pathway (Figure 11B). To further explore the biological function of CCL17 in THCA, GSEA analysis was conducted and the results showed that CCL17 was enriched in apoptosis, JAK-STAT signaling pathway, and ECM receptor interaction (Figure 11C-11E).

Figure 11 Enrichment analysis. (A) The DEGs between low- and high-CCL17 groups visualized by the volcano plot. (B) The KEGG pathway enrichment analysis of the DEGs between low- and high-CCL17 expression groups. (C-E) The enriched signal pathway of CCL17. DEGs, differentially expressed genes; KEGG, Kyoto Encyclopedia of Genes and Genomes.

Combined with the GSEA results, the JAK-STAT signaling pathway was the shared pathway. Therefore, CCL17 expression might promote lymph node metastasis and tumor progression by regulating the JAK-STAT pathway.

CCL17 affects the malignant phenotype of THCA cells via the JAK-STAT signaling pathway

The role of CCL17 in regulating the malignant phenotype of THCA cells and its influence on the JAK-STAT signaling pathway were explored. We first established TPC-1 cell lines with CCL17 knockdown (si-CCL17) and overexpression (oe-CCL17). Cells transfected with si-NC and oe-NC served as corresponding negative controls. WB analysis confirmed the successful construction of both si-CCL17 and oe-CCL17 cell lines (Figure 12A). In addition, we also verified that overexpression or silencing of CCL17 was successful at the level of mRNA (Figure 12A). Transwell assay results showed that the number of invading cells was significantly decreased in the si-CCL17 group compared to the si-NC group (P=0.01), whereas it was significantly increased in the oe-CCL17 group compared to the oe-NC group (P=0.019, Figure 12B). Colony formation assay results showed that the number of proliferative colonies was significantly reduced in the si-CCL17 group compared to the si-NC group (P<0.001), while it was significantly increased in the oe-CCL17 group compared to the oe-NC group (P=0.01, Figure 12C). Wound healing assay results demonstrated that the migration rate was significantly decreased in the si-CCL17 group relative to the si-NC group (P=0.001), whereas it was significantly increased in the oe-CCL17 group compared to the oe-NC group (P=0.003, Figure 12D). WB analysis revealed that the phosphorylation levels of JAK2 and STAT3 were significantly reduced in the si-CCL17 group compared to the si-NC group (P<0.001). In contrast, phosphorylation of JAK2 and STAT3 was significantly increased in the oe-CCL17 group compared to the oe-NC group (P<0.001, Figure 12E). These findings indicate that CCL17 overexpression significantly promotes the malignant phenotype of THCA tumor cells. Furthermore, high CCL17 expression may activate the JAK-STAT signaling pathway.

Figure 12 The effect of CCL17 on the proliferation, migration, invasion and JAK-STAT signaling pathway in PTC cells. (A) The efficiency of CCL17 protein knockout or overexpression. (B) Transwell analysis of the invasion. (C) Clone formation assay of cell proliferation. (D)Wound healing assay of cell migration. (E) WB analysis of JAK and STAT3’s phosphorylation level. Crystal violet staining was used in both the Transwell cell invasion assay and the plate cloning assay. The scales of the figure of the Transwell cell invasion assay and wound healing assay are 1:200 µm and 1:500 µm, respectively. oe-NC, overexpression negative control; PTC, papillary thyroid cancer; si-NC, siRNA negative control; WB, western blot.

To investigate whether CCL17 enhances tumor cell malignancy through JAK-STAT pathway activation, we treated oe-CCL17 TPC-1 cells with JAK inhibitors (oe-CCL17 + JAK inhibitor), using oe-CCL17 alone as the control. WB analysis revealed that JAK inhibition substantially reduced JAK and STAT3 phosphorylation compared to the control (P<0.001, Figure 13A). The Transwell assay demonstrated significantly fewer invasive cells in the JAK inhibitor-treated group versus controls (P<0.001, Figure 13B). Colony formation assays showed markedly reduced proliferative capacity in inhibitor-treated cells (P=0.01, Figure 13C). Wound healing assay indicated slower migration in oe-CCL17 + JAK inhibitor cells than in oe-CCL17 controls (P=0.006, Figure 13D).

Figure 13 CCL17 promotes malignant characteristics of THCA cells by activating the JAK-STAT signaling pathway. (A) The expression of proteins in the JAK-STAT pathway using WB. (B) Transwell analysis of the invasion. (C) Clone formation assay of cell proliferation. (D) Wound healing assay of cell migration. Crystal violet staining was used in both the Transwell cell invasion assay and the plate cloning assay. The scales of the figure of the Transwell cell invasion assay and wound healing assay are 1:200 µm and 1:500 µm, respectively. THCA, thyroid cancer; WB, western blot.

Screening of small-molecule drugs

A CMap analysis was conducted to predict potential therapeutic targets for THCA by uploading the top 50 CCL17-associated up-regulated and down-regulated genes. Five target drugs (irsogladine, rupatadine, PYM-50028, AZD-3514, and TG-101348) with the top five norm_connectivity scores were screened and listed in Table 4. Then, molecular docking was performed to analyze the binding effect of CCL17 and the five target drugs. Affinity <−5.0 kcal/mol indicates that the molecule compound has a good docking result with the protein. Table 5 shows the binding effect between CCL17 and five target drugs <−5.0 kcal/mol, suggesting good affinity. The optimal docking posture of the receptor and ligand is presented in Figure 14A-14E.

Figure 14 Molecular docking analysis. The optimal docking posture of CCL17 and isogladine (A), Rupatadine (B), PYM-50028 (C), AZD-3514 (D), and TG-101348 (E).

Moreover, we established three experimental groups to assess TG-101348’s impact on THCA: a TPC-1 control group treated with 1 µM dimethyl sulfoxide, a TPC-1 group treated with 1 µM TG-101384 (TG-101348), and an oe-CCL17 TPC-1 group treated with 1 µM TG-101384 (TG-101348 + oe-CCL17). WB results showed that the expression of CCL17 mRNA and protein in the TG-101348 group was significantly lower than that in the control group (P<0.001). Compared with the TG-101348 group, the mRNA and protein levels of CCL17 were increased in the TG-101348 + oe-CCL17 group (P=0.005) (Figure 15A). The results of the Transwell assay showed that the number of invasive cells in the TG-101348 group was significantly lower than that in the control group (P<0.001). Compared with the TG-101348 group, the number of invasive cells in the TG-101348 + oe-CCL17 group was significantly increased (P<0.001, Figure 15B). Besides, the number of cloned proliferating cells in the TG-101348 group was significantly lower than that in the control group (P<0.001). Compared with the TG-101348 group, the number of cloned proliferating cells in the TG-101348 + oe-CCL17 group was significantly increased (P<0.001, Figure 15C). The results of the wound healing assay showed that the migration rate of TG-101348 group was significantly lower than that of the control group (P<0.001). Compared with the TG-101348 group, the migration rate of the TG-101348 + oe-CCL17 group was significantly increased (P<0.001, Figure 15D). WB results showed that the phosphorylation of JAK2 and STAT3 in the TG-101348 group was significantly lower than that in the control group (P<0.001). Compared with the TG-101348 group, the phosphorylation of JAK2 (P=0.001) and STAT3 (P<0.001) in the TG-101348 + oe-CCL17 group was significantly increased (Figure 15E).

Figure 15 Molecular docking experimental verification (taking TG-101348 as an example). (A) WB analysis of CCL17 protein expression. (B) Transwell analysis of the invasion. (C) Clone formation assay of cell proliferation. (D) Wound healing assay of cell migration. Crystal violet staining was used in both the Transwell cell invasion assay and the plate cloning assay. The scales of the figure of the Transwell cell invasion assay and wound healing assay are 1:200 µm and 1:500 µm, respectively. WB, western blot.

Discussion

Herein, we aimed to identify DECRGs related to THCA prognosis and developed a risk model. The high-risk group was related to poor prognosis and worse immunotherapy response. Notably, CCL17 among the three hub genes stood out due to its predictive ability comparable to the risk model. Therefore, the clinical relevance of CCL17 in THCA and the underlying mechanism were also explored.

First, we constructed a risk model based on three hub genes (ACKR3, CCL2, CCL17) and found that a high-risk group exhibited an unfavorable survival. To gain insights into the underlying mechanisms of the poor prognosis in the high-risk group, GSVA was performed, and the results revealed that the high-risk group was closely associated with angiogenesis, EMT, and various immune-related signaling pathways. Additionally, our investigation unveiled significant decreases in the immune infiltration levels of CD8 T cells and M2 macrophages, while showing significant increases in TIDE score, Tfh cells, and M1 macrophages within the high-risk group. These findings imply that the poor prognosis observed in the high-risk group might be attributed to the activation of immune escape mechanisms and the suppression of immune responses.

After conducting the expression and survival analyses of the three hub genes, we compared the value of the hub genes and risk model in predicting THCA prognosis. K-M curves showed that only ACKR3 and CCL17 groups were linked to patient prognosis. Surprisingly, there were no significant differences in ACKR3, CCL17, and the risk model for predicting survival. ACKR3 is an atypical chemokine receptor and is a member of class A of G protein-coupled receptors (29). ACKR3 upregulation was related to the poor prognosis of various cancers, such as glioblastoma and colorectal cancer (30,31). Some experiments showed that overexpressed ACKR3 promoted the growth of lung cancer cells and cervical cancer (32,33). However, our results exhibited that ACKR3 was downregulated in tumor samples, which indicated that ACKR3 may not be involved in the growth of THCA. Besides, Li et al. identified that ACKR3 expression in primary tumors of patients with lymphatic metastasis is higher than that of patients without metastasis (34). Further research determined that ACKR3 promoted the EMT to enhance tumor migration by activating the TGF-β signal pathway (35). Herein, ACKR3 was not closely related to immune cell infiltrations, except for DC resting and ICG expression. It suggested that ACKR3 may promote THCA progression by activating the TGF-β signal pathway to promote the EMT instead of the immune response.

CCL17 is a member of the CC chemokine family, and it is a ligand for CCR4 (36,37). CCL17 was upregulated in some tumors and promoted the development and metastasis of some tumors (36,38). Consistent with the previous literature, we also found that the expression level of CCL17 was higher in THCA and its upregulation was related to shorter survival. The results of immune infiltration analysis showed that CCL17 was negatively correlated with CD8 T cells and M2 macrophages, positively correlated with follicular helper T cells and DCs, and positively correlated with a variety of human leukocyte antigen molecules. This may be related to the chemotaxis dominated by the CCL17-CCR4 axis. CCL17 efficiently chemotacticizes Tfh expressing CCR4 to the tumor site, promoting B cell activation and lymphoid follicle formation (39); at the same time, it recruits DCs, enhances their HLA molecules (HLA-DRB1/DMA/DQB1, etc.) expression, and improves antigen presentation ability (40). Moreover, CCL17-driven CCR4⁺ regulatory T cells are enriched in tumors, directly inhibiting CD8⁺ T cell function by secreting TGF-β/interleukin (IL)-10 and highly expressing CTLA4 (41), resulting in reduced CD8⁺ T cell infiltration. Treg and Th1 cells recruited by CCL17 secrete interferon-γ (IFN-γ), antagonizing IL-4/IL-13-driven M2 polarization (42). These results indicate that CCL17 is a core driver of immune escape in THCA and can be used as a therapeutic target.

To further explore the mechanism by which high expression of CCL17 promotes the poor prognosis of THCA, we found that CCL17 was enriched in the JAK-STAT signaling pathway through GSEA analysis. The experimental results showed that CCL17 promoted the phosphorylation of the JAK2 and STAT3 in this pathway, thereby activating this pathway and promoting the malignant phenotype development of tumor cells. This is mainly due to the fact that the STAT3 gene induces the expression of c-Myc and Cyclin D1, promotes the cell cycle process, and thus promotes the proliferation of tumor cells (43,44). In addition, STAT3 activation can antagonize the expression of anti-tumor T helper cell 1 cytokines (such as IL-12 and IFN-γ) mediated by NF-κB and STAT1, thereby inhibiting the anti-tumor immune response and contributing to the immune escape of tumor cells (45). STAT3 can also activate MMP-2/9 and VEGF genes to promote cell matrix degradation, participate in angiogenesis, and contribute to tumor cell invasion and metastasis (46). In many studies, inhibition of STAT3 can inhibit tumor development (47,48). In our study, we found that inhibiting the high expression of CCL17 by drugs can significantly inhibit the phosphorylation of STAT3, thereby inhibiting the development of the malignant phenotype of the tumor. These results also suggest that CCL17 can be used as a potential target for the treatment of THCA.

There are some limitations in our study. Firstly, the samples were only obtained from the TCGA database due to a lack of prognosis information in other datasets. In follow-up studies, we should supplement the additional sample based on our own cohorts to strengthen the scientific evidence. Secondly, although we have conducted cell experiments for expression validation, the potential function and mechanism, especially how CCL17 regulates the JAK-STAT signal pathway to activate immune escape to promote THCA deterioration, need to be verified by in vitro experiments. Our results paved the way for future mechanism exploration and offered a novel insight into THCA treatment.

Conclusions

In conclusion, we identified three hub CRGs (ACKR3, CCL2, CCL17) and constructed a risk model, which could steadily and effectively predict the prognosis of patients with THCA. The high-risk group exhibited a worse prognosis, which may be due to an impaired immune response. Moreover, CCL17 may serve as a core biomarker in THCA deterioration and metastasis by regulating the JAK-STAT signal pathway.

Acknowledgments

None.

Footnote

Reporting Checklist: The authors have completed the TRIPOD reporting checklist. Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2024-2654/rc

Data Sharing Statement: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2024-2654/dss

Peer Review File: Available at https://tcr.amegroups.com/article/view/10.21037/tcr-2024-2654/prf

Funding: This work was supported by the following grants: Xuzhou Medical University (Natural Science) Project (No. 2018KJ16).

Conflicts of Interest: All authors have completed the ICMJE uniform disclosure form (available at https://tcr.amegroups.com/article/view/10.21037/tcr-2024-2654/coif). The authors have no conflicts of interest to declare.

Ethical Statement: The authors are accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. This study was conducted in accordance with the Declaration of Helsinki and its subsequent amendments. The Ethics Committee of The Second Affiliated Hospital of Xuzhou Medical University deemed that this research is based on open-source data, so the need for ethics approval was waived.

Open Access Statement: This is an Open Access article distributed in accordance with the Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International License (CC BY-NC-ND 4.0), which permits the non-commercial replication and distribution of the article with the strict proviso that no changes or edits are made and the original work is properly cited (including links to both the formal publication through the relevant DOI and the license). See: https://creativecommons.org/licenses/by-nc-nd/4.0/.

References

1. Bai M, Ke S, Yu H, et al. Key molecules associated with thyroid carcinoma prognosis: A study based on transcriptome sequencing and GEO datasets. Front Immunol 2022;13:964891. [Crossref] [PubMed]
2. Xie F, Li L, Luo Y, et al. Long non-coding RNA LINC00488 facilitates thyroid cancer cell progression through miR-376a-3p/PON2. Biosci Rep 2021;41:BSR20201603. [Crossref] [PubMed]
3. Santhanam P, Solnes L, Nath T, et al. Real-time quantitation of thyroidal radioiodine uptake in thyroid disease with monitoring by a collar detection device. Sci Rep 2021;11:18479. [Crossref] [PubMed]
4. Ren X, Du H, Cheng W, et al. Construction of a ferroptosis-related eight gene signature for predicting the prognosis and immune infiltration of thyroid cancer. Front Endocrinol (Lausanne) 2022;13:997873. [Crossref] [PubMed]
5. Kee JY, Ito A, Hojo S, et al. CXCL16 suppresses liver metastasis of colorectal cancer by promoting TNF-α-induced apoptosis by tumor-associated macrophages. BMC Cancer 2014;14:949. [Crossref] [PubMed]
6. Gilchrist A. Chemokines and Bone. Handb Exp Pharmacol 2020;262:231-58. [Crossref] [PubMed]
7. Entschladen F, Drell TL 4th, Lang K, et al. Neurotransmitters and chemokines regulate tumor cell migration: potential for a new pharmacological approach to inhibit invasion and metastasis development. Curr Pharm Des 2005;11:403-11. [Crossref] [PubMed]
8. Liu Y, Gao Y, Lin T. Expression of interleukin-1 (IL-1), IL-6, and tumor necrosis factor-α (TNF-α) in non-small cell lung cancer and its relationship with the occurrence and prognosis of cancer pain. Ann Palliat Med 2021;10:12759-66. [Crossref] [PubMed]
9. Li YF, Han CC, Wang Y, et al. Combined PGE2 with TNF-α promotes laryngeal carcinoma progression by enhancing GRK2 and TRAF2 interaction. Neoplasma 2020;67:354-63. [Crossref] [PubMed]
10. Qiu X, Ye H, Li X, et al. IL-6/JAK2-dependent G6PD phosphorylation promotes nucleotide synthesis and supports tumor growth. Mol Metab 2023;78:101836. [Crossref] [PubMed]
11. Leelawat S, Leelawat K, Yimsoo T, et al. Antitumor Effects of Delta (9)-Tetrahydrocannabinol and Cannabinol on Cholangiocarcinoma Cells and Xenograft Mouse Models. Evid Based Complement Alternat Med 2022;2022:6477132. [Crossref] [PubMed]
12. Stephens OW, Zhang Q, Qu P, et al. An intermediate-risk multiple myeloma subgroup is defined by sIL-6r: levels synergistically increase with incidence of SNP rs2228145 and 1q21 amplification. Blood 2012;119:503-12. [Crossref] [PubMed]
13. Ciobanu DA, Poenariu IS, Crînguș LI, et al. JAK/STAT pathway in pathology of rheumatoid arthritis Exp Ther Med 2020;20:3498-503. (Review). [Crossref] [PubMed]
14. Moser B. Chemokine Receptor-Targeted Therapies: Special Case for CCR8. Cancers (Basel) 2022;14:511. [Crossref] [PubMed]
15. Wang C, Xiang R, Zhang X, et al. Doxycycline inhibits leukemic cell migration via inhibition of matrix metalloproteinases and phosphorylation of focal adhesion kinase. Mol Med Rep 2015;12:3374-80. [Crossref] [PubMed]
16. Zhang Y, Li Y, Thompson KN, et al. Polarized NHE1 and SWELL1 regulate migration direction, efficiency and metastasis. Nat Commun 2022;13:6128. [Crossref] [PubMed]
17. Bai H, Yang J, Meng S, et al. Oral Microbiota-Driven Cell Migration in Carcinogenesis and Metastasis. Front Cell Infect Microbiol 2022;12:864479. [Crossref] [PubMed]
18. Jeng KS, Jeng CJ, Jeng WJ, et al. Role of C-X-C chemokine ligand 12/C-X-C chemokine receptor 4 in the progression of hepatocellular carcinoma. Oncol Lett 2017;14:1905-10. [Crossref] [PubMed]
19. Gunes-Bayir A, Guler EM, Bilgin MG, et al. Anti-Inflammatory and Antioxidant Effects of Carvacrol on N-Methyl-N'-Nitro-N-Nitrosoguanidine (MNNG) Induced Gastric Carcinogenesis in Wistar Rats. Nutrients 2022;14:2848. [Crossref] [PubMed]
20. Chen J, Zhu H, Liu Q, et al. DEPTOR induces a partial epithelial-to-mesenchymal transition and metastasis via autocrine TGFβ1 signaling and is associated with poor prognosis in hepatocellular carcinoma. J Exp Clin Cancer Res 2019;38:273. [Crossref] [PubMed]
21. Neves KB, Rios FJ, Jones R, et al. Microparticles from vascular endothelial growth factor pathway inhibitor-treated cancer patients mediate endothelial cell injury. Cardiovasc Res 2019;115:978-88. [Crossref] [PubMed]
22. Fan W, Wang D, Li G, et al. A novel chemokine-based signature for prediction of prognosis and therapeutic response in glioma. CNS Neurosci Ther 2022;28:2090-103. [Crossref] [PubMed]
23. He C, He L, Lu Q, et al. The functions and prognostic values of chemokine and chemokine receptors in gastric cancer. Am J Cancer Res 2022;12:3034-50.
24. Zhang Z, Mao M, Wang F, et al. Comprehensive analysis and immune landscape of chemokines- and chemokine receptors-based signature in hepatocellular carcinoma. Front Immunol 2023;14:1164669. [Crossref] [PubMed]
25. Dubrot J, Du PP, Lane-Reticker SK, et al. In vivo CRISPR screens reveal the landscape of immune evasion pathways across cancer. Nat Immunol 2022;23:1495-506. [Crossref] [PubMed]
26. Xu J, Luo Z, Xu D, et al. Immunotyping of thyroid cancer for clinical outcomes and implications. Cancer Immunol Immunother 2025;74:221. [Crossref] [PubMed]
27. Gao B, Wang Y, Li C, et al. Estrogen-related genes influence immune cell infiltration and immunotherapy response in Hepatocellular Carcinoma. Front Immunol 2023;14:1114717. [Crossref] [PubMed]
28. Laha D, Nilubol N, Boufraqech M. New Therapies for Advanced Thyroid Cancer. Front Endocrinol (Lausanne) 2020;11:82. [Crossref] [PubMed]
29. Neves M, Fumagalli A, van den Bor J, et al. The Role of ACKR3 in Breast, Lung, and Brain Cancer. Mol Pharmacol 2019;96:819-25. [Crossref] [PubMed]
30. Salazar N, Carlson JC, Huang K, et al. A Chimeric Antibody against ACKR3/CXCR7 in Combination with TMZ Activates Immune Responses and Extends Survival in Mouse GBM Models. Mol Ther 2018;26:1354-65. [Crossref] [PubMed]
31. Yang J, Miao RR, Li YN, et al. Atypical chemokine receptor 3 induces colorectal tumorigenesis in mice by promoting β-arrestin-NOLC1-fibrillarin-dependent rRNA biogenesis. Acta Pharmacol Sin 2022;43:2967-76. [Crossref] [PubMed]
32. Liu H, Cheng Q, Xu DS, et al. Overexpression of CXCR7 accelerates tumor growth and metastasis of lung cancer cells. Respir Res 2020;21:287. [Crossref] [PubMed]
33. Xu L, Li C, Hua F, et al. The CXCL12/CXCR7 signalling axis promotes proliferation and metastasis in cervical cancer. Med Oncol 2021;38:58. [Crossref] [PubMed]
34. Li JT, Jia LT, Liu NN, et al. MiRNA-101 inhibits breast cancer growth and metastasis by targeting CX chemokine receptor 7. Oncotarget 2015;6:30818-30. [Crossref] [PubMed]
35. Kim N, Ryu H, Kim S, et al. CXCR7 promotes migration and invasion in head and neck squamous cell carcinoma by upregulating TGF-β1/Smad2/3 signaling. Sci Rep 2019;9:18100. [Crossref] [PubMed]
36. Rasheed K, Abdulsalam I, Fismen S, et al. CCL17/TARC and CCR4 expression in Merkel cell carcinoma. Oncotarget 2018;9:31432-47. [Crossref] [PubMed]
37. Yuan H, Nishikori M, Otsuka Y, et al. The EZH2 inhibitor tazemetostat upregulates the expression of CCL17/TARC in B-cell lymphoma and enhances T-cell recruitment. Cancer Sci 2021;112:4604-16. [Crossref] [PubMed]
38. Morimura S, Sugaya M, Oka T, et al. Increased Regulatory T Cells and Decreased Myeloid-Derived Suppressor Cells Induced by High CCL17 Levels May Account for Normal Incidence of Cancers among Patients with Atopic Dermatitis. Int J Mol Sci 2021;22:2025. [Crossref] [PubMed]
39. Maeda S, Murakami K, Inoue A, et al. CCR4 Blockade Depletes Regulatory T Cells and Prolongs Survival in a Canine Model of Bladder Cancer. Cancer Immunol Res 2019;7:1175-87. [Crossref] [PubMed]
40. Villani AC, Satija R, Reynolds G, et al. Single-cell RNA-seq reveals new types of human blood dendritic cells, monocytes, and progenitors. Science 2017;356:eaah4573. [Crossref] [PubMed]
41. Curiel TJ, Coukos G, Zou L, et al. Specific recruitment of regulatory T cells in ovarian carcinoma fosters immune privilege and predicts reduced survival. Nat Med 2004;10:942-9. [Crossref] [PubMed]
42. Barclay AN, Van den Berg TK. The interaction between signal regulatory protein alpha (SIRPα) and CD47: structure, function, and therapeutic target. Annu Rev Immunol 2014;32:25-50. [Crossref] [PubMed]
43. Yu H, Pardoll D, Jove R. STATs in cancer inflammation and immunity: a leading role for STAT3. Nat Rev Cancer 2009;9:798-809. [Crossref] [PubMed]
44. Bollrath J, Phesse TJ, von Burstin VA, et al. gp130-mediated Stat3 activation in enterocytes regulates cell survival and cell-cycle progression during colitis-associated tumorigenesis. Cancer Cell 2009;15:91-102. [Crossref] [PubMed]
45. Kujawski M, Kortylewski M, Lee H, et al. Stat3 mediates myeloid cell-dependent tumor angiogenesis in mice. J Clin Invest 2008;118:3367-77. [Crossref] [PubMed]
46. Han H, Yuan Y, Li C, et al. RNA-binding motif protein 28 enhances angiogenesis by improving STAT3 translation in hepatocellular carcinoma. Cancer Lett 2024;604:217191. [Crossref] [PubMed]
47. Deswal B, Bagchi U, Santra MK, et al. Inhibition of STAT3 by 2-Methoxyestradiol suppresses M2 polarization and protumoral functions of macrophages in breast cancer. BMC Cancer 2024;24:1129. [Crossref] [PubMed]
48. Kortylewski M, Kujawski M, Wang T, et al. Inhibiting Stat3 signaling in the hematopoietic system elicits multicomponent antitumor immunity. Nat Med 2005;11:1314-21. [Crossref] [PubMed]