A TP53-associated gene signature for prediction of prognosis and therapeutic responses in lung squamous cell carcinoma

ABSTRACT

The tumor-suppressor gene tumor protein p53 (TP53) is one of the most commonly mutated genes in human lung cancer, and TP53 mutations are associated with a worsened prognosis and causes resistance to cancer therapy. RNA sequencing and TP53 mutation data were downloaded to determine specific TP53-associated signature based on differentially expressed genes between patients with lung squamous cell carcinoma (LUSC) with wild type (TP53 ^WT) and mutated (TP53^MUT) TP53. We investigated the predictive value of this signature on the immune microenvironment, tumor mutational burden (TMB), and likelihood of response to immunotherapy and chemotherapy. In total, 1,556 differentially expressed genes were identified based on TP53 mutation status. Three genes (KLK6, MUC22 and CSN1S1) identified by univariate and multivariate Cox regression analyses, comprised the prognostic signature which was an independent and specific prognostic marker of overall survival in patients with LUSC. A nomogram was also established to validate this signature for clinical use. Moreover, the high-risk group was characterized by increased infiltration of monocytes and macrophages M1, and decreased T cells CD8 and T cells follicular helper. High-risk group exhibited a higher TMB, and was much more sensitive to immunotherapy and chemotherapy. KLK6 and CSN1S1 expression and the prognostic prediction values were further validated in clinical samples. The derived TP53-associated signature is a specific and independent prognostic biomarker for LUSC patients, and could provide potential prognostic biomarker or therapeutic targets for the development of novel immunotherapies and chemotherapies.

KEYWORDS:

• TP53
• lung squamous cell carcinoma
• signature
• prognosis
• immune/chemotherapy

Introduction

Lung cancer is one of the common human malignancies and is the leading cause of cancer-associated deaths world-wide.¹ There are three different subtypes of non-small cell lung cancer (NSCLC): (1) lung adenocarcinoma (LUAD); (2) lung squamous cell carcinoma (LUSC); and (3) large cell lung carcinoma (LCLC). LUSC constitutes approximately 30% of all lung cancer cases and globally, results in approximately 400,000 deaths per year.^2,3 Though great progress has been made toward its prevention, diagnosis and targeted treatment, the clinical outcome of lung cancer remains unsatisfactory. Growing evidence indicates that the malignant phenotypes of cancers are influenced by the tumor microenvironment.^4–6 Lung cancer, an immune-sensitive malignancy, is infiltrated by numerous immune cells, including macrophages, eosinophils, neutrophils, dendritic cells, mast cells, natural killer cells, B cells, and T cells. However, few studies have systematically explored the relationship between the immune phenotype of the lung cancer microenvironment and its prognosis.

Tumor protein p53 (TP53), located on chromosome 17p13.1, encodes the tumor suppressor protein p53.⁷ It is the most commonly mutated gene in humans with more than 50% of all human cancers indicating alterations in p53 signaling.⁸ Functionally, p53 binds directly to chromatin, and plays important roles in the regulation of the cell cycle, apoptosis, autophagy, and DNA repair in response to damaging agents.⁹ However, when mutated, loss of these functions results in abnormal cell proliferation and promotes cancer.¹⁰ TP53 mutations, mainly missense mutations, are also the most common mutations in NSCLC, and are more prevalent in LUSC than in LUAD.^11–13 Epidemiological studies have reported that TP53 mutations are closely related to smoking, and it occurring more frequently in patients with tobacco-associated lung cancer than in never-smokers.^14,15 An increasing body of literature suggests that TP53 mutations in lung cancer are associated with increased resistance to cancer therapies and poorer survival prognosis.^16–18 In addition, TP53 mutations are associated with higher vascular endothelial growth factor (VEGF) synthesis and angiogenesis.¹⁹ Recently, TP53 mutation status was associated with cancer-related microenvironment.^20,21 We hypothesized that the overall survival of patients with LUSC harboring TP53 mutations might be particularly influenced by the lung cancer microenvironment. Therefore, we identified genes affected by TP53 mutation status, and established a three-gene gene signature that is a robust prognostic biomarker and predictive factor that can be used in the clinic.

Materials and methods

Data sources

VarScan 2-based somatic mutation data from patients with LUSC and LUAD, combined with gene expression data and corresponding clinical features, were accessed from the Cancer Genome Atlas (TCGA) website. This study meets TCGA’s publication guidelines. All LUSC gene expression, clinical, and somatic mutation data were downloaded through the Data Coordinating Center. We also downloaded somatic mutation data from the International Cancer Genome Consortium (ICGC) to estimate the somatic mutations of patients with LUSC.

Screening of differentially expressed genes (DEGs)

First, the raw counts of gene expression data from TCGA were normalized using a weighted trimmed mean of log ratios-based method.²² To obtain DEGs between patients with (n = 388) and without (n = 100) TP53 mutations in the TCGA LUSC cohort, the R package “edgeR” was used in the standard comparison mode.²³ The DEG threshold was set at a |log2 fold change| ≥ 1 and a false discovery rate < 0.05.

Gene set enrichment analysis (GSEA)

To identify potential differences in biological functions between LUSC patients with and without TP53 mutations, GSEA annotation was performed using the R package “clusterProfiler”.^24,25 The GSEA threshold for significantly enriched functional annotations was set at a p-value < 0.05 and an enrichment score > 1.0.

Identification of TP53-related prognostic signature

Univariate Cox regression analysis was performed using the R package “survival” to evaluate correlations between the DEG expression levels and the overall survival of patients with LUSC. DEGs with p-value < 0.05 by univariate Cox regression analysis were identified as prognostic. Multivariate Cox regression analysis was then used to determine the prognostic values of specific gene signatures. A three-gene-based prognosis risk score was calculated to assess each patient’s risk was calculated according to the following formula:

$\mathrm{R}\mathrm{i}\mathrm{s}\mathrm{k}\mathrm{s}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e}=\sum _{i=1}^{\mathrm{N}}\left({\mathrm{E}\mathrm{x}\mathrm{p}}_i\times {\mathrm{C}}_i\right)$

in which N represents the number of prognostic genes, Exp_i represents the expression of gene_i profile and C_i represents the estimated regression coefficient of gene_i determined by multivariable Cox regression analysis.^26–28 Patients with LUSC with available survival data were separated into high- and low-risk groups using the median score as a cutoff. Survival analysis was performed using the Kaplan-Meier method and the log-rank test was applied to evaluate the statistical significance of the differences.

Estimation of tumor-infiltrating immune cells

We uploaded normalized gene expression data with standard annotation files to the Cell type Identification by Estimating Relative Subsets of RNA Transcripts (CIBERSORT) web portal, and the fractions of 22 immune cell types using 1,000 permutations and the LM22 gene signature as previously described.²⁹ The R package “Genefilter” was used to screen each sample, and the threshold was set at p-value < 0.05. The final CIBERSORT output was subsequently analyzed.

Immunotherapeutic and chemotherapeutic response prediction

The programmed cell death 1 (PCDC1, also known as PD-1)/CD274 molecule (CD274, also known as PD-L1) and cytotoxic T-lymphocyte associated protein 4 (CTLA-4) pathways in are implicated in tumor immune evasion, and therefore, immune checkpoint inhibitors targeting PD-1 and CTLA-4 thereby enhance antitumor immunity.³⁰ We used the Tumor Immune Dysfunction and Exclusion (TIDE) algorithm and subclass mapping to predict clinical responses to immune checkpoint inhibitors as previously described.^31,32 As chemotherapy is commonly used to treat NSCLC, we used the R package “pRRophetic” to estimate the chemotherapeutic response determined by the half maximal inhibitory concentration (IC50) of each LUSC patient on the Genomics of Drug Sensitivity in Cancer (GDSC) website.³³

Estimation of tumor mutational burden (TMB)

The TMB (in mutations per megabase), is an emerging therapeutic measure of sensitivity to immunotherapy. The TMB score of each patient with LUSC was calculated as previously described.³⁴

Independence of the gene-related prognostic signature from other clinical features

To determine whether the predictions of the prognostic signature were independent of the clinical characteristics (including age, gender, tumor/node/metastasis (TNM) stage, T stage, N stage, and M stage) of patients with LUSC, univariate and multivariate Cox regression analyses were performed, and nomograms were constructed to assess the probability of 1-, 3-, and 5-year overall survival for LUSC patients based on the signature.

Immunohistochemistry

Tissue microarray chips containing 90 samples of LUSC and 90 samples of paired normal lung tissue were purchased from Outdo Biotech (Shanghai, China). Immunohistochemistry (IHC) staining was performed as previously described.^35,36 Microarray chips were stained with PD-L1 (Gene Tech), KLK6 (R&D) and CSN1S1 (Bioss). Quantitative analysis of the staining was independently assessed based on the percentage of positive cells and staining density by two experienced pathologists.³⁷ The final score was determined by adding the staining intensity score and the average proportion of positive cells score.

Statistical analysis

All statistical analyses were performed using R version 3.6.1, and the data from different groups were compared by Mann-Whitney-Wilcoxon Test. Pearson’s chi-square test was performed to measure the level of significance for association amongst variables. All reported P values were two-tailed, and p < .05 was considered statistically significant.

Results

Mutations in LUSC and LUAD

Traditionally, lung cancer treatment decisions have been based on histological considerations. In the last few years, novel insights in tumor biology and the opportunity to identify genetic alterations have rapidly changed the process of therapeutic selection. We initially sought to identify somatic mutations in patients with LUSC and LUAD. According to TCGA, TP53 mutations were the most frequent, and were more prevalent in LUSC than LUAD (77% vs. 47%; Figure 1). We also identified LUSC mutations in the ICGA database. Consistently, TP53 was also the most frequently mutated gene (ranked second), which was consistent with its high frequency in the TCGA database (Supplemental Figure 1).

Figure 1. Mutations in LUSC and LUAD samples (a) Overview of somatic mutations in all samples in the (A) LUSC and (b) LUAD TCGA cohorts

Identification of DEGs between LUSC patients with and without TP53 mutations

The high rate of TP53 mutations indicated that TP53 status was closely linked to LUSC. TP53 mutation status is a well-known clinically relevant molecular marker in lung cancer.³⁸ Therefore, we separated LUSC patients into TP53 mutated and wild-type groups and explored DEGs between them. In total, 773 upregulated genes and 783 downregulated genes were identified (Figure 2(a,b)). To gain insight into DEG functions, we performed gene ontology (GO) enrichment analysis based on GSEA analysis. As a result, TP53 mutation status genes were clustered most enriched for terms related to immune functions, such as major histocompatibility complex (MHC) class II protein complex, establishment of T cell polarity, immunoglobulin complex, and circulating and immunoglobulin complex indicating that genes related by TP53 mutation status mainly probably played an important roles in immune-related processes in LUSC.

Figure 2. Identification of DEGs in patients with LUSC with and without TP53 mutations (a) Heatmap and (b) volcano plot of identified DEGs. (c) GSEA analysis of samples with and without TP53 mutations

Construction of the DEG-based prognostic signature

To identify DEGs that correlated with the overall survival patients with LUSC, univariate Cox regression analysis was performed. At p < .05, 75 genes were identified as prognostic of overall survival in patients with LUSC (Supplemental Figure 2). Multivariate Cox regression analysis was then applied to inspect the interrelated relationship among genes with overall survival, and only three genes showed significant prognostic signatures for LUSC, containing kallikrein related peptidase 6 (KLK6), mucin 22 (MUC22), and casein alpha s1 (CSN1S1; Figure 3(a)). A risk score to predicting prognostic value was calculated, and patients with LUSC were separated into low- or high-risk groups (Figure 3(b–d)). Kaplan-Meier analysis indicated that high-risk patients had significantly worse overall survival than low-risk patients (p < .0001; Figure 3(e)).

Figure 3. Construction of the DEG-based prognostic signature (a) Prognostic values of DEGs by multivariate Cox regression analysis. (b) The distribution of the three-gene-based risk score. (c) Vital statuses of patients in the high‐ and low‐risk groups. (C) Heatmap of three-gene expression profiles in the high‐ and low‐risk groups; (d) Kaplan-Meier survival curves of the relative overall survival of high- and low-risk patients

Kaplan-Meier analyses of overall survival according to TP53 mutation status

Consistent with the prognostic capacity of the three-gene signature, TP53 mutation status was also significantly correlated with the prognosis of patients with LUSC (Figure 4(a)). To investigate whether the three-gene signature was independent of TP53 mutation status, patients with LUSC were divided into high- and low-risk groups based on TP53 mutation status. Kaplan-Meier overall survival curves of the two groups based on the three-gene signature were significantly different in the TP53 ^WT and TP53^MUT LUSC cohorts (Figure 4(b,c)). Different types of TP53 mutations can occur and we found that the TP53 mutation type affects the prognosis of LUSC patients (Figure 4(d)). To explore whether the three-gene signature was independent of TP53 mutation type, we performed prognostic analysis of the largest subgroup, which contained TP53 missense mutations. Interestingly, the three-gene signature played a negative role, with high expression of these genes shortening the overall survival of patients with LUSC (Figure 4(e)).

Figure 4. Kaplan-Meier analysis of overall survival according to TP53 mutation status Kaplan-Meier survival by TP53 status (a); in the TP53 mutation subgroup (b); in the TP53 wild-type subgroup (c); with different types of TP53 mutations (d); and in the TP53 missense mutation subgroup (e)

Immune cell infiltration landscapes of high- and low-risk patients with LUSC

We next investigated differences in immune infiltration between high- and low-risk patients with LUSC. As shown in Figure 5(a,b) the proportions of 22 tumor-immune cell types were significantly different between high- and low-risk patients with LUSC. Moreover, immune cell proportions were weakly to moderately correlated (Figure 5(c)). High-risk patients with LUSC had significantly higher proportions of monocytes and macrophages M1 and lower proportions of T cells CD8 and T cells follicular helper. We also calculated TMB scores based on TGCA somatic mutation data. The TMB was significantly higher in patients with LUSC (Supplemental Figure 3). The heterogeneity of immune infiltration observed in LUSC may provide prognostic indicators and targets for immunotherapy and could have significant clinical implications.

Figure 5. Immune cell infiltration landscapes in high- and low-risk patients with LUSC (a) Relative proportions of immune cell infiltration in high- and low-risk patients. (b) Correlation matrix and (c) heatmap of the 22 immune cell proportions. (d) Differences in immune cell infiltration abundances between high- and low-risk patients

Immunotherapeutic and chemotherapeutic responses of high- and low-risk patients with LUSC

Immune checkpoint blockade using immunotherapies, targeting CTLA-4 and PD-1, has emerged as a promising strategy to treat a variety of malignancies.³⁰ Thus, we estimated the clinical response to immune checkpoint blockade (targeting CTLA-4 and PD-1 in high- and low-risk patients with LUSC). Interestingly, high-risk patients with LUSC showed more promise in response to anti-PD-1 therapy (Bonferroni-corrected p = .005; Figure 6(a)). We also investigated the response to chemotherapy in high- and low-risk patients with LUSC, and found that 29 chemotherapeutic drugs displayed significant differences in estimated IC50 between high- and low-risk patients with LUSC, and that high-risk patients with LUSC showed increased sensitivity to all 29 chemotherapies (Figure 6(b)).

Figure 6. Immunotherapeutic and chemotherapeutic responses in high- and low-risk patients with LUSC (a) Immunotherapeutic responses to anti-CTLA-4 and -PD-1 treatments in high- and low-risk patients. (b) Differential chemotherapeutic responses in high- and low-risk patients

Correlations between the three-gene signature and clinical characteristics

We next investigated whether the three-gene signature-derived risk score was an independent biomarker with regards with conventional clinical information. Univariate Cox regression analysis indicated that the TNM stage, T stage, and risk score were correlated with poorer survival in patients with LUSC, and multivariate Cox regression analysis indicated

that the risk score could be used as a specific prognostic indicator in LUSC (p < .001; Figure 7(a,b)). A prognostic nomogram was established according to the factors in the multivariate analysis. Compared with the TNM stage and T stage, the risk score displayed superior predictive performance (Figure 7(c)).

Figure 7. Correlation between the three-gene signature and clinical characteristics (a, b) Univariate and multivariate Cox regression analyses of correlations between the three-gene signature and clinical characteristics with overall survival. (c) Nomogram for predicting the 1-, 3-, and 5-year overall survival of patients with LUSC

Validation of the gene signature in clinical tissue samples

To confirm the reliability of the identified gene signature, we examined KLK6 and CSN1S1 expression levels by immunohistochemistry using a tissue microarray containing 90 pairs of LUSC tissues and adjacent normal tissues. The results showed that KLK6 and CSN1S1 proteins were significantly overexpressed in tumor tissues when compared with those in normal tissues (Figure 8(a,b)). Moreover, to determine whether KLK6 and CSN1S1 could be prognostic biomarkers for patients with LUSC, survival analysis was performed, and the results showed that KLK6 and CSN1S1 expression levels were critically associated with overall survival in patients with LUSC (Figure 8(c)). We also assessed the correlations between the expression of gene signature (KLK6 and CSN1S1) and PD-L1, and found that KLK6 and CSN1S1 were significantly related to the expression of PD-L1 (Figure 8(d)).

Figure 8. Validation of the gene signature in clinical tissue samples (a) Representative images of IHC staining for CSN1S1 and KLK6 expression in LUSC tissues and adjacent non-tumor tissues. (b) Protein expression scores in LUSC tissues and normal lung tissues. (c) Kaplan–Meier analysis of overall survival according to CSN1S1 and KLK6 expression levels. (d) Correlation of PD-L1 expression with KLK6 and CSN1S1 expression

Discussion

Lung cancer is the most common cause global cancer mortality,³⁹ and TP53 mutated lung cancer displays an increased mutational burden, increased expression of immune checkpoint proteins, increased T cell infiltration, and PD-L1 amplification and derives remarkable clinical benefit from PD-1 inhibitors.²¹ However, the mechanisms by which TP53 mutation affects the microenvironment and lung cancer prognosis. Therefore, it is vital to further elucidate the immune-related effects of TP53 mutation status.

For the first time, this study has identified immune-related genes affected by TP53 mutation status, providing novel prognostic biomarkers and therapeutic targets for LUSC. Generally, patients in the TP53 mutation group had worse clinical outcomes. Genes modulated according to TP53 status were specifically enriched for GO terms related to immune response. We established a three-gene signature to predict LUSC prognosis, and separated patients with LUSC into two groups based on a related risk score. High-risk patients had worse overall survival.

According to the tumor immunoediting hypothesis, less immunogenic cancer cells are selected for during tumor development in immune-competent hosts, to evade antitumor immune responses.⁴⁰ This may result in increased immunosuppressive cells (e.g., regulatory T cells and tumor-associated macrophages), decreased immunoreactive cells (e.g., T cells follicular helper), and increased expression of immunosuppressive molecules (e.g., CTLA-4 and PD-1) in the tumor.⁴¹ PD-1 is a central regulator of T cells CD8 exhaustion and blockade of this inhibitory pathway enhances T cell immunity in several different types of cancers.^42,43 Thus, we hypothesized that patients in different groups would have different immunotherapeutic responses. As expected, we found that high-risk LUSC patients generally had higher monocytes and macrophages M1 infiltration and lower fractions of T cells CD8 and T cells follicular helper than low-risk patients. Moreover, we found that high-risk LUSC patients had higher TMB values than low-risk patients. Higher TMB values were correlated with increased probability of a favorable immunotherapy.^44,45 Consistently, TIDE prediction, indicated that high-risk LUSC patients were more likely to respond to anti-PD-1 therapy. Chemotherapy is the most common method used to treat lung cancer, and high-risk patients with LUSC were more sensitive to 29 chemotherapies than low-risk patients. The above results suggest that the poorer prognosis for high-risk patients with LUSC is due to higher immunosuppression and lower immunoreactivity in the tumor microenvironment, and that these differences promote tumor growth, progression, invasion, and metastasis. Importantly, due to these differences, high-risk patients with LUSC may derive greater benefit from immunotherapy and chemotherapy.

We found that the clinical TNM stage, T stage, and risk score significantly affected the overall survival of patients with LUSC, and determined that the three-gene signature was an effective independent prognostic model for LUSC. In routine clinical practice, the pathologic stage is a vital prognostic determinant of lung cancer. However, clinical outcomes differ among patients at the same stage, indicating that current staging systems are insufficient for effective prognosis, and cannot fully reflect the biological heterogeneity of patients with lung cancer. Therefore, it is essential to find potential biomarkers to use as prognostic and therapeutic indicators. To the best of our knowledge, this is the first prognostic gene signature related to TP53 mutation status identified. This gene signature provides a novel method to evaluate patients with LUSC, and guide prognostic prediction and therapeutic decisions. The three-gene signature could even distinguish the prognoses of patients with LUSC with different types of TP53 mutations. Importantly, we provide a nomogram that can be used to combine the risk score with LUSC clinical data to predict patient outcomes, which demonstrated that the risk score was a better signature for predicting the short- and long-term survival of patients with LUSC.

The expression of KLK6 and CSN1S1 and the prognostic prediction value was further validated in clinical samples. KLK6 and CSN1S1 proteins were significantly overexpressed in tumor tissues when compared with those in normal tissues. Survival analysis showed that high expression of these proteins was related to poor prognosis. The results might further prove the reliability of our analysis.

Overall, for the first time, this study identifies a three-gene signature associated with TP53 mutation status that can independently predict the survival of patients with LUSC. High-risk patients with LUSC may derive greater benefit from immunotherapy and chemotherapy, and using the signature could have significant impact of the clinical treatment of these patients.

Disclosure of Potential Conflicts of Interest

No potential conflicts of interest were disclosed.

Supplementay materials

Supplemental data for this article can be accessed on the publisher’s website.

Additional information

Funding

This study was supported by grants from the National Natural Science Foundation of China (81672640), the Grant for Key Disciplinary Project of Clinical Medicine under the Guangdong High-level University Development Program (002-18119101), the Project of Innovating and Strengthening Universities in Guangdong Province Supported by Department of Education of Guangdong Province (2018KTSCX066), and the Special Funds for Innovation Strategy of Science and Education in Guangdong Province (2018-157), the Science and Technology Planning Project of Shantou City ([2019]77), and Supporting Program of the First Affiliated Hospital of Shantou University Medical College ([2019]70).