Up-regulation of Mitotic Nuclear Division-related Genes and Down-regulation of TGF-β-related Genes Correlate with Survival in Malignant Pleural Mesothelioma

Malignant pleural mesothelioma (MPM) is a highly aggressive and lethal malignancy; however, its molecular origins remain largely unknown. The purpose of this study was to use transcriptomics to explore the molecular mechanisms of MPM tumorigenesis toward gaining a better understanding of disease development. The transcriptomes from MPM tissues and paired normal pleural tissues were compared to identify signicantly differentially expressed genes. These genes were then subject to Gene Ontology analysis to explore pathways that are dysregulated in MPM, and The Cancer Genome Atlas (TCGA) gene expression prole data were analyzed to assess the correlation between the signicantly dysregulated pathways and patient survival. Moreover, three independent transcriptomic datasets were used to validate the association between signicantly dysregulated genes and the clinical features of MPM. We identied 136 up-regulated and 599 down-regulated genes in MPM tissues as compared to normal pleural tissues. Functional enrichment analysis showed that the up-regulated genes were mainly associated with mitotic nuclear division, whereas the down-regulated genes were mainly associated with transforming growth factor beta (TGF-β) signaling. The 14 most signicantly up-regulated mitotic nuclear division-related genes were more likely to be highly expressed in pathological subtypes with higher malignancy, whereas the down-regulated TGF-β pathway-related gene PPARGC1A displayed the opposite trend. In the TCGA dataset, up-regulation of mitotic nuclear division-related genes was associated with a poor prognosis, whereas down-regulation of TGF-β pathway-related genes was associated with a positive prognosis. The down-regulated PPARGC1A were further validated by immunohistochemical analysis. In summary, the correlation between up-regulation of mitotic nuclear division-related genes and downregulation of the TGF-β pathway-related gene PPARGC1A with overall survival indicates an important role for these genes in MPM development and progression.


Introduction
Malignant pleural mesothelioma (MPM) is a rare, highly lethal, chest malignancy, of which there were over 30,000 new cases and 25,000 deaths worldwide in 2018. Median survival following diagnosis is generally less than 1 year [1,2]. MPM can be divided into three subtypes: epithelioid, biphasic or mixed, and sarcomatoid, of which epithelioid MPM is the most common subtype and has the best prognosis [3].
Treatment options available for MPM include surgery, chemotherapy and radiotherapy. Unfortunately, it is usually diagnosed at an advanced stage when surgery has little chance of success and has a high recurrence rate of 35-65% [4]. Therapy involving pemetrexed in combination with cisplatin has not been improved in the last 15 years [5]. In addition, currently, no strong guidance for second-line treatment options exists [6]. Encouragingly, recent ndings from phase I/II trials suggest that immune checkpoint inhibitor therapy might improve this situation, but further evidence is still needed [7,8]. Given the high malignancy of MPM and the lack of effective treatment strategies, there is still an urgent need to explore the mechanism underlying the development of MPM.
Omics studies are powerful tools for exploring the mechanisms of tumorigenesis and discovering new therapeutic targets for malignancies. Previous genomic studies have shown that numerous genes are frequently mutated in MPM patients, including BAP1, NF2, TP53, SETD2, DDX3X, ULK2, RYR2, CFAP45, SETDB1 and DDX51 [9]. Interestingly, inactivation of the CDKN2A and NF genes is common in MPM patients. Deletions of 9p21, which includes the CDKN2A gene, occur in over 70% of MPM patients [10], whereas mutations in NF2 are seen in over 40% of MPM patients [11]. Other chromosomal deletions have also been observed in MPM patients. For example, the loss of 3p21.1, which includes the BAP1 gene, has been found in over 30% of MPM patients [12]. Mutations in the TP53 and SETDB1 genes and a loss of heterozygosity have been suggested to de ne a genomic subtype of MPM, which is independent of histological subtype and has a prognostic value [13]. Gene expression analyses can also identify distinct MPM subtypes which are signi cantly correlated with patient survival. Aberrant Hippo, mTOR, histone methylation, RNA helicase and p53 signaling pathways have also been reported in MPM patients [9]. The mitotic spindle assembly checkpoint pathway and the cytoskeleton/spindle microtubule network are both signi cantly altered in MPM tumor tissue. In fact, in vitro studies have shown that targeting microtubules might be a potential treatment option for MPM patients. High levels of expression of MAD2L1 have also been seen in MPM tissues [14]. Another expression pro ling study has identi ed candidate oncogenes such as NME2, CRI1 and PDGFC in up-regulated genes in MPMs, and tumor suppressor genes such as GSN which is underexpressed in MPMs [15]. Expression of the immune checkpoint gene VISTA has also been shown to be higher in epithelioid MPM patients than in patients with other pathological subtypes [13]. Several of these overexpressed genes identi ed in MPM patients have the potential to be useful as diagnostic and prognostic biomarkers. These genes include MSLN, PDGFRB, AURKA, and BIRC5 [10,16]. However, despite these reports, there remain no established markers for MPM.
The low incidence, high mortality rate, and lack of effective treatment options make it urgent to explore the molecular mechanisms in MPM. In this study, we used publicly available independent expression pro ling datasets in MPM to explore the underlying molecular mechanisms behind MPM tumorigenesis and to mine for new prognostic factors in MPM.

Datasets
All microarray expression datasets were downloaded from the GEO database (https://www.ncbi.nlm.nih.gov/geo/). Data obtained from 38 pairs of MPM tissue and normal lung tissue from GSE51024 were assessed using the paired Student's t-test to identify genes that were signi cantly up-regulated or down-regulated in MPM tissues. Expression pro le datasets that included pathological subtype information (GSE2549, GSE29354, and GSE42977) were also used to validate the association between signi cantly abnormal pathways and the clinical features of MPM. MPM RNA sequencing data from TCGA were obtained for survival analysis using the R package "RTCGAToolbox" [17]. The characteristics of these datasets are summarized in Table I.

Functional analysis and visualization
Gene Ontology (GO) analysis was used to perform functional and pathway enrichment analysis for assessing the effect of the differentially expressed genes on cellular function. The R package "clusterPro ler" was used to perform the GO analysis and visualize the results [18]. The R package "ggplot2" was used to visualize the other data [19].

Immunohistochemistry (IHC) analysis
This study was approved by the ethics committee of Shandong Provincial Hospital (SWYX: No.2019 − 213). Para n-embedded tissue sections from 20 MPM patients and 11 tuberculous pleurisy patients were obtained for IHC analysis. The tissue sections were dewaxed and hydrated, followed by antigen retrieval in EDTA solutions. After blocking the endogenous peroxidase, sections were blocked with 3% BSA and incubated with anti-PGC1α antibody (ab54481, Abcam, USA; 1:250) overnight. HRP-labeled goat anti-rabbit IgG antibody (GB23303, Servicebio, China; 1:200) was added and incubated for 50 minutes.
Diaminobenzidine was used as the chromogen and the counterstaining was performed with hematoxylin. The sections were then dehydrated and prepared for staining intensity and positivity analysis. Semiquantitative analysis of the ratio of tissue sections to cells used pale yellow to tan as markers of positive cells. The IHC results were reviewed again by two trained pathologists. Staining intensity was scored as follows: Negative (0), weak (1), moderate (2) and strong (3). The percentage of positive cells was scored as follows: 0% (0), 1-25% (1), 26-50% (2), 51-75% (3) and 76-100% (4). Three 200-fold elds were randomly taken from each slide. The nal PGC1α IHC score was the multiplication of these two scores for the three elds. All methods were carried out in accordance with the Declaration of Helsinki principles and approved guidelines in Shandong Provincial Hospital. Written informed consent was obtained from all subjects.

Statistical analysis
All statistical analyses were performed using R software version 3.4.2 (http://www.r-project.org). The R packages "survival", "survminer" and "pheatmap" were used to perform hierarchical clustering and survival analysis. Student's t-test was used to compare relative gene expression levels between two groups. Analysis of variance (ANOVA) was used to compare relative gene expression levels among multiple histological groups. Tukey's HSD method was used to perform post hoc all-pairs comparisons. Pearson correlation analysis was used to assess the correlation among the targeted genes. All statistical tests were two-sided and a p value of < 0.05 was considered statistically signi cant [20][21][22].

Dysregulated genes and functions in MPM
The paired Student t-test at the probe level using the GSE51024 dataset identi ed 735 differentially expressed genes in MPM, including 136 up-regulated genes and 599 down-regulated genes ( Fig. 1A; Supplementary Table I). Unsupervised hierarchical clustering of these differentially expressed genes revealed that four normal samples clustered with the malignant samples (Fig. 1B). GO analysis of the upregulated genes revealed dysregulation in functions related to extracellular matrix and chromosome segregation ( Fig. 1C; Supplementary Table II), whereas the down-regulated genes were associated with cilium organization, microtubule-based movement, urogenital system development, and response to TGFβ processes ( Fig. 1D; Supplementary Table III).
Mitotic nuclear division and TGF-β-related genes are associated with survival in MPM Based on previous reports, we focused on 14 mitotic nuclear division-related genes in the up-regulated gene set (IGF1, CCNB1, NCAPG, NDC80, BIRC5, CDC20, BUB1B, NUSAP1, CENPF, KIF18B, PRC1, AURKA, ANLN, and TTK) and 19 TGF-β-related genes in the down-regulated gene set (CDKN2B, ID1, F11R, PPARGC1A, PRDM16, PRKCZ, CAV1, CAV2, CFLAR, CITED2, FOLR1, NKX2-1, EDN1, HPGD, ITGA8, TGFBR3, EPB41L5, SMAD6, and SMAD7) to evaluate their relationship with the clinical features of MPM patients. We next used an unsupervised clustering method with each gene set to separate the TCGA dataset into two groups. Neither group correlated with pathological MPM subtypes ( Fig. 2A and B). With the exception of IGF, the other 13 mitotic nuclear division genes clustered in one class, and patients with high levels of expression for the genes in this group had a worse prognosis (Fig. 2A). The 19 TGF-β pathway-related genes were also divided into two groups and were found to be related to patient survival (Fig. 2B).
Speci cally, patients with a better prognosis had relatively high levels of expression of one group of TGFβ-related genes, consisting of CDKN2B, ID1, F11R, PPARGC1A, PRDM16, and PRKCZ. Subsequent analysis of each single gene revealed that patients with relatively high levels of expression of PPARGC1A, CDKN2B, and PRKCZ had better prognoses (Fig. 2C), whereas relatively high levels of expression of all of the mitotic nuclear division-related genes (except for IGF1) were indicators of a poor prognosis ( Supplementary Fig. 1).
Correlation analysis using the mitotic nuclear division and TGF-β pathway-related genes To further explore the nature of gene dysregulation in MPM, we next examined the correlation between the mitotic nuclear division-related genes and TGF-β pathway-related genes. Using a correlation analysis of 33 genes in four different datasets, we found that the mitotic nuclear division-related genes were highly correlated with each other (Fig. 3). Despite a relatively poor correlation within the TGF-β-related gene set, PPARGC1A showed a negative correlation with the mitotic nuclear division-related genes among the four different datasets (Fig. 3).

Differential expression of the mitotic nuclear division and TGF-β pathway-related genes among different pathological subtypes
Because different pathological subtypes of MPM have varying degrees of malignancy, we assumed that the genes differentially expressed between tumor and tissues would be related to the malignant transformation potential of pleural cells, and thus might be associated with degree of malignancy.
Accordingly, we explored whether or not highly correlated genes were differentially expressed among different pathological subtypes of MPM. In the TCGA dataset, the expression levels of several mitotic nuclear division-related genes were relatively low in epithelial MPM ( Fig. 4A; Supplementary Fig. 2), whereas a high level of PPARGC1A expression was observed in epithelial MPM (Fig. 4B). This opposite trend suggests that mitotic nuclear division and TGF-β-related genes play different roles in the malignant transformation of MPM.

IHC validation of PGC1α downregulation in MPM tissues
Based on the negative correlation between the levels of PPARGC1A and mitotic nuclear division-related genes, we next evaluated the expression pattern of PGC1α, which is encoded by PPARGC1A, in the MPM tissues by IHC analysis. Compared to the tuberculous pleurisy tissues, the PGC1α staining score was signi cantly lower in MPM tissues (p = 0.023; Fig. 5).

Discussion
In this study, we used the expression pro les of paired MPM and normal tissues to identify signi cantly differentially expressed genes and dysregulated pathways in MPM patients. We found that the expression of mitotic nuclear division-related genes was up-regulated, whereas TGF-β-related genes were down-regulated in the tumor tissues compared with the paired normal lung tissues from MPM patients. Correlation analysis revealed that the expression levels of PPARGC1A, as a TGF-β-responsive gene, were negatively correlated with the levels of mitotic nuclear division genes. We also validated the down-regulated expression of the encoded product of PPARGC1A, PGC1α, in MPM tissues by immunohistochemistry. In addition, single-gene analysis revealed that high PPARGC1A expression levels are a positive prognostic predictor for MPM and that PPARGC1A tended to be overexpressed in epithelioid MPM, whereas most mitotic nuclear division-related genes such as BRIC5 and CCNB1 showed the reverse trend.
Previous expression pro ling studies have identi ed several prognosis-related genes in MPM patients. Early studies using the GSE51024 dataset found that the mitotic spindle assembly checkpoint pathway is signi cantly up-regulated in MPM tumors, and that low MAD2L1 protein expression levels were associated with a better prognosis in MPM patients [14]. In the present study, we also found that genes associated with the mitotic nuclear division process were signi cantly up-regulated in MPM tumors, and that high expression levels of this gene set, and most of the single genes contained within it, were related with a poor prognosis of MPM patients. Other studies have identi ed several events as independent adverse prognostic factors in MPM, such as AURKA overexpression and CDKN2A deletion [10]. Elevated serum and pleural uid levels of the well-known apoptosis inhibitor survivin, which is encoded by the BIRC5 gene, have also been shown to be negative prognostic factors in MPM patients [23,24]. Positive IHC expression of survivin could effectively distinguish MPM from reactive mesothelial hyperplasia with 100% accuracy [25]. This consistency of gene expression pro le studies and single-gene studies indicates that disruption of mitotic nuclear division is a major abnormal process in MPM, and may serve as a diagnostic and prognostic biomarker, or potentially a therapeutic target.
Unlike the relatively concentrated functional enrichment of overexpressed genes, the 599 down-regulated genes were enriched in 15 biological process terms. Three of these terms were related to the cellular response to TGF-β. The TGF-β signaling pathway plays important roles in carcinogenesis, and is known to have dual functions during the development of cancer, since it is involved in the early stages of tumor suppression but also promotes the later stages of cancer [26]. Here, we focused on genes related to TGFβ signaling, and found that high expression levels of three genes are associated with a positive prognosis in MPM patients. In particular, the peroxisome proliferator-activated receptor gamma coactivator 1-alpha (PPARGC1A) gene showed the tightest association with the mitotic nuclear division genes, which are wellestablished to correlate with MPM malignancy. The PPARGC1A gene encodes the protein PGC1α, which is a transcriptional coactivator that regulates the activities of key mitochondrial genes, including cAMP response element binding protein (CREB) and nuclear respiratory factor (NRF) genes. As a critical regulator of mitochondrial biogenesis, PGC1α has been widely studied in the context of metabolic diseases such as diabetes [27]. Interestingly, PGC1α plays a controversial role in cancer among different cancer types: it is suppressed in VHL-de cient clear cell renal cell carcinoma, whereas its expression suppresses tumor growth and improves the sensitivity to cytotoxic therapies [28]. However, in vitro studies have indicated that PGC1α is essential in maintaining the survival and progression of PGC1αpositive melanoma cells, since knock-down of PGC1α activated the mitochondrial apoptotic pathway and decreased cell numbers [29]. Our results show that PPARGC1A expression is signi cantly down-regulated in MPM patients and also tended to be overexpressed in epithelioid MPM in contrast to the more malignant histological phenotypes. Survival analysis also revealed that high expression levels of PPARGC1A are related to a positive prognosis for MPM. Since the expression levels of PPARGC1A were found to be negatively correlated with those of mitotic nuclear division-related genes, a disorder in mitochondrial energy metabolism may be a key mechanism in MPM carcinogenesis and may be related to the mitotic nuclear division process.
Mechanistic studies have revealed that activated TGF-β1 signaling can down-regulate the otherwise high levels of PGC1α expression in the human skeletal muscle, contributing to a failure to improve insulin sensitivity after exercise in type 2 diabetes [27]. An in vivo study of thoracic aortic aneurysms in a mouse model showed that increased TGF-β levels reduced PGC1α levels [30]. However, the detailed mechanisms by which TGF-β causes the down-regulation of PPARGC1A remain largely unknown (Fig. 6). We also examined the role of another TGF-β pathway-related gene, CDKN2B, in MPM. Previous studies have identi ed the loss of chromosome 9p21, which contains CDKN2A, as a common genetic event in MPM patients (69% in epithelial mesotheliomas and 100% in sarcomoid mesotheliomas) [12,31,32]. Murine mesothelioma models have shown the homozygous deletion of the CDKN2A/Arf locus and the adjacent CDKN2B [33]. In vitro analysis of malignant mesothelioma cell lines revealed a deletion of chromosome 9 and the absence of expression of CDKN2A RNA. However, the minimal region of overlap of the deletions was determined to be 24 kb and did not include CDKN2B, which indicates that CDKN2B expression might be regulated by other mechanisms [34]. TGF-β can regulate the transcription of CDKN2B through the canonical SMAD-dependent pathway and thus induce cell cycle arrest (Fig. 6) [35]. Overall survival analysis of mesothelioma patients after surgery according to p16/CDKN2A homozygous deletion status showed a signi cant survival advantage for the patients without the deletion [10]. Our study provides new evidence that the expression levels of CDKN2B are down-regulated in MPM patients and that a high expression level is associated with a positive prognosis.
Another gene identi ed to be potentially associated with MPM in our study was PRKCZ, which encodes the protein PKCζ that is known as a critical metabolic tumor suppressor in mouse and human cancers. PKCζ de ciency can promote the reprogramming of cancer cell metabolism to utilize glutamine for adapting to variable nutrient circumstances. Patients with lower PKCζ levels also tend to have a poor prognosis [36]. Unlike CDKN2B, PRKCZ is regulated by TGF-β through non-canonical SMAD-independent signaling. TRAF6 interacts with the activated TGF-β receptors, which induces activation of the PI3K/AKT pathway [37]. Studies on islet glucose metabolism suggested that activation of the PI3K pathway can then induce the phosphorylation of PKCζ (Fig. 6) [38]. Thus, the down-regulated genes identi ed as prognosis-related genes in MPM may be valuable for patient strati cation in clinical practice.

Conclusions
We adopted a bioinformatics approach to explore the key dysregulated pathways in MPM patients, and found that genes associated with the mitotic nuclear division process and response to TGF-β were signi cantly up-and down-regulated, respectively. Three of the TGF-β pathway-related genes were identi ed as negative prognostic factors, and one in particular, PPARGC1A, was negatively correlated with mitotic nuclear division genes and showed reduced expression levels in the less malignant phenotype.
These gene sets have potential prognostic value in MPM and might provide further information for the diagnostic and prognostic evaluation of MPM patients.

Data availability
The data included in this study are publicly available at the TCGA database and https://www.ncbi.nlm.nih.gov/geo/ website by searching the GEO accession number displayed in the Table I.  Correlation plots of the mitotic nuclear division and TGF-β-related genes among the four different datasets.