{'fraction_non_alphanumeric': 0.027649769585253458, 'fraction_numerical': 0.017629407351837958, 'mean_word_length': 5.376153057738299, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 3, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Integrated Data Analysis Implicates PLAU, SERPINE1, SPP1, and MMP1 as Prognostic Factors in Head and Neck Squamous Cell Carcinoma

Aim: Head and neck squamous cell carcinoma (HNSCC) ranks seventh in incidence among the most common types of cancer worldwide. Although smoking and drinking are risk factors for HNSCC, the specific causes of, and molecular mechanisms underlying, of HNSCC have not been identified. Despite improvements in tumor treatment technology over the last 40 years, the prognosis of patients with HNSCC has not changed significantly. Therefore, to improve HNSCC treatment strategies and for early diagnosis, it is important to study its pathogenesis, and to identify prognostic markers for it. Materials and methods: We used an integrated bioinformatics approach to identify key pathogenic and prognostic genes involved in HNSCC and to reveal the potential underlying molecular mechanisms. The expression profiles of the GSE6631 and GSE107591 datasets were downloaded from the Gene Expression Omnibus (GEO) database, and the tertiary RNA-sequencing dataset of HNSCC and corresponding clinical information were downloaded from The Cancer Genome Atlas (TCGA). These three datasets were integrated to identify differentially expressed genes (DEGs), and DEGs were analyzed using bioinformatic tools, including R packages. Results: We identified 83 DEGs among these datasets. Gene Ontology analysis showed that the biological functions of the identified DEGs are primarily associated with regulating extracellular signal cascades, epidermis development, adhesion, and other tumor cytology behaviors. Kyoto Encyclopedia of Genes and Genomes pathway analysis showed that these DEGs were mainly involved in the PI3K-Akt, human papillomavirus infection, and IL-17 signaling pathways. A protein–protein interaction network was constructed to reveal the 20 most closely related genes among the DEGs. Conclusion: Four key genes in the network (PLAU, SERPINE1, SPP1, and MMP1) were demonstrated to have prognostic relevance based on their significant associations with survival in patients with HNSCC. We verified the high expression of SERRPINE1 in hypopharyngeal carcinoma cell lines by qRT-PCR and Western blotting. Clinical significance: Our findings would help elucidate the molecular mechanisms underlying the development and provide the possibility to improve the prognosis of HNSCC.


IntroductIon
Head and neck squamous cell carcinomas (HNSCCs) encompass tumors of various subsites, including those of the oral cavity, oropharynx, larynx, and hypopharynx. HNSCC is an aggressive epithelial malignancy, ranking seventh in incidence among the most common cancers worldwide, at approximately 36/100,000 men and 7/100,000 women. 1 Its incidence has increased annually, 2 with 600,000 new cases and >300,000 deaths caused by HNSCC globally each year. 3 Over 90% of HNSCCs are squamous cell carcinomas that are associated with a poor prognosis, despite diverse treatments being available. 4 The 5-year survival rate is approximately 35%, 5 and nearly 50% of patients have locoregionally advanced disease. 6 However, the causes, mechanisms, genetic markers, and underlying molecular events associated with HNSCC are not well defined. The use of alcohol and tobacco, immunosuppression, and exposure to high-risk human papilloma virus (HPV) are the predominant environmental risk factors for developing HNSCC. Nevertheless, not all smokers and drinkers develop cancer. Only 10-15% of smokers develop lung cancer, and fewer develop HNSCC, suggesting that there are other environmental or genetic factors at play. 7 Several studies have indicated that the occurrence of oral and oropharyngeal HNSCCs, but not cancer of the nasopharynx or hypopharynx, is linked to HPV infection. 8,9 Despite improvements in tumor treatment technology over the last 40 years, the prognosis of patients with HNSCC has not changed significantly. 10 The molecular pathogenesis of HNSCC is gradually being revealed through genomic, proteomic, and transcriptomic approaches. Epidermal growth factor receptor (EGFR), nuclear factor kappalight-chain enhancer of activated B cells (NF-κB), signal transducers and activators of transcription (STATs), and cyclin-dependent kinase inhibitor 2A (CDKN2A) have all been identified to play a key role in HNSCC. The roles of these proteins in clinical applications and prognosis are also being continuously studied. Bioinformatics, which combines molecular biology and information technology, has emerged as a valuable tool to study the underlying molecular mechanisms of various diseases and to identify tumor biomarkers. 11 With advances in molecular biology techniques, abundant gene expression data are publicly available for a diverse set of neoplasms. The Gene Expression Omnibus (GEO) and The Cancer Genome Atlas (TCGA), two important databases, provide information regarding the expression of coding or non-coding RNAs and clinical information submitted by research groups using different microarray platforms for analysis. 12,13 They are widely used for the discovery of tumor biomarkers and mechanisms. Here, to identify the key genes involved in HNSCC, we applied an integrated bioinformatics approach using mRNA expression data from the GEO and TCGA datasets. These findings would help reveal the molecular mechanisms underlying the development of HNSCC and associated tumor biomarkers.

Microarray Datasets
The keyword "HNSCC" was used to search the GEO DataSets database (https://www.ncbi.nlm.nih.gov/geo/), and gene expression profiles from the GSE107591 14 and GSE6631 15 datasets were obtained. The GPL6244 [HuGene-1_0-st] Affymetrix Human Gene 1.0 ST Array platform was used for the GSE107591 dataset, which includes 23 normal and 24 tumor tissues. The GPL8300 [HG_ U95Av2] Affymetrix Human Genome U95 Version 2 Array platform was used for the GSE6631 dataset, which includes 22 normal and 22 tumor tissues. Detailed information about the GEO datasets is provided in Table 1. We also downloaded an RNA-sequencing dataset of mRNA level 3 and corresponding clinical information for 544 patients from the TCGA, including detailed biological data from 500 tumor tissue samples and 44 samples of adjacent normal tissue. R v. 3.6.0 was used to process the downloaded files for conversion and filtering of unqualified data with the GEOquery, dplyr, limma, heatmap, ggplot2, and clusterProfiler packages, along with analysis tools from the STRING, Cytoscape, and Gene Expression Profiling Interactive Analysis (GEPIA) databases.

Data Processing and Screening of Differentially Expressed Genes Associated with HNSCC
We employed R language software and annotation packages to handle the downloaded platform data and series of matrix files. Specifically, GEO data were calibrated, standardized, and log 2 -transformed to obtain expression levels using the R packages limma, dplyr, and tibble. The platform IDs corresponding to the probe names were converted to standard recognized gene symbols and preserved in an R data file. TCGA data were calibrated and standardized, and differential gene expression analysis and annotation were performed using the R packages dplyr, lapply, tidyr, DESeq2, and edgeR. DEGs (genes that were differentially expressed between HNSCC and normal tissue) were identified according to the criteria |log 2 fold change (FC)| ≥ 1 and adjusted p value < 0.05. Heatmaps and volcano plots were created using the R package ggplot2. Subsequently, we used the R venn package to identify DEGs common to the three datasets.

Gene Ontology and Kyoto Encyclopedia of Genes and Genomes Pathway Enrichment Analyses of Differentially Expressed Genes
Gene ontology (GO) annotation includes three categories: biological process, cellular component, and molecular function. Therefore, GO enrichment analysis can offer a crude understanding of the biological functions, pathways, or cell localization of DEGs. GO analysis was performed using the R package clusterProfiler. 16 Kyoto encyclopedia of genes and genomes is a database for associating molecular-level information with the high-level functions and utilities of a biological system at the cell, organism, and ecosystem levels and is widely used for analyzing large-scale molecular datasets generated by genome sequencing and other high-throughput techniques. We also used the R clusterProfiler package to search the key pathways associated with DEGs via the KEGG database. A p value <0.05 was considered to indicate statistically significant enrichment of a gene with a given KEGG pathway or GO term.

Protein-protein Interaction Network Integration
A protein-protein interaction (PPI) refers to the process in which two or more proteins form complexes through noncovalent bonds; thus, a PPI network can be used to identify interactions between known and predicted proteins. We used the STRING database (http:// string-db.org) to identify PPIs for the DEGs found to be associated with HNSCC. The hub genes in the network represent candidate key or core genes with important physiological regulatory functions. We exported the results as a simple tab-separated value (TSV) table, which were inputted to Cytoscape v. 3.7.2 for network visualization analysis to select the hub genes.

Association of Patient Survival with Hub Genes
We used GEPIA (http://gepia.cancer-pku.cn/) to analyze the associations of the identified hub genes with overall survival of patients with HNSCC. The GEPIA server has been running for the last 2 years and has processed 280,000 analysis requests for 110,000 users from 42 countries. For this analysis, TCGA data were classified into two groups of high and low expression based on the median value. Genes with expression levels significantly (p value < 0.05) associated with overall survival were considered as candidate markers of HNSCC prognosis.

Western Blotting
Cells were lysed in RIPA protein extraction reagent (Sangon Biotech, Shanghai, China) with a protease inhibitor. Total protein concentration was quantified using the Bio-Rad Protein Assay Kit. Approximately 50 μg of each protein extract was separated by SDS-PAGE on a 7.5% acrylamide gel; the resultant protein bands were subsequently transferred to PVDF membranes. Immunoblotting was performed using rabbit anti-SERPINE1 antibodies (1 μg/mL; ab66705, Abcam, China) and anti-GAPDH antibodies (2 × 10 − 4 μg/ mL; ab181602, Abcam, China).

Identification of Differentially Expressed Genes
The standardized data of the HNSCC expression microarray datasets (GSE107591 and GSE6631) are summarized in Figure 1. DEGs (400), comprising 254 upregulated and 146 downregulated genes, were identified using the GSE107591 dataset, which contains 18,832 genes. Under the same criteria, 162 DEGs (64 downregulated and 98 upregulated) were identified using the GSE6631 dataset. In were identified using the TCGA dataset. The DEGs identified in the three microarray datasets are shown in Figure 2. For simplicity, gene symbols are shown only for those with |log 2 FC| > 3 (Figs 2A and B) or |log 2 FC| > 9 (Fig. 2C). Hierarchical clustering of DEGs is shown using heatmaps (Fig. 3). A Venn diagram revealed 83 intersecting DEGs from three datasets ( Fig. 4 and Supplementary Table S1).

GO Term Enrichment Analysis of DEGs
The top 10 enriched GO terms associated with the identified DEGs are presented in Figure 5. In the biological process group, the DEGs were mainly enriched in extracellular matrix organization, extracellular structure organization, epidermis development, skin development, cell-substrate adhesion, and regulation of cell-substrate adhesion, which are all processes related to the occurrence and development of HNSCC. In the cellular component group, the DEGs were mainly enriched in extracellular matrix, collagen-containing extracellular matrix, endoplasmic reticulum lumen, extracellular matrix component, basement membrane, and collagen trimer, which are mainly associated with cell proliferation. In the molecular function group, the DEGs were mainly enriched in the extracellular matrix structural constituent, endopeptidase activity, metallopeptidase activity, serine-type peptidase activity, serine hydrolase activity, and structural constituent conferring tensile strength, which are largely focused on regulating enzyme activity, cell proliferation, migration, and invasion.

KEGG Pathway Analysis of DEGs
KEGG pathway analysis demonstrated that the DEGs were mainly involved in the extracellular matrix (ECM)-receptor interaction, PI3K-Akt signaling pathway, HPV infection, interleukin (IL)-17 signaling pathway, amoebiasis, AGE-RAGE signaling pathway in diabetic complications, focal adhesion, relaxin signaling pathway, and protein digestion and absorption (Fig. 6). The specific results for the enrichment of DEGs in the ECM-receptor interaction pathway are shown in Supplementary Figure S1.

PPI Analysis of DEGs
The PPI network based on the DEGs has 89 edges, with 83 nodes and an average node degree of 2.14 ( Fig. 7).

Higher Expression of SERPINE1 in FaDu Cell Lines than in DOK Cell Lines
SERPINE1 was selected from the four key genes associated with overall survival, and was screened for cytological verification of its upregulation. We used a head and neck tumor cell line (FaDu) for the analysis; a human immortalized normal mucosal cell line (DOK) was used as a control. qPCR and Western blot results showed that SERPINE1 mRNA and protein expression levels were higher in the FaDu tumor cell line than in the DOK control cell line (Fig. 10).

dIscussIon
We used an integrated bioinformatics approach to identify key pathogenic and prognostic genes, which is an effective method to explore the underlying molecular mechanisms of HNSCC. Alterations in three types of genes-oncogenes, tumor-suppressor genes, and stability genes-are responsible for tumorigenesis. 17 Mammalian cells have several safeguards to protect against the potentially lethal effects of cancer gene mutations. Therefore, mutations in several genes are generally required for the development of invasive cancer. 17 Although all genes can undergo changes in mutation rates, only mutations in oncogenes and tumor-suppressor genes affect net cell growth and thereby confer a selective growth advantage on the mutant cells. Thus, the DEGs screened from the three databases are likely to be key factors in the cell proliferation network, and thus in the occurrence and development of HNSCC. Three studies have sought to elucidate the mechanisms underlying the occurrence of HNSCC and two of these studies focused on a single gene or single cohort. 18,19 These studies were unsuccessful in establishing clear and effective treatment options or targets for HNSCC, and the 5-year survival rate remains relatively low. Therefore, additional genes must be researched for their roles in the tumorigenesis, development, and prognosis of HNSCC. A good starting point is to methodically explore the contents of an existing database. Using three separate databases, the present study included the gene expression profiles from 635 samples, which were integrated for bioinformatic analysis, revealing 83 DEGs as candidate markers associated with HNSCC. Functional annotation showed that most of the DEGs associated with HNSCC are significantly enriched in functions related to cell proliferation and transcriptional activity.
Several matrix metalloproteinase (MMP) genes, which have been widely linked to promoting tumor invasion and metastasis when overexpressed, 20,21 were upregulated in HNSCC. MMP1 was previously shown to play an important role in initiating tumor growth and promoting HNSCC cell invasion. 22 However, another study revealed that MMP1 was not clinically useful as a diagnostic marker for HNSCC. 23 The PPI network of these DEGs identified 20 hub genes, four of which (PLAU, SERPINE1, SPP1, MMP1) were identified by GEPIA as significantly associated with survival and are candidate prognostic markers in HNSCC patients. Urokinasetype plasminogen activator (PLAU) and serpin family E member 1 (SERPINE1) were previously identified as candidate biomarkers for HNSCC using an integrated bioinformatic analysis. 24 We verified our findings for SERPINE1 by evaluating its mRNA expression using qPCR and its protein expression using Western blotting. SERPINE1 encodes a member of the serine proteinase inhibitor (serpin) superfamily, which is the principal inhibitor of tissue plasminogen activator and urokinase and hence an inhibitor of fibrinolysis. PLAU and SERPINE1 were identified as independent prognostic factors for predicting overall survival in HNSCC patients from the TCGA cohort; however, their combination had better prognostic value than either gene alone. [25][26][27] Secreted phosphoprotein 1 (SPP1) is involved in the attachment of osteoclasts to the mineralized bone matrix and binds to hydroxyapatite with high affinity. Moreover, SPP1 may also interact with the vitronectin receptor located in the cell membrane. Another study suggested that SPP1 protein and mRNA expression levels were upregulated in metastatic castration-resistant prostate cancer, and that SPP1 may regulate the androgen receptor signaling pathway. Therefore, SPP1 may be a new target for effective therapeutic approaches for the prevention of, and intervention in, drug resistance and metastasis, and also as a potential biomarker and target for the diagnosis and treatment of metastatic castrate-resistant prostate cancer. 28 Although our integrated bioinformatics approach revealed some key genes associated with HNSCC, the interpretation of their clinical value and prospects is limited without additional cytological and gene signal pathway research. Thus, further cytological verification is needed to better understand the role of these DEGs in HNSCC. However, the genes identified herein can serve as a starting point for investigating the causes and underlying mechanisms of HNSCC. These DEGs may show clinical relevance in the early diagnosis and prevention of HNSCC and may provide effective targets for treatment.

Author contrIbutIons
Pin Dong and Xinwei Chen designed the study and provided technology and software. GLW, Xinjiang Ying and Bin Shen interpreted the data. Lixiao Chen drafted the initial manuscript and performed the experiments. All authors read and approved the final manuscript. addition, we would like to acknowledge the helpful comments and suggestions on this article received from our reviewers and editors.

dAtA stAteMent
The datasets used and analyzed during the current study are available from the corresponding author on reasonable request.