{'fraction_non_alphanumeric': 0.03860438353298112, 'fraction_numerical': 0.021237641889417794, 'mean_word_length': 4.978111319574734, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 0, 'https://': 3, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A CD8+ T cell-associated immune gene panel for prediction of the prognosis and immunotherapeutic effect of melanoma

Background Skin cutaneous melanoma (SKCM) is the most frequently encountered tumor of the skin. Immunotherapy has opened a new horizon in melanoma treatment. We aimed to construct a CD8+ T cell-associated immune gene prognostic model (CDIGPM) for SKCM and unravel the immunologic features and the benefits of immunotherapy in CDIGPM-defined SKCM groups. Method Single-cell SKCM transcriptomes were utilized in conjunction with immune genes for the screening of CD8+ T cell-associated immune genes (CDIGs) for succeeding assessment. Thereafter, through protein-protein interaction (PPI) networks analysis, univariate COX analysis, and multivariate Cox analysis, six genes (MX1, RSAD2, IRF2, GBP2, IFITM1, and OAS2) were identified to construct a CDIGPM. We detected cell proliferation of SKCM cells transfected with IRF2 siRNA. Then, we analyzed the immunologic features and the benefits of immunotherapy in CDIGPM-defined groups. Results The overall survival (OS) was much better in low-CDIGPM group versus high CDIGPM group in TCGA dataset and GSE65904 dataset. On the whole, the results unfolded that a low CDIGPM showed relevance to immune response-correlated pathways, high expressions of CTLA4 and PD-L1, a high infiltration rate of CD8+ T cells, and more benefits from immunotherapy. Conclusion CDIGPM is an good model to predict the prognosis, the potential immune escape from immunotherapy for SKCM, and define immunologic and molecular features.


Introduction
Pigment cells originate from the neural crest locate in the epidermis, in which they primarily serve as protectors for keratinocytes against UV-elicited DNA injury (1,2). Skin cutaneous melanoma (SKCM), developing from vicious transformation of melanocytes, turns out to be the fatal form of skin cancer (3). The proportion of new SKCM cases is mounting in spite of the decreasing proportions of most cancers. Moreover, SKCM causes about 72% of deaths in skin cancer on account of its strong potentials in metabolism and metastasis (4,5).
Recently, immunotherapy has revolutionized the therapy for cancer. In particular, on the condition of the approval of the use of CTLA-4, PD-L1, and PD-1-specific immune checkpoint inhibitors, immunotherapy performs better than conventional therapies in lengthening the overall survival (OS) of cases of heterogeneous tumors (6)(7)(8)(9). SKCM is a class of tumor displaying the most sensitive to immunotherapeutic methods (10,11). According to the latest clinical study, after SKCM patients underwent nivolumab + ipilimumab combination therapy, the 5-year OS rate was 52% (12). In addition, it appears that T cells are promotors for immune responses and immunotherapy, among which CD8 + T cells occupy a dominant position (13,14). It is worth noting that CD8+ T cells play an important role in the prognosis of melanoma. It has been reported that the oxidative phosphorylation CD8+ T cell subset is predictive of immunotherapy resistance in melanoma patients (15), and CD8+ T-cell infiltration could influence patient survival in cutaneous melanoma directly (16). Therefore, tumor-infiltrating CD8 + T cells-associated immune genes (CDIGs) are probably targets for the identification of SKCM patients with sensitivity to immunotherapy.
Investigating prognostic markers for SKCM is central to this study, which can be conducive to the prediction of traditional therapeutic outcomes and the suggestion of immunotherapeutic value. For a detailed assessment of CD8+ T cell-related genes (the differential genes in CD8+ T cells) in SKCM, we explored the single-cell RNA sequencing (scRNA-seq) dataset using the Tumor Immune Single-Cell Hub (TISCH). Here, we constructed the CDIGs prognostic model (CDIGPM). Furthermore, we described the immunological characteristics of CDIGPM defined groups. Finally, we detected the ability of CDIGPM to predict the prognosis and immunotherapeutic efficacy in SKCM patients. Our research show that CDIGPM is an encouraging prognostic model.

Recognition of CDIGs in SKCM
The 708 differentially expressed genes (DEGs) in CD8 + T cells in SKCM (Table S1) were obtained from TISCH (GSE120575, http://tisch.comp-genomics.org/). As a scRNAseq database concentrating on the tumor microenvironment (TME), TISCH enables the exploration of TME. The filtering of DEGs in CD8 + T cells was based on the threshold of P < 0.05 and | log2 FC| ≥ 0.5. Meanwhile, the updated immune genes were retrieved from ImmPort and InnateDB. Later, DEGs of CD8 + T cells and immune genes were intersected to obtain CDIGs (Table S2), which would be examined later in this study.

Recognition of hub CDIGs
Through the online database STRING (https://string-db.org/), we generated a PPI network of CDIGs. Then, 32 hub CDIGs were filtered by the number of adjacent nodes ≥ 30.

Construction of the CDIGPM
From 32 hub CDIGs, using univariate and multivariate Cox regression analyses, six genes (MX1, RSAD2, IRF2, GBP2, IFITM1, and OAS2) were screened to construct a CDIGPM (Table S3). In the model, we calculated the CDIGPM score using the formula: CDIGPM = [Expression value × gene coefficient]. After that, we drawn Kaplan-Meier survival curves and conducted a log-rank test to explore the performance of the CDIGPM on the TCGA and GEO cohorts. Additionally, by application of univariate and multivariate Cox regression analyses, we clarified independent prognostic value of the CDIGPM.

Cell proliferation detection
The proliferation of SKCM cells was detected by Cell Counting Kit-8 kit (Beyotime Biotechnology, Shanghai, China). Approximately 1×10 3 cells were incubated in triplicate in 96-well plates. At 48h, the Cell Counting Kit -8 reagent (10mL) was added to each well and incubated at 37°C for 2h. Absorbance at 450 nm was used.
Immunologic features and immunotherapeutic effect in the two CDIGPM groups Limma package of R was used for differential expression analysis of all genes. ClusterProfiler package of R was used for gene set enrichment analysis (GSEA) on GO and KEGG gene sets. Genetic alteration data were downloaded from TCGA for gene mutation analysis. Correlation analyses were performed between CDIGPM and CTLA4 and CD274 expressions.
To determine immune features of SKCM samples in different CDIGPM subgroups, we used CIBERSORT (https:// cibersort.stanford.edu/) to calculate the relative proportion of 22 types of immune cells. The relative proportions of 22 types of immune cells were then compared between the two CDIGPM subgroups, and the results were showed in a landscape map.

Statistical analysis
The continuous variables were analyzed by t-test. The categorical data were analyzed by c 2 test. Kaplan-Meier survival analysis and the log-rank test were used for univariate survival analysis. Cox regression model was used for multivariate survival analysis. P value < 0.05 was considered to be significant differences.

CDIGPM
As unveiled by the univariate Cox regression analysis, totaling 29 from 32 hub CDIGs (P < 0.05) were clearly connected with OS in TCGA cohort ( Figure 2A; Figure S1). Then, the multivariate Cox regression analysis disclosed that six genes (MX1, RSAD2, IRF2, GBP2, IFITM1, and OAS2) were prognostic hallmarks, which were employed to build a CDIGPM (Table S3). In Cox model, the CDIGPM score of all samples were calculated using the formula: CDIGPM = [Expression value × gene coefficient] (Table  S3). We then explored the expression of these genes in SKCM using GEPIA (http://gepia.cancer-pku.cn/). Based on TCGA and GTEx data, the expression levels of MX1, RSAD2, IRF2, GBP2, IFITM1, and OAS2 were showed in the Figure S2. The performance of the CDIGPM in prognosis prediction was verified by the Kaplan-Meier survival curve and log-rank test on the TCGA and GEO cohorts. With the cutoff value of the median CDIGPM, samples fell into low-CDIGPM group and high CDIGPM group. It was disclosed that low-CDIGPM patients had better OS vs. high-CDIGPM patients ( Figures 2B,  C). We then evaluated the independent prognostic value of CDIGPM via univariate and multivariate Cox regression analyses. Results showed that score of CDIGPM, tumor stage and age (P < 0.05) were independent prognostic factors (Figures 2D, E; Table S4). Since the absolute value of IRF2 coefficient is the largest (-0.419071422539304), we wonder whether it affects the prognosis by regulating the malignant behavior of tumor cells. Results showed that down-regulation of IRF2 promoted SKCM cell proliferation ( Figure 2F). keratinocyte differentiation in high-CDIGPM group ( Figure 3A). Further, low-CDIGPM was connected immune responses related pathways ( Figure 3B). Enriched KEGG gene sets in the two CDIGPM groups were showed in Figure S3. Next, gene mutation analysis was performed and the gene mutation in low-CDIGPM score group and high-CDIGPM score group was showed in the Figures 3C, D. Subsequently, the relation between CDIGPM score and immune checkpoint genes were explored. It followed that the CDIGPM score displayed negative relevance to CTLA4 expression ( Figures 4A, B), CD274 expression ( Figures 4C, D), and PDCD1 (PD1) expression ( Figures 4E, F).

Immunologic features of the two CDIGPM groups
We detected the constituents of immune cells in the two CDIGPM groups. The results uncovered more activated memory CD4 + T cells, M1 macrophages (anti-tumor phenotype), and CD8 + T cells in low-CDIGPM group, but more M2 macrophages (pro-tumor phenotype) in high-CDIGPM group (Figures 5A, B). Then, we defined the immune and molecular function between the two groups by certain gene signatures. As a result, the immune and molecular function were more active in low-CDIGPM group ( Figure 5C).  Figure 6A shows the clinical features in high-and low-CDIGPM groups. We could find from Figures 6B, C that CDIGPM was related to tumor stage and size (P = 0.001, c 2 test). Specifically, there were more Stage I samples and fewer Stage II samples in the low-CDIGPM group versus the high-CDIGPM group (P = 0.001). In Figure 6C, more T0-T2 samples were belonged to the low-CDIGPM group and more T3-T4 samples were belonged to the high-CDIGPM group (P = 0.001).

Relationship between CDIGPM and immunotherapy
To explore the role of CDIGPM in immunotherapeutic effect, we analyzed the expression profile in SKCM patients (GSE35640). The analysis yielded that the CDIGPM in SKCM patients who responded to immunotherapy was lower than it in SKCM patients who did not respond to immunotherapy ( Figure 7A). Furthermore, receiver operating characteristic (ROC) analysis was performed to determine the performance of CDIGPM score in immunotherapeutic efficacy prediction ( Figure 7B). Our findings mean that CDIGPM could predict whether there is a response to immunotherapy in SKCM patients. The graphical abstract of this research is displayed in the Figure S4.

Discussion
Immunotherapy is the most promising therapy for several tumors, including melanoma. Its mechanism is to activate autologous immune responses by interfering with the tolerance of human cancer and re-inducing the tumor-resistant impacts on the TME (8,(17)(18)(19). Nevertheless, some patients couldn't get satisfactory efficacies due to the complex mechanisms underlying tumor immunity (20). Melanoma treatment has recently made headway with the advent of immunotherapies (anti-CTLA4 and anti-PD-1 antibodies). Whereas sustained responses may be observed with anti-PD-1 antibody therapy, around 60% of patients still develop resistance (12,21). Cancer cell phenotype plasticity (22, 23), tumor microenvironment (24), the expression of immune checkpoint genes (25) and other factors may be related to immune escape in melanoma. What's more, the reported genomic and immune biomarkers are not accurate enough in evaluating therapeutic effects (26). Hence, though it is challenging, finding a better predictor is in sore need for the accurate assessment of clinical outcomes before immunotherapy.
T cells are tumor-resistant effector cells of paramount significance owing to their direct attacks on cancer cells. It has been recently presented that the outcome of immune checkpoint therapy (ICT) targeting T cells is promising in melanoma cases. The efficacy of ICT is only favorable in certain tumor cases, which appears to be affected by the extent of the activation or infiltration of immune cells, especially CD8 + T lymphocytes (16,27,28). From GSE120575, which was deposited in the TISCH, 708 DEGs of CD8 + T cell from SKCM were obtained. Then, the overlapped genes of immune genes and differential genes in CD8 + T cells were regarded as CDIGs. Through PPI, 32 hub CDIGs were screened for subsequent analysis. We screened 29 CDIGs pertaining to the prognosis of SKCM from these hub genes by univariate Cox analysis. By application of multivariate Cox regression analysis, we built the CDIGPM. Using TCGA and GEO arrays, the CDIGPM was proven to be an effective model for the prognosis of SKCM.
CDIGPM is composed of six genomes: MX1, RSAD2, IRF2, GBP2, IFITM1, and OAS2. MX1 is an interferon-inducible dynamin GTPase that is essential for the suppression of replication of multifold viruses (29). It impedes the early stage in the replication cycle of discrepant viruses to exert inhibitory effect on these viruses (30). RSAD2, a gene stimulated by Relationship between CDIGPM and immunotherapy. (A) CDIGPM in patients who did not respond to immunotherapy is higher than it in patients who responded to immunotherapy. (B) ROC curve. (*p < 0.05). CDIGPM: CD8 + T cells-associated immune genes prognostic model.
interferons, is engaged in congenital immunity and primarily accountable to antiviral responses. It is reported that knockdown of RSAD2 makes mature dendritic cells unable to stimulate the production of proinflammatory cytokines and T cell proliferation (31). Besides the participation in antiviral immune responses, RSAD2 is also a potent driver for adaptive immune responses mediated by mature dendritic cells (31). IRF2, a constitutive transcription factor pertaining to cancer development, could exert anti-oncogenic activities by regulating tumor cell apoptosis, growth, and drug resistance (32,33). GBP2 is indispensable for the protective immunity against microorganisms (34). Moreover, up-regulation of GBP2 expression corresponds to a better prognosis of breast cancer patients, and might participate in T-cell defense against breast cancer (35). IFITM1, belonging to the IFN-induced transmembrane protein family, exhibits high expressions in tumor tissues and cells, and it is an independent prognostic biomarker for patients suffering from tumors including gallbladder carcinoma, esophageal adenocarcinoma, colorectal cancer, and gastric cancer (36). According to a report, OAS2 exists in the patients with malignant diseases (37). It appears that the high mRNA expression of OAS2 represents more favorable outcomes in breast cancer patients (38). In our CDIGPM, the coefficients of MX1 and OAS2 were positive numbers, while those of RSAD2, IRF2, GBP2 and IFITM1 were negative numbers. Therefore, CDIGPM is negatively correlated with RSAD2, IRF2, GBP2 and IFITM1 but positively correlated with MX1 and OAS2. Since the absolute value of IRF2 coefficient is the largest (-0.419071422539304), we wonder whether it affects the prognosis by regulating the malignant behavior of tumor cells. Results showed that down-regulation of IRF2 promoted SKCM cell proliferation. It has been reported that IRF2 inhibits cell proliferation by inducing CLDN7 upregulation in oral squamous cell carcinoma (39), and inhibits cancer proliferation by promoting AMER-1 transcription in human gastric cancer (40). What's means that CDIGPM might affects the prognosis of SKCM patients by regulating the malignant behavior of tumor cells. In summary, CDIGPM is a model that is related to the prognosis and tumor immunotherapy. Then, we explored the correlation of CDIGPM with existing predictive markers including CD274 (PD-L1) and CTLA4 for immunotherapy. CD274 + and CTLA4 + tumors tend to respond better to ICT than negative tumors (41)(42)(43). Here, we found a negative correlation between CDIGPM score and CTLA4 and CD274.
Exploring the TME might be helpful for finding new methods for the immunotherapy of SKCM. Between the two CDIGPM subgroups, there were differences in the activity of immune functions and the infiltration of certain immune cells. More CD8 + T cells, activated CD4 + memory T cells, and M1 macrophages were existed in low-CDIGPM group, while more M2 macrophages were existed in high-CDIGPM group. It has been reported that more infiltration of CD8 + T cells is related to a good prognosis in cancers (44-46). Activated M1 macrophages can trigger adaptive immune responses. M2 macrophages play an immunosuppressive role and exert a tumor growthpromoting effect (47, 48). All these mean that low-CDIGPM group has better tumor immunity potential, while high-CDIGPM group has immunosuppressive characteristics.
Aimed to explore its predictive value in cancer immunotherapy, we analyzed immunotherapy data GSE35640. We found that the CDIGPM score in patients who respond to immune therapy was lower than it in patients who did not respond to immune therapy. These results indicate that CDIGPM might be a prediction model for the effect of cancer immunotherapy.
Nevertheless, the study still has shortcomings. Most importantly, prospective studies are needed to further confirm the value of this prognostic model.
In total, CDIGPM is an encouraging prognostic model. It may help identify immunologic features and predict the prognosis of SKCM patients. Meanwhile, CDIGPM might have predictive value for immune escape.

Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi. nlm.nih.gov/geo/, GSE120575.

Author contributions
QL and YS made important contributions to the study conception and design. SS, JS and ZZ conducted data analysis and interpretation. All authors contributed to the article and approved the submitted version.

Funding
This research was supported by Natural Science Research Plan of Huaian (No. HAB202022).

Publisher's note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

SUPPLEMENTARY FIGURE 2
The expression levels of CDIGPM genes in SKCM. The expression levels of MX1, RSAD2, IRF2, GBP2, IFITM1, and OAS2 in SKCM were explored using GEPIA. Based on TCGA and GTEx data, the expression levels of MX1 and OAS2 are higher in SKCM cells. CDIGPM: CD8 + T cells-associated immune genes prognostic model.

SUPPLEMENTARY FIGURE 3
Enriched KEGG gene sets in high-and low-CDIGPM groups. Graphical Abstract.