The Immune-Related PRF1 May Interact with Component 3 to Facilitate a Favorable Prognosis of Skin Cutaneous Melanoma

Background Previous studies have revealed the immune cytolytic activity (CYT) mediated by Perforin 1 (PRF1) is associated with a favorable prognosis of skin cutaneous melanoma (SKCM) and displays signicant differences between primary and metastatic SKCM. However, the role of PRF1 in SKCM is poorly understood. Methods By comparing the PRF1 co-expression network between metastatic and primary SKCM, we identied the potential hub genes interacting with PRF1 and then veried their correlation in lung metastatic melanoma cell lines. Survival analysis was performed for the four groups of SKCM dened by the median expression values of PRF1 and C3. Meanwhile, the possible relationship of PRF1 with immunotherapy response is investigated in the melanoma expression prole of an adoptive T cell therapy dataset. Results Almost all primary co-expression genes (n= 502, 96.91%) and more than half of metastasis-specic co-expression genes (n= 830) were identied in metastatic SKCMs. Functional enrichment analysis revealed that extensive CYT-related biological processes were disturbed in the metastatic SKCM, including immune response activation and inammation response. The hub gene C3 showed a positive correlation with PRF1 in multiple datasets. Prognosis analysis indicated that the group with high expression of both PRF1 and C3 showed the best prognosis. Moreover, high expression of PRF1 and C3 is associated with positive immunotherapy response. Besides, we identied ve human leucocyte antigen (HLA) genes as the central nodes of the PPI network. Conclusions Remarkable alterations in the CYT-related biological processes were observed between primary and metastatic SKCMs. Preliminary results demonstrate that PRF1 may act synergistically with C3 to facilitate the favorable prognosis of SKCM.


Background
Melanoma is a malignant form of skin cancer that develops from melanocytes, of which skin cutaneous melanoma (SKCM) is the most common type [1]. Recently, the incidence of SKCM presents a steady increase and has become one of the top 20 most common cancers, causing approximately 1.6% of all new cancers and 0.6% of all deaths each year even though the approved application of new treatment techniques has reduced its mortality [2]. Numerous complex factors contribute to the development of SKCM, whereas UV radiation exposure is generally considered to be the key risk factor [3]. Age is the other important factor related to the prognosis of SKCM, with the disease tending to occur more aggressively in older patients (> 70 years), having lower survival rates [4]. Juvenile or pediatric patients are usually diagnosed at an advanced stage, while the prognosis is better in comparison with adult patients. In addition, sex is also associated with SKCM incidence and men are more likely to develop melanoma than women [5]. Perforin 1 (PRF1), known as pore-forming protein, encodes a protein that forms pores in the cytoplasmic membrane, thus allowing cytotoxic lymphocytes or neutral killer (NK) cells to release granzyme (GZMA) in order to mediate the lysis of target cells [6]. PRF1 acts with GZMA synergistically to perform the immune cytolytic activity (CYT), an essential component of adaptive immunity, thus allowing cytotoxic T cells and NK cells to induce cytolysis through programmed cell death, aiming to remove the infected or pre-cancerous transformed abnormal cells [7]. Mutated or expression defected PRF1 can cause dysfunction of the immune system, leading to serious immune disorders such as familial hemophagocytic lymphohistiocytosis and macrophage activation syndrome [8]. Additionally, PRF1 also exerts an important role in provoking skin in ammation, as evidenced by its increased expression in chronic in ammatory skin diseases such as atopic dermatitis (AD) and psoriasis [9]. The relationship of altered PRF1 with primary cancers is not well-known. In a genotyping study, it was found that the genotype of PRF1 was not associated with increased susceptibility to malignant tumors such as melanoma, lymphoma, colorectal cancer and ovarian cancer [10]. In head and neck squamous cell carcinoma (HNSCC), PRF1 is recognized as a good prognostic predictor with an elevated expression, especially in HPV positive HNSCC [11].
Accumulated evidence has demonstrated that PRF1 plays an essential role in the control of viral infections and tumor surveillance. However, the role of PRF1 in SKCM is poorly studied. Extensive studies have focused on the CYT, which is denominated by the expression of the toxins PRF1 and GZMA [12]. It has been reported that CYT is associated with a favorable prognosis of SKCM, and the CYT score, de ned by average expression of PRF1 and GZMA, is signi cantly higher in metastatic SKCMs than that in primary SKCMs [13,14]. In the TCGA SKCM dataset, metastatic patients show a better prognosis than primary patients. This prompted us to propose a hypothesis that the different CYT between primary and metastatic SKCMs may have an important impact on the patients' prognosis, which can be mediated by shaping the adaptive immunity of tumor cells. By comparing the PRF1 co-expression network between metastatic and primary SKCM, we identi ed the potential key genes that interact with PRF1 in the network, and then veri ed their correlations in lung metastatic melanoma cell lines. Meanwhile, the possible relationship of PRF1 with immunotherapy response is also investigated in the melanoma expression pro le of an adoptive T cell therapy dataset.

Data collection and preparation
The gene expression pro les of SKCM were obtained from the Cancer Genome Atlas (TCGA) database (https://cancergenome.nih.gov/), which contains 456 primary (n=98) and metastatic (n=358) samples (Supplementary table 1). The overall survival time was shorter in primary patients than that of metastatic patients, and there was no signi cant difference between male and female proportions (P = 0.55). The other four expression pro les, GSE22155 (n=29) [15], GSE54467 (n=70) [16], GSE65904 (n=206) [17] and GSE100797 (n=25) [18] were retrieved from the Gene Expression Omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/). We merged GSE22155, GSE54467, and GSE65904 as an independent validation set, as they had small sample sizes and survival time, and were generated by Illumina HumanHT-12 platform. The Z-score method was used to normalize the expression pro les before merging to remove batch effects. GSE100797 contains 25 samples that derived from melanoma patients after adoptive T cell therapy, including 5 complete responders (CR), 5 partial responders (PR), 10 stable disease (SD) and 5 progressive disease (PD).

Co-expression analysis
Pearson correlation coe cients between PRF1 and all other genes were calculated in primary and metastatic SKCMs, respectively, meanwhile signi cant P values were estimated. Low-expression genes, de ned as the FPKM values below 1 in more than 50% samples, were removed before the correlation analysis to avoid impacting the assessment of correlation coe cient. Multiplex test was conducted to evaluate the false positive discovery rate using the BH method, and genes with FDR < 0.05 were screened as co-expression genes. Gene ontology (GO) enrichment and gene set enrichment analysis (GSEA) were performed for the co-expression genes.

Identi cation of hub genes
The co-expression genes were used to construct a protein-protein interaction (PPI) network by integrating the interaction information from STRING (Search Tool for the Retrieval of Interacting Genes) database [19] to identify the potential hub genes. The connections with score >=0.8 were selected as reliable protein-protein interactions. The Cytoscape software [20] was used to show the PPI network and the topological parameters of this network were estimated by the NetworkAnalyzer plugin [21]. The potential hub genes were determined based on the number of directed edges and betweenness centrality of nodes. Besides, the other four parameters, closeness centrality, clustering coe cient, average shortest path length and stress were also estimated.

Construction of recombinant vector
Since the PRF1 gene is barely expressed on the MALME-3M melanoma cell line, the gene sequence could not be cloned by a cDNA library. The vitro recombination technique were used to clone the full length of PRF1 cDNA. Brie y, two pairs of primers were designed to amplify two exons of PRF1 respectively (Supplementary table 2), and then the two ampli ed fragments were ligated to a complete PRF1 sequence. Then the complete PRF1 fragment was inserted into the pcDNA3.0-3xFlag vector to construct the recombinant vector.

Cell culture and transfection
Human Malme-3M melanoma cells was purchased from ATCC (https://www.atcc.org/products/htb-64) in May 2019 and authenticated by STR pro ling. Mycoplasma contamination was rst tested for the cells and then they were revived on Dulbecco's modi ed Eagle's complete medium (DMEM), and incubated in a cell incubator at 37°C and 5% CO 2 for 4-6 hours. The revived cells were then digested by trypsin for 3-5 min at 37°C. 10ul of the cell suspension were added to the blood cell counting plate for cell counting (2.5 * 10 ^ 6 / ml). Then, 80ul of cell suspension was added to a six-well plate, incubating in a cell incubator at 37°C and 5% CO 2 for 24 hours. Solution 1 was prepared by mixing 100ul of Opti-MEM medium (without antibiotics and serum) with 3ul of transfection reagent (Lipo3000). Solution 2 and solution 3 were prepared by mixing 100ul of Opti-MEM medium with 0.1ug and 1ug recombinant vector. Solution 1 was mixed with solution 2 and solution 3 respectively, and stand for 15 min. Then, the mixed solutions and 800 ul of preheated DMEM medium (serum free and double antibody free) were transferred to the six-well plate that contained revived cells. The cells were rstly incubated in a 37°C cell culture incubator for 4-6h. The medium was subsequently replaced by complete medium (containing 10% fetal bovine serum and 1% double antibodies) and after 48h, the cells were collected for RNA extraction.
Quantitative PCR Real-time PCR was used to quantify the expression of PRF1 and C3 in transfected cells. Total RNA was extracted from the collected cells using Trizol reagent according to the instructions. The HiScript Q Select RT SuperMix for qPCR kit (vazyme, Nanjing, China) and SYBR Green Master Mix kit (vazyme, Nanjing, China) were utilized for cDNA synthesis and quantitative PCR, respectively. Real-time qPCR was performed in the BioRad CFX96 PCR system. The primers used in the study were listed in supplementary table 3. GAPDH was selected as an internal reference, and 2 -ΔΔCt was calculated to quantify the relative expressions of PRF1 and C3.

Survival analysis
Firstly, samples with survival time less than 30 days were removed. For categorical variables, the log rank test (R package 'survival') was used to evaluate the prognostic survival difference among groups. For continuous variables, the samples were divided into high-expression and low-expression groups according to the median expression value, and then the prognostic survival difference between the two groups were assessed. The survival curve and the corresponding 95% con dence interval was plotted in GraphPadPrism 8 (v8.0.0) software. Each point in the survival curve represents a censored value.

Statistical analysis
The data preparation, processing and statistical analysis were performed on R software (v3.6.2). Survival curves were plotted on GraphPadPrism software (v9.0.0). For continuous variables, wilcox rank sum test and krusk Wilma test was conducted for the comparisons between unpaired two groups and multiple groups, respectively. The chi-square test was used for the comparisons of category variables between groups. Spearman's rho statistic is used to estimate the correlation of two continuous variables.

Results
The association of PRF1 with SKCM prognosis The primary SKCM patients had shorter overall survival time the metastatic SKCM patients (P < 0.001, Supplementary table 1) ( Figure 1A) with PRF1 expression upregulated in metastatic samples ( Figure 1B). Patients in the PRF1 high expression group showed a better prognosis than that in the low expression group ( Figure 1C). As one of the key marker genes of immune cytolytic activity (CYT), the PRF1 expression exhibited a strong positive correlation with lymphocyte in ltration score ( Figure 1D). These results re ected a tight association between PRF1 and the metastasis and prognosis of SKCM, implying the remarkable varied CYT between primary and metastatic SKCMs.
The CYT-related biological processes in primary and metastatic SKCM We identi ed 518 and 1332 co-expression genes in primary and metastatic SKCMs, of which 502 genes were shared by both (Supplementary gure 1), accounted 96.91% and 37.69% of the total co-expression genes respectively. Functional enrichment analysis revealed that both the primary and metastatic coexpression genes were enriched in the processes of T cell activation, lymphocyte activation, etc. ( Figure   2A). For primary speci c co-expression genes, they were mainly related to tumor necrosis, while immune response activation and in ammation response were the major enriched processes for metastatic speci c co-expression genes ( Figure 2B) Figure 2C).

Identi cation of C3 as the hub gene in PPI network
Since extremely few primary speci c co-expression genes (n=16) were identi ed, the metastatic speci c co-expression genes (n=830) were selected to construct a protein-protein interaction network based on the StringDB database. The network consists of 948 edges and 257 nodes ( Figure 3A). The topological parameters of PPI network estimated by the NetworkAnalyzer plugin in Cytoscape indicated that C3 gene had the largest degree, the number of node connections ( Figure 3B). The betweenness centrality of C3, which re ects the amount of control the node exerts over the interactions of other nodes in the network, ranked second ( Figure 3C). In addition, the other four parameters including closeness centrality, clustering e ciency, average shortest path length and stress, also suggested the key role of C3 in the network (Supplementary gure 2).

The association of PRF1 with C3
The expression of the potential hub gene, C3, was found upregulated in the metastatic SKCMs ( Figure  4A), and its high expression group exhibited a better prognosis than the low expression group (Figure 4B), suggesting that C3 might be a favorable prognostic predictor. Besides, the expression of PRF1 and C3 presented a moderate correlation ( Figure 4C). Then, we divided SKCM samples into four groups, PRF1high/C3-high, PRF1-high/C3-low, PRF1-low/C3-high, and PRF1-low/C3-low according to the median values of PRF1 and C3 expressions. Survival analysis indicated that PRF1-high/C3-high group showed the best prognosis and PRF1-L/C3-L group had the worst prognosis ( Figure 4D). Additionally, most of the PRF1-L/C3-L group samples were primary SKCMs, while the PRF1-H/C3-H group samples were dominated by metastatic SKCMs (Supplementary Figure 3).

Validation of PRF1 and C3 in independent set
In the TCGA SKCM dataset, the potential synergy between PRF1 and C3 was observed that may contribute to a good prognosis of this disease. Then, we veri ed this relationship in an independent GEO dataset. Due to their small sample sizes, three GEO datasets including GSE22155 (n= 29), GSE54467 (n= 70) and GSE65904 (n= 206) were integrated into one set by employing Z-score method. The combined dataset showed no batch effect after Z-score normalization (Supplementary gure 4). The expression of PRF1 was positively associated with C3 in the validation set ( Figure 5A), which was in line with that in TCGA dataset. Four groups divided by median values of PRF1 and C3 expressions also presented a signi cant difference in prognosis, with PRF1-H/C3-H group having the best prognosis and PRF1-L/C3-L group the worst prognosis ( Figure 5B).

Over-expression of PRF1 in melanoma cell lines
The expression pro les of Cancer Cell Line Encyclopedia (CCLE) revealed that PRF1 and C3 were almost not expressed in SKCM cell lines. As a result, we attempted to overexpress PRF1 on the lung metastasis cell line, MALME-3M, and quantify the C3 expression to investigate their possible relationship. A recombinant vector containing the full fragment of PRF1 was rst constructed based on the recombination clone technology (Supplementary gure 5). Then, the empty vector, 0.1ug and 1ug recombinant vector were transfected to MALME-3M cells. The relative expressions of PRF1 and C3 were quanti ed by qPCR with GAPDH as the control gene. Results indicated that the PRF1 expression was signi cantly up-regulated in cells transfected by recombinant vectors, and the expression of C3 was also increased ( Figure 6A). The relative expression values of both were associated with the concentration of transfected recombinant vectors (Supplementary Table 5). Moreover, the relative expressions of PRF1 and C3 also showed a strong positive correlation ( Figure 6B).
The association of PRF1/C3 with immunotherapy PRF1 is well-known as the marker of immune cytolytic activity, and C3 is related to the in ammatory response of complement system. In the expression pro les of melanoma patients received adoptive T cell therapy, we investigated the expressions of PRF1 and C3 in four immune response groups including complete responders (CR, n=5), partial responders (PR, n=5), stable disease (SD, n=10), and progressive disease (PD, n=5). The PD group showed the lowest expressions of PRF1 and C3 in comparison with the other three groups ( Figure 7A). A strong correlation between PRF1 and C3 was observed among the four group samples ( Figure 7B).

Discussion
Although studies of CYT de ned by PRF1 and GZMA have been reported in SKCM, the related biological pathways remain unclear, especially the relationship between CYT and metastatic SKCM prognosis is rarely involved. Most of the previous studies are focused on the impact of dysfunctional PRF1 caused by genomic variations on diseases [22]. In the current study, we rstly investigated the alterations of PRF1related biological processes between primary and metastatic SKCMs, and then identi ed C3 as the potential hub gene of the PPI network constructed by PRF1 co-expression genes. Further investigations indicated that the PRF1/C3 axis was associated with favorable immunotherapy response and prognosis. Finally, the relationship between PRF1 and C3 was veri ed in an independent validation set and malignant melanoma cell lines.
Almost all of the co-expression genes (96.91%) in primary SKCMs were also found in metastatic SKCMs and more metastasis-speci c co-expression genes (62.31%) were identi ed, indicating that a few of new biology processes were stimulated in metastatic SKCMs, except for the pathways shared by both [12]. These results were con rmed by functional enrichment analysis. Both of the primary and metastatic coexpression genes were enriched in the regulation of lymphocyte activation, T cell activation, and leukocyte activation. However, the speci c co-expression genes of both sides showed dramatic differences in BP terms, re ecting the extensive CYT-related biological processes involved in the metastasis of melanoma cells from local to distal. Although the metastasis of melanoma is a highly complex event such as the genomic mutations, the interaction between extracellular matrix and tumor cells, etc. [23,24], our ndings provide a new insight into the perspective of CYT.
Component 3 that was recognized as the hub gene of the PPI network, showed a positive correlation with PRF1, suggesting a possible synergistic interaction between them in metastatic SKCMs. The subsequent prognostic analysis revealed that PRF1-H/C3-H group had the best prognosis and the PRF1-L/C3-L group the worst prognosis, which was veri ed in a validation set. Overexpression of PRF1 in the malignant melanoma cell lines increased the expression of C3. These results demonstrated the potential interaction between PRF1 and C3. Indeed, C3 is one of the well-known pro-in ammatory molecules that are upregulated in most in ammatory diseases [25,26]. The CYT-related modules are positively correlated with in ammatory response, suggesting a stronger in ammatory response in tumors compared to normal tissues. Besides, it has also been reported that C3 is required for the metastasis of lung and breast cancer cell lines to the soft meningeal space [27], which is line with our observations that C3 was up-regulated in metastatic SKCMs.
The other ve genes with highest degrees were all human leucocyte antigen (HLA) group genes (HLA-DRB1, HLA-DPA1, HLA-DRB5, HLA-DQA2 and HLA-DQA1). HLA group genes are encoded by human major histocompatibility complex (MHC) genes characterized by high polymorphism [28]. These molecules exert an important role in regulating the host immune response, mainly participated in the presentation of antigens to T-cell receptors (TCRs) [29]. Our results revealed that the ve HLA genes were the central nodes of the PPI network, indicating their extensive interactions with other PRF1-related processes. These ndings suggested a high immune activity in the immune microenvironment of metastatic SKCMs, which might contribute to a favorable prognosis of this disease.
In melanoma patients received the adoptive T cell therapy, we observed that high expression of PRF1 and C3 was associated with a good therapy response, and the two genes also showed a strong correlation, re ecting a synergistic interaction between them, which may bene t a positive response to adoptive T cell therapy. The relationship between T cell densities at the margins and centers of invasive tumor was investigated in colorectal cancer, and high granzyme density was found to be associated with a lower probability of tumor recurrence and prolong overall survival, implying that high activity of CYT can be a negative regulator of tumor recurrence [30]. The acquisition of cytolytic property is a key step for activated T cells (CD4 + T cell and CD8 + T cell) to exert cytotoxic function [31]. Cytolytic property plays an important role in the prognosis of melanoma patients, for patients with high cytolytic score having a better prognosis. Therefore, it is recognized as one of the functional effectors of toxic T cells [13]. In addition, the upregulated PRF1 was also observed in anti-CTLA-4 and anti-PD-L1 immunotherapy response groups [32]. Previous studies and the current study imply that PRF1/C3 may serve as a predictor of immunotherapy in melanoma patients.
However, some issues still need to be addressed in the future. Firstly, it remains unclear whether there is a causal relationship between PRF1 and C3, because the current results only revealed a positive correlation between them. Although overexpression of PRF1 increased the C3 expression, more thorough molecular experiments are required to investigate their potential relationship. Second, the synergistic effect of PRF1 and C3 cannot tell us if there exists a direct or indirect interaction between them. Third, the prognosis of melanoma is affected by numerous complex factors, and some other pathways may interact with PRF1/C3 axis to regulate the patient prognosis.

Conclusions
In summary, we observed vast variations between primary and metastatic SKCMs in the CYT-related biological processes. C3 was identi ed as the hub gene in the PPI network by comparing the coexpression genes of PRF1 in primary and metastatic SKCMs. The subsequent prognostic survival analysis and cytological experiments provide preliminarily evidence that C3 may act synergistically with PRF1 and thus contributed to the favorable prognosis of SKCM. We further observed a correlation between the PRF1/C3 expression and positive immunotherapy response, indicating their potential to serve as the predictors of immunotherapy.

Declarations
Ethics approval and consent to participate All data used in this study are publicly accessible to any researcher, our research did not require the approval of an ethics committee.

Consent for publication
Not applicable.

Availability of data and materials
The gene expression pro les are publicly accessible in TCGA and GEO database. All the intermediate data and results are presented in the manuscript and supplementary materials.     Overexpression of PRF1 in melanoma cells. A: Relative expression of PRF1 and C3 in MALME-3M cells. B: The correlation of PRF1 and C3 relative expressions.