Identification and validation of Birc5 as a novel activated cell cycle program biomarker associated with infiltration of immunosuppressive myeloid‐derived suppressor cells in hepatocellular carcinoma

Abstract Background Preclinical studies and clinical trials have demonstrated that tumor‐intrinsic activation of the cell cycle program impedes anticancer immunotherapy. Identification of cell cycle‐related biomarkers may provide novel therapeutic targets to augment the efficacy of immunotherapy in hepatocellular carcinoma (HCC). Method and Results Based on the genes related to cell cycle program, two clusters (Cluster 1 and Cluster 2) were detected in HCC patients via non‐negative matrix factorization algorithm. Multivariable‐adjusted Cox regression analysis indicated that the cell cycle gene‐based classification was a significant prognostic factor for predicting the clinical outcome of HCC patients. Cluster 1 showed shorter overall survival time and progression‐free interval time was associated with activated cell cycle program, higher infiltration of myeloid‐derived suppressor cells (MDSCs) and less sensitivity to immunotherapy. A three‐gene prognostic model, including BIRC5, C8G, and SPP1, was constructed to characterize the cell cycle‐based classification of HCC, which had strong robustness and a stable predictive performance. Notably, Birc5 was positively correlated with CD11b expression (a MDSC marker) in HCC tissue. Concordant high expression of Birc5 and intratumor infiltration level of MDSCs were correlated with worse prognosis of HCC patients. In vitro, hepatocellular Birc5 overexpression promoted immunosuppressive CD11b+CD33+HLA‐DR−MDSC expansion from human peripheral blood mononuclear cells. Genetically modified animal model of liver cancer revealed that Birc5 depletion upregulated the genes related to lymphocyte‐mediated immunity, natural killer cell‐mediated immunity, interferon‐gamma production, T‐cell activation, and T‐cell‐mediated cytotoxicity. These results suggest an immunosuppressive function of Birc5 in HCC. Conclusion Birc5 was a potential biomarker and inducer of intratumor infiltration of MDSCs, which led to T cell exclusion or dysfunction in tumor immune microenvironment, consequently resulting in reduced response to ICIs in HCC.


| INTRODUCTION
Hepatocellular carcinoma (HCC) is the most common form of liver cancer, and therapeutic interventions, including surgical resection, loco-regional percutaneous ablation, and chemotherapy, have limited effects on HCC with a high rate of recurrence within 5 years. [1][2][3] For patients with HCC who are diagnosed at advanced stages and not suitable for surgery, their median survival time is generally around 6 months. 4,5 Therefore, the overall survival (OS) of HCC patients is still far from satisfactory. During the last decade, immunotherapy has become a powerful weapon for cancer treatment and revolutionized the systemic management of advanced cancers (e.g., melanoma and renal cell carcinoma). 6,7 However, the response rate of immune checkpoint inhibitors (ICIs) targeting CTLA-4 or PD-1/PD-L1 in HCC is only 15%-20%, highlighting the urgent need to identify biomarkers that can robustly predict the effectiveness of ICI therapy. 8 Beyond that, clinical researches have uncovered that the immune cell composition of HCC is markedly associated with survival outcomes and response to immunotherapy. [8][9][10][11] Hence, it is proposed that biomarker-directed therapy with good patient stratification and prudent selection of therapeutic agent is important for enhancing the efficacy of immunotherapy against HCC. 6,8,10 Genes related to cell cycle program are frequently mutated in human neoplasms. Currently, increasing evidence has revealed that genes related to cell cycle program are capable of remodeling tumor microenvironment (TME) and can regulate tumor immune responses, underscoring the potential strategies of targeting cell cycle-associated genes for immunotherapy of human cancers. 12,13 For instance, inhibiting the expression of CDKN2A, a member of cyclin-dependent kinase (CDK) inhibitors, could lead to invasion, progression, and immunotherapy resistance in melanoma. [14][15][16] CDKN2A deletion repressed T-cell infiltration by suppressing chemokine expression in a cell cycledependent manner and might account for immunotherapy failure. 15 Besides, TP53, associated with STK11 and EGFR, has been proved to be a major determinant of tumor immune profile. Previous research showed that patients with TP53 mutation and without alteration in STK11 and EGFR had the highest levels of CD8 + T-cell and PD-L1 expression and exhibited best responses to PD-1 blockade. 17 In another study, TP53 mutation was related to activated immune cell infiltration and was regarded as a potential biomarker of immunotherapy in lung adenocarcinoma. 18 Moreover, retinoblastoma 1 (RB1), the downstream target of CDK4/6, was associated with immunological features of tumor tissues. [19][20][21][22][23] RB1-dependent cell cycle inhibition with CDK4/6 inhibitors increased cytotoxic T cell influx and limited intrinsic immunological exclusion in tumors. 19,20 However, the role of cell cycle-associated genes in remodeling immune profile and their impacts on immunotherapy efficiency in HCC have been rarely investigated. Hence, further exploiting the impact of cell cycle program on immunological features may provide novel approaches to enhance sensitivity to immunotherapy and develop more durable combinatorial therapies in HCC.
In this study, transcriptome analysis was performed to classify HCC into two clusters based on the genes related to cell cycle program via non-negative matrix factorization (NMF) algorithm. Next, the associations of immunosuppressive CD11b + CD33 + HLA-DR − MDSC expansion from human peripheral blood mononuclear cells. Genetically modified animal model of liver cancer revealed that Birc5 depletion upregulated the genes related to lymphocytemediated immunity, natural killer cell-mediated immunity, interferon-gamma production, T-cell activation, and T-cell-mediated cytotoxicity. These results suggest an immunosuppressive function of Birc5 in HCC.

Conclusion:
Birc5 was a potential biomarker and inducer of intratumor infiltration of MDSCs, which led to T cell exclusion or dysfunction in tumor immune microenvironment, consequently resulting in reduced response to ICIs in HCC.

K E Y W O R D S
Birc5, cell cycle, hepatocellular carcinoma (HCC), immunotherapy, myeloid-derived suppressor cells (MDSCs), tumor microenvironment the clusters with the prognosis and clinical features of HCC patients were evaluated. Integrated analysis was further conducted to explore the differences in immunogenomic profile and benefits of immunotherapy between the cell cycle gene-based subtypes. Then, a 3-gene prognostic model, including BIRC5, C8G, and SPP1, was constructed and validated in different datasets, which could better predict the prognosis of liver cancer. Finally, we identified and validated that Birc5 expression might serve as a prognostic biomarker associated with the infiltration of immune-suppressive myeloidderived suppressor cells (MDSCs) and a promising therapeutic target for HCC patients.

| Workflow and data collection
The overall workflow of this research is presented in Figure S1. Clinical information and RNA sequencing data of HCC were retrieved from Gene Expression Omnibus (GEO) repository (https://www.ncbi.nlm.nih.gov/geo/), International Cancer Genome Consortium (https:// dcc.icgc.org/) and The Cancer Genome Atlas (TCGA) data portal (https://portal.gdc.cancer.gov/), patients from which were named as GSE14520 cohort (n = 221), ICGC (LIRI-JP) cohort (n = 231) and TCGA-LIHC cohort (n = 342), respectively. Gene expression profiles of mouse liver cancers with or without BIRC5 depletion were downloaded from the GEO repository (GSE64804). 24 The TCGA-LIHC cohort was employed for data training and internal validation, while the ICGC (LIRI-JP) and GSE14520 cohorts were used for external validation. 693 reactome cell cycle-related genes were compiled from MsigDB (https://www.gsea-msigdb.org/gsea/msigd b/). 25 The somatic mutation data of patients from the TCGA-LIHC cohort were downloaded from the TCGA database. The frequencies of genes and mutation types were determined using the maftools package in R. 26

| Consensus molecular subtyping with NMF
NMF was performed to cluster TCGA-LIHC samples into cell cycle-related subtypes. 27 Consensus clustering was conducted using NMF package in R, with a setting of standard "brunet" and 50 iterations. The optimal number of subtypes was determined based on the cophenetic, dispersion, and silhouette coefficients, and the optimal clustering number was selected as 2. Then, the patients in TCGA-LIHC cohort were categorized into 2 molecular subtypes, termed Cluster 1 (166 patients, 48.5%) and Cluster 2 (176 patients, 51.5%).

| Development and validation of a prognostic cell cycle-related gene signature
1793 immune-associated genes were downloaded from the ImmPort database (https://www.immpo rt.org/share d/genel ists). 28 Veen diagram revealed that 82 genes were common between these immune-related genes and the 1202 DEGs. Finally, 24 immune-related genes with log 2 (count+1) value not less than 1 in all TCGA-LIHC samples were applied for univariable survival analysis using Cox proportional hazards model and were further integrated into the least absolute shrinkage and selection operator (LASSO) regression analysis. Multivariable survival analysis was conducted to investigate the overall effect of the immune-related genes on prognosis using the survival package in R software. The prognostic risk score (RS) model was constructed according to the multivariable co-efficiency multiplied by expression data. The RS model was further evaluated and validated with ICGC (LIRI-JP) and GSE14520. The median RS of all patients from each dataset was used as a cutoff value to categorize the whole group into low-and high-risk groups. To verify the role of prognostic features, the area under the curve (AUC) was constructed using the timeROC package in R software.

| Human specimens
Paired tumor and non-tumor tissues from 28 HCC patients were recruited from the First Affiliated Hospital of USTC between October 2020 and August 2021 for immunohistochemical analyses. Inclusion criteria for HCC patients were histopathological diagnosis of HCC. The clinical information of the enrolled HCC patients is presented in Table S1. Blood specimens were additionally collected from 6 healthy people without any history of infectious diseases within last year and without any concurrent diseases. The study protocol involving human subjects was approved by the Ethics Committees of the First Affiliated Hospital of USTC (2021-RE-113). All participation gave written informed consent in accordance with the Declaration of Helsinki.

| Immunohistochemistry analysis
Immunohistochemical staining was conducted according to a previous method. 29 Five-millimeter serial sections were prepared from formalin-fixed and paraffinembedded tissues, followed by deparaffinization and rehydration. Antigen retrieval was conducted with sodium citrate buffer (0.01 M, pH 6.0) for 10 min at 98°C. The endogenous peroxidase activity was suppressed by 3% hydrogen peroxide in dH 2 O for 10 min, and then incubated with 5% normal goat serum for 30 min. The serial sections were further incubated with polyclonal rabbit anti-Birc5 (1:100 dilution; 10508, Proteintech) or anti-CD11b antibody (1:100 dilution; 21851, Proteintech) at room temperature for 4 h. The sections were exposed to HRP-conjugated anti-rabbit secondary antibody for 1 h. After staining with 3,3′-diaminobenzidine, the sections were counterstained with hematoxylin. Five randomly selected visual fields (×200 magnification) were counted per each section and ImageJ software 30 was utilized to evaluate the area percentage of positive staining area of Birc5 and CD11b.

| Cell culture and transfection
The HCC cell line Huh7 and the immortal human liver cell line THLE-3 were obtained from the American Type Culture Collection and cultured in high-glucose DMEM (Invitrogen/Thermo-Fisher-Scientific) containing 10% FBS (Gibco). All cells were maintained at 37°C in a humidified incubator with 5% CO 2 . Lentiviral vectors carrying the BIRC5 gene and control lentiviral vector (Hunan Fenghui Biotechnology Co., Ltd) were added into the THLE-3 cells. The transfection efficiency was evaluated by western blotting. Fresh peripheral blood mononuclear cells (PBMCs) were cultured in DMEM medium containing 10% FBS (Gibco), 1 mM sodium pyruvate (Gibco) and 2 mM l-glutamine (Gibco).

| Flow cytometry analysis and evaluation of MDSC expansion
PBMCs were isolated from the buffy coats of anonymous human healthy donors using Ficoll-Paque Premium (GE Healthcare). PBMCs were treated with the supernatant collected from THLE-3 cells expressing Birc5 or control vector for 5 days. For flow cytometry analysis, cells were collected and washed and incubated with mouse serum (BioLegend) for 10 min at 4°C to block Fcγ receptors. Then, the cells were stained by anti-human myeloid markers CD11b, CD14, CD33 and HLA-DR, followed by flow cytometry analysis. Flow cytometry data were acquired using the CytoFLEX (BECKMAN COULTER), and the percentages of proliferating cells were analyzed by FlowJo software 10 (Tree Star, Inc.).

| Western blot analysis
Western blotting was conducted according to a previous method. 31 Total cell protein lysates were extracted on ice using the protease inhibitor cocktail-containing radioimmunoprecipitation assay buffer. Protein concentrations were detected with a Bio-Rad protein assay. Then, 10-30 μg of total cell protein extracts were subjected to 10% SDS-PAGE and electroblotted onto PVDF membranes (Millipore). The membranes were exposed to primary antibodies at 4°C overnight, followed by HRP-linked anti-mouse and anti-rabbit antibodies for 2 h. Primary antibodies used were rabbit anti-Birc5 polyclonal antibody (10508, Proteintech), mouse anti-PD-L1 monoclonal antibody (66248, Proteintech), mouse anti-GAPDH (60004, Proteintech) and alpha-Tubulin (66031, Proteintech).

| Statistical analysis
Continuous data were presented as mean ± SD. The differences between groups were compared using the two-tailed unpaired Student t test. The χ 2 test or Fisher's exact test was employed to assess the correlation between the cell cycle gene-based subtypes and clinical characteristics of HCC patients. The correlation between normally distributed variables was evaluated using Pearson's correlation test. p < 0.05 and |correlation coefficient (R)| > 0.3 indicate a significant correlation. All statistical tests were conducted with R (https://www.r-proje ct.org) and GraphPad Prism software (version 5.03; GraphPad Software Inc.). The Kaplan-Meier method was applied, and the survival of groups was compared using the log-rank test. Findings with two-sided p < 0.05 were deemed statistically significant.
Identification of differentially expressed genes (DEGs), DEGs enrichment analysis, Gene Set Variation Analysis (GSVA), Tumor Immune Dysfunction and Exclusion (TIDE) algorithm, and TIMER2.0 database analysis are detailed in Appendix S1.

| Identification of two clusters with different clinical outcomes using the NMF algorithm
The NMF algorithm was applied to identify the molecular subtypes of HCC patients based on the 578 cell cyclerelated genes. HCC patients in the TCGA-LIHC cohort were categorized into Cluster 1 and Cluster 2 molecular subtypes, according to the dispersion, cophenetic, and silhouette indicators ( Figure 1A). As displayed in Figure 1B, Kaplan-Meier survival analysis indicated that Cluster 1 was markedly related to worse OS and shorter progressionfree interval (PFI). The median OS time of HCC patients in Cluster 1 was remarkably shorter compared with Cluster 2 (46 vs. 83 months, p < 0.001), and the median PFI time of Cluster 1 was markedly shorter than that of Cluster 2 (13 vs. 33 months, p < 0.001), as shown in Table S2. Next, the clinical characteristics of the two clusters were compared. Significant differences were observed for tumor histological grade and pTNM stages, as shown in Table S3. The prognostic significance of the two clusters in the context of various clinicopathological variables was evaluated by multivariable-adjusted Cox regression analyses. The results revealed that HCC patients in Cluster 2 had significant better OS (HR = 0.51, p = 0.002) and PFI (HR = 0.62, p = 0.004) when compared to Cluster 1 ( Table 1).

| Cluster 1 exhibits hyperactivated cell cycle program and frequent TP53 mutation
To explore the underlying mechanisms related to the cell cycle gene-based classification of HCC, we first identified 1202 differentially expressed genes (DEGs), including  Figure S2). Then, GSVA analysis was performed and differentially enriched molecular pathways were identified, including 14 and 56 pathways positively correlated with Cluster 1 and Cluster 2, respectively ( Figure 1C; Table S4). Briefly, Cluster 1 was shown to be correlated with p53, DNA replication, and cell cycle pathways, whereas Cluster 2 was mainly correlated with molecular pathways related to metabolism, such as primary bile acid biosynthesis, fatty acid metabolism, and citrate cycle. Somatic mutation analysis revealed that each subtype possessed specific mutated genes. As shown in Figure 1D, TP53 (48%) and CTNNB1 (34%) were the most frequently mutated genes in Cluster 1 and Cluster 2, respectively. Taken together, these findings uncovered that activated cell cycle program and frequently TP53 mutation were observed in Cluster 1 and might be associated with unfavorable clinical outcomes of HCC.

| Construction and validation of immune gene signature for prognostic prediction of HCC
Among the 1202 DEGs, 82 genes were immune-related genes extracted from the Immport database (Figure 2A). A total of 24 immune genes were selected for the construction of RS model in all enrolled TCGA-LIHC samples. A prognostic RS model composing of BIRC5, SPP1, and C8G was constructed via the LASSO algorithm. The results showed that the higher expression of BIRC5 or SPP1 implied a worsen prognosis of HCC ( Figure 2B). Moreover, the number of deaths in the TCGA-LIHC cohort was increased with increasing RS values ( Figure 2C,D). To confirm and validate the performance of the prognostic RS model, Kaplan-Meier analysis was conducted to evaluate the association between OS and the 3-gene-signature risk model in TCGA-LIHC and the validation datasets including GSE14520 and ICGC (LIRI-JP). It was found that patients with a higher RS had shorter OS time ( Figure 2E). Furthermore, the reliability of the prognostic RS model was evaluated by time-dependent ROC curves. As shown in Figure 2F, the AUCs were 0.785, 0.700 and 0.671 for 1-, 3-, and 5-year survival, respectively, in the TCGA-LIHC training set. In the ICGC (LIRI-JP) validation cohort, the AUCs were 0.693, 0.767, and 0.810 for 1-, 3-, and 4-year survival, respectively. In the GSE14520 cohort, the AUCs were 0.648, 0.664, and 0.595 for 1-, 3-, and 5-year survival, respectively. These results demonstrated a great potential of the constructed RS in estimating the overall survival of HCC patients. T A B L E 1 Associations of molecular subtypes and clinicopathologic variables with overall survival (OS) and progression-free survival (PFS), assessed using multivariable-adjusted Cox regression.

MDSCs and is not sensitive to ICI therapies
Genomic alterations have been proposed to play potential roles in modulating tumor immunity and immune infiltration patterns. 32 Then, the differences in immune profiles between the two clusters were further investigated. The M2 subtype of tumor-associated macrophages, cancer-associated fibroblasts (CAFs) and MDSCs are three cell types that suppress T cell infiltration in tumors. 33 By utilizing the Tumor Immune Dysfunction and Exclusion (TIDE) algorithm, 34 MDSC, CAF and M2 intratumor infiltration scores were calculated in patients from TCGA-LIHC cohort. Multivariable-adjusted Cox regression analyses revealed that among the three immune suppressive cells, the MDSC infiltration score was markedly related to OS and PFI like pTNM stages (Table 2), suggesting intratumor MDSC infiltration level was a significant prognostic factor for predicting the clinical outcomes of HCC patients. Notably, our results also revealed that patients in Cluster 1 had a significantly higher TIDE prediction score accompanied by a higher MDSCs infiltration score when compared to those in Cluster 2 ( Figure 3A). Additionally, MDSC infiltration score was positively correlated with T-cell exclusion score in C1 ( Figure 3B) and survival analysis revealed that patients in Cluster 1 with a higher MDSCs score had shorter OS time compared to those with a lower MDSCs score (p = 0.0146) ( Figure 3C). Collectively, these results revealed that intratumor infiltration of MDSCs level was an independent prognostic factor for HCC and patients in Cluster 1 possessed higher TIDE and MDSCs scores, which might account for the insensitivity of C1 to ICI therapies.

T A B L E 2 Associations of TIDE scores and clinicopathologic variables with overall survival (OS) and progression-free survival (PFS)
, assessed using multivariable-adjusted Cox regression.

Birc5 and higher MDSC infiltration is correlated with poor prognosis of HCC
Among the 3 prognostic genes, BIRC5 is the only gene belonging to cell cycle-related and immune-related genes ( Figure 2A). Immunohistochemical staining confirmed that Birc5 was overexpressed in the nucleus of HCC tumor cells ( Figure 4A). Survival analysis revealed that patients with higher Birc5 expression had shorter OS time and DFI in the TCGA-LIHC cohort ( Figure 4B). We further created a nomogram combining Birc5 and clinicopathological information using the TCGA-LIHC cohort (C index, 0.71), and the calibration curves illustrated that the predicted and actual survival had a good consistency ( Figure S3). In addition, Birc5 expression was positively correlated with MDSC infiltration score in TCGA-LIHC cohort ( Figure 4C). In HCC tumor tissue samples, CD11b was used to indicate MDSC infiltration level, and a positive correlation between CD11b and Birc5 staining was observed ( Figure 4D,E). Survival analysis indicated that the patients with Birc5 overexpression and MDSC infiltration score had shorter OS time when compared to those with HCC. By utilizing TIMER dataset, it was observed that BIRC5 expression was positively correlated with MDSC infiltration score in multiple cancers ( Figure 5A,B), and survival analysis demonstrated that patients with higher BIRC5 expression and MDSC infiltration score had shorter OS time in sarcoma, lung adenocarcinoma, kidney renal clear cell carcinoma, and adrenocortical carcinoma ( Figure 5C). Therefore, these findings implied that BIRC5 might be used as a biomarker associated with MDSC infiltration in multiple cancers, including HCC, and concordant overexpression of Birc5 and higher MDSC infiltration was correlated with poor prognosis of several tumors.

| Birc5 regulates the immune microenvironment of HCC
Reportedly, Birc5 inhibits apoptosis and promotes proliferation, invasion, migration, and drug resistance in HCC cells. [35][36][37] To further investigate whether Birc5 contributes to MDSC infiltration, the expression of Birc5 was evaluated in the normal human hepatocyte cells THLE-3 and human liver cancer cells Huh7. It was found that Birc5 was highly expressed in Huh7 when compared to THLE-3 ( Figure S4). Subsequently, THLE-3 cells with stable Birc5 overexpression were established using lentivirus, and the efficacy of overexpression was detected by western blotting ( Figure 6A). Next, PBMCs from healthy human donors were cultured in the conditioned medium collected from the supernatant of THLE-3 cells with stable Birc5 overexpression. The flowchart of cultivation strategy is shown in Figure 6B. Flow cytometry showed that the supernatant of Birc5-overexpressing THLE-3 cells significantly promoted the amplification of MDSCs ( Figure 6C). By analyzing the whole-genome expression profile of mouse liver cancer  with or without BIRC5 depletion (GSE64804), 343 DEGs were obtained ( Figure 6D). KEGG pathway enrichment analysis revealed that DEGs related to BIRC5 depletion were enriched in chemokine, cytokine-cytokine receptor interaction, T-cell receptor, JAK-STAT, and NF-kappa B pathways ( Figure 6E). Heatmap showed that 30 immune-related genes enriched in these pathways were dysregulated upon BIRC5 depletion ( Figure 6F). Gene ontology enrichment analysis indicated that these 30 genes were remarkably enriched in adaptive immune response, chemokine-mediated signaling pathway, lymphocyte-mediated immunity, natural killer cell-mediated immunity, interferon-gamma production, T-cell activation, and T-cell-mediated cytotoxicity ( Figure 6G). Taken together, these findings demonstrate that Birc5 plays a vital role in modulating immune profile and response in HCC. The underlying mechanism might be associated with the function of Birc5 in recruiting MDSCs and regulating the expression of cytokine or chemokine, as well as inhibiting the activation of antitumor immune cells.

| DISCUSSION
Growing evidence from preclinical studies and clinical trials indicates that the hyperactivation of cell cycle program in tumor cells might confer protection from immune surveillance. 38 Therefore, identification and validation of new biomarkers related to cell cycle program may offer novel mechanistic insights and therapeutic targets to improve the immune responses of HCC to immunotherapy. In this work, based on the expression profiles of genes related to cell cycle program, we stratified HCC from the TCGA-LIHC cohort into two clusters (Cluster 1 and Cluster 2) using the NMF method. Multivariable-adjusted Cox regression analysis showed that the cell cycle genebased classification like pTNM stage was a significant prognostic factor for predicting the clinical outcome of HCC patients. Briefly, HCC patients in Cluster 1 presented worse OS and PFI compared with those in Cluster 2. GSVA analysis revealed that genes related to p53, DNA replication, and cell cycle pathways were significantly enriched in the Cluster 1. Genomic variation analysis showed that TP53 mutation was more frequently observed in Cluster 1, whereas Cluster 2 was characterized by CTNNB1 mutation. Previously, TP53 mutation has been reported to be a potential prognostic marker of several cancers, including HCC, lung cancer, clear cell renal cell carcinoma, and acute myeloid leukemia. 39 Mechanistically, TP53 mutation has been reported to play vital roles in determining tumor immune profile and PD-L1 expression in lung adenocarcinoma. 17 Mutant TP53 suppressed innate immune signaling and promoted immune evasion. 40 Therefore, these results imply that activated cell cycle program may account for unfavorable clinical outcomes and affect the tumor immune profile of HCC patients.
Next, we further investigated the differences in immune cell infiltration patterns between Cluster 1 and Cluster 2. TIDE algorithm was used to predict the therapeutic responses to ICIs, and the results demonstrated that HCC patients in Cluster 1 exhibited a higher TIDE score compared with those in Cluster 2, implying a worse response of patients in Cluster 1 to ICI therapy. Among three immunosuppressive cells (MDSCs, CAF, and M2), MDSCs was predominant in Cluster 1 and multivariableadjusted Cox regression analysis revealed that the MDSC infiltration score was a significant prognostic factor associated with OS and PFI. It has been well known that MDSCs are a pivotal cell population that contributes to the formation of tumor-induced immunosuppression through restricting T-cell accumulation in the vicinity of tumor cells and promoting tumor cell evasion from the immune system. 33,41 Moreover, the positive correlation between MDSC infiltration and T-cell exclusion has been proven to be a potential mechanism corresponded to the failure of antitumor therapy including ICI in several tumors. 20,33,34 Zhou et al 42 showed the CD11b + CD33 + HLA-DR − MDSCs was upregulated in both peripheral blood and intratumor tissues and could inhibit CD8 + T-cell proliferation in patients with HCC. Consistent with their findings, we found that MDSCs were positively correlated with T-cell exclusion score and survival analysis revealed that HCC patients with a higher MDSC score had shorter OS time compared to those with a lower MDSC score in Cluster 1. Taken together, these findings imply that a higher level of MDSC infiltration may account for poor prognosis and immune therapy failure in HCC patients in Cluster 1, probably through inhibiting the intratumor accumulation of cytotoxic T cells.
To characterize the cell cycle-based classification of HCC, we constructed a three-gene prognostic model, including BIRC5, C8G, and SPP1, which had strong robustness and a stable predictive performance in the datasets from two different cohorts. Notably, among the three risk genes, BIRC5 was the only gene positively correlated with the infiltration of MDSCs in HCC. Birc5, also known as Survivin, is significantly upregulated in several tumors. 43,44 It has been shown that Birc5 is a poor prognostic marker of HCC by promoting cancer cell proliferation, inducing angiogenesis, and attenuating the sensitivity of tumor cells to radiotherapy and chemotherapy. 45 Recently, the relationships among Birc5, immune cell infiltration, and the clinical outcome of multiple cancers have been reported. 43 It was found that Birc5 expression was related to immune cell infiltration in the TME of several tumors. However, the underlying mechanism is less reported in HCC. Our results showed that Birc5 expression was positively correlated with MDSC infiltration level and CD11b expression (a MDSC marker) in HCC as well as in multiple cancers. Of note, patients with upregulated expression of Birc5 and higher level of MDSC infiltration exhibited worse clinical outcomes in HCC and other cancers such as KIRC, LUAD, and SARC. Taken together, we speculate that tumor-intrinsic expression of Birc5 may contribute to the infiltration of MDSCs in tumor immune microenvironment.
To further investigate the underlying mechanism, Birc5 was overexpressed by lentivirus in human hepatocytes in vitro. By treating human PBMC with the supernatant of hepatocytes with Birc5 overexpression, our findings demonstrated that hepatocellular Birc5 overexpression promoted immunosuppressive CD11b + CD33 + HLA-DR − MDSC expansion. RNA-seq revealed that Birc5 overexpression could affect the expression levels of several inflammatory factors or chemokines, such as IL1B, TGFBR2, MME, CCL2, and CXCR4, in hepatocellular cells, and the genes related to Birc5 overexpression were enriched in pathways related to hematopoietic cell lineage, TGF-β, cytokine-cytokine receptor interaction, p53, and Wnt signaling (data not shown). Among these genes, IL1B, which encodes IL-1β, has been shown to promote the infiltration of polymorphonuclear MDSCs in renal cell carcinoma. 46 Furthermore, we analyzed the microarray data of an animal liver cancer model with or without BIRC5 depletion. 24 The findings revealed that the DEGs upon BIRC5 depletion were enriched in chemokine-mediated signaling pathway, lymphocyte-mediated immunity, natural killer cell-mediated immunity, interferon-gamma production, T-cell activation, and T-cell-mediated cytotoxicity. Moreover, Li et al 24 showed that Birc5 deficiency could lead to increased CD3 + T lymphocytes in tumor sites and remarkable HCC suppression. Altogether, these findings suggest that Birc5 may be a regulator of antitumor immune responses by inhibiting the activation of natural killer cells and T cells.
Birc5 expression was related to immune cell infiltration in the tumor microenvironment of several tumors. However, the underlying mechanism is less reported. In this study, through bioinformatics analyses and in vitro experiments, we showed that activated cell cycle program was associated with immunosuppression, weakened immunotherapy efficiency, and worse prognosis in HCC patients, and that Birc5 was a potential biomarker and inducer of intratumor infiltration of MDSCs. MDSCs are a heterogeneous group of pathologically activated immature cells and can be further identified into two different subsets: polymorphonuclear-(PMN-) and monocytic-(M-) MDSCs. 47 In HCC, the roles of M-MDSC and PMN-MDSCs in suppressing antitumor immunity remain controversial. For instance, increased accumulation of M-MDSCs in mice liver tissues with liver fibrosis was significantly correlated with decreased infiltration of lymphocytes and increased tumorigenicity, instead of PMN-MDSCs. 48 Nevertheless, Zhou et al 42 revealed that PMN-MDSCs were able to inhibit tumor infiltration of IFN-γ+TNF-α+CD8 + T lymphocytes through the expression of Arg-1. Therefore, to clearly illustrate the immune-relevant function of Birc5 in regulating HCC tumor immune microenvironment, we still need to further answer which type of MDSCs contribute to the immunosuppressive tumor microenvironment and by which signaling pathways, cytokines, and chemokines, Birc5 impacts the amplification and infiltration of MDSCs in tumor tissue, as well as the clinical significance of targeting Birc5 in immunotherapy of HCC. Cancer therapeutics targeting Birc5/survivin include survivin-partner protein interaction inhibitors, survivin homodimerization inhibitors, survivin gene transcription inhibitors, survivin mRNA inhibitors, and survivin immunotherapy. 44 YM155 is a small molecule strongly inhibiting Birc5 expression at both the protein and mRNA levels, and we plan to use YM155 for Birc5 inhibitory experiments in our future research.
There are some other HCC-related markers of note as revealed by some important studies. Wang et al 9 identified a novel eight-gene (PAM, NUP155, GOT2, KDELR3, PKM, NSDHL, ENO1, and SRD5A3) signature associated with cell glycolysis to predict survival of HCC patients. Feng et al 11 identified ARID1A as a prognostic biomarker associated with immune infiltrates in HCC: high ARID1A expression was associated with worse prognosis and low immune cell infiltration; decreasing ARID1A expression inhibited cell cycle, proliferation, migration, and invasion in HCC.
The main limitation of this study is the lack of further in vivo experiments to validate the effect of biomarkerdirected therapy. Further explorations should be focused on the underlying mechanisms by which Bric5 regulates the infiltration of different types of MDSCs (PMN-MDSCs or M-MDSCs) in the TME and on how they contribute to T-cell exclusion or dysfunction in HCC.

| CONCLUSIONS
Birc5 was a potential biomarker and inducer of intratumor infiltration of MDSCs, which led to T-cell exclusion or dysfunction in tumor immune microenvironment, consequently resulting in reduced response to ICIs in HCC.