Genomic Network-Based Analysis Reveals Pancreatic Adenocarcinoma Up-Regulating Factor-Related Prognostic Markers in Cervical Carcinoma

We previously showed that PAUF is involved in tumor development and metastases in cervical cancer. This study was conducted to discover novel molecular markers linked with PAUF in cervical cancer using genomic network analysis and to assess their prognostic value in cervical cancer. Three PAUF-related genes were identified using in-silico network-based analysis of the open genome datasets. To assess the expression of these genes and their relationship to the outcome of cervical cancer, immunohistochemical analysis was performed using cervical cancer TMA. The associations of the identified proteins with clinicopathologic characteristics and prognosis were examined. AGR2, BRD7, and POM121 were identified as interconnected with PAUF through in-silico network-based analysis. AGR2 (r = 0.213, p < 0.001) and POM121 (r = 0.135, p = 0.013) protein expression were positively correlated with PAUF. BRD7High and AGR2Low were significantly associated with favorable disease-free survival (DFS) (p = 0.009 and p < 0.001, respectively), and in combination with PAUFHigh, even more significantly favorable DFS observed (p < 0.001 for both). In multivariate analysis, AGR2High (HR = 3.16, p = 0.01) and BRD7High (HR = 0.5, p = 0.025) showed independent prognostic value for DFS. In a random survival forest (RSF) model, the combined clinical and molecular variable model predicted DFS with significantly improved power compared with that of the clinical variable model (C-index of 0.79 vs. 0.75, p < 0.001). In conclusion, AGR2 and BRD7 expression have prognostic significance in cervical cancer and provide opportunities for improved treatment options. Genomic network-based approaches using the cBioPortal may facilitate the discovery of additional biomarkers for the prognosis of cervical cancer and may provide new insights into the biology of cervical carcinogenesis.


INTRODUCTION
Cervical cancer is the second most commonly diagnosed cancer and is the third leading cause of cancer death among females in less-developed countries (1). Although screening for precancerous lesions and human papilloma virus (HPV) infection and preventive HPV vaccinations are good options to prevent cervical cancer, once invasive cancer occurs, recurrence remains a major problem despite improvements in treatment (2).
To improve survival in cervical cancer patient, proper treatment should be adopted. If tumor behavior could be reliably predicted at the initial diagnosis, tailored treatments could be implemented to improve survival. Although clinical factors such as International Federation of Gynecology and Obstetrics (FIGO) stage, lymph node metastasis, and tumor size may serve as markers for overall prognosis, they have limitations in accurately predicting survival (3). Thus, novel biomarkers including molecular markers are needed for accurate survival estimates. Recent advances in molecular biology and wellcurated, publicly available data, such as The Cancer Genome Atlas (TCGA), have made it possible to identify multiple valuable prognostic biomarkers for cancer studies.
We previously reported that pancreatic adenocarcinoma upregulated factor (PAUF) can be used as a prognostic molecular marker in patients with cervical cancer (4). PAUF expression is upregulated in cancers, particularly in glandular cells or adenocarcinoma, and cytoplasmic expression is independently associated with poor survival. However, the specific receptor involved, its downstream signaling pathways, and the molecular mechanisms of PAUF-dependent gene regulation need to be elucidated. Recently, there have been several network-based studies utilizing interaction information between genes in cancer (5,6). Considering the complex nature of the interaction between genes, a single genetic abnormality can spread along the links of the complex intracellular network to alter the functions of related genes (7)(8)(9). In addition, we hypothesized that multigene markers could have more accuracy in prognostication than a single molecular marker (10). To this end, the markers interconnected with PAUF could be identified with in silico analysis using publicly available datasets. The purpose of this study was to find novel proteins connected with PAUF in cervical cancer using in silico analysis and to investigate the clinical significance of those proteins in a well-defined cohort of cervical cancers using immunohistochemistry.

In silico Analysis to Select Candidate Genes Connected With PAUF
To identify the novel prognostic molecular markers related to PAUF protein expression, data from the Gene Expression Omnibus (GSE44001) and TCGA datasets were analyzed, and a detailed schematic is depicted in Figure 1A. For the details, we downloaded the GES44001 dataset from GEO website (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi? acc=GSE44001). The obtained data were dichotomized by the median expression of ZG16b (Entrez Gene ID: 124220, ID: ILMN_1753139), and each 150 samples were allocated as high or low PAUF expression (File S1). Subsequently, the differential expression analysis was performed with the "Analyze with GEO2R" in the GEO website with options of Benjamini & Hochberg (False discovery rate), Auto-detect, and Submitter supplied (https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc= GSE44001). The detailed instruction of GEO2R is available in https://www.ncbi.nlm.nih.gov/geo/info/geo2r.html#how_to_ use. The genes selected after the above steps are described in Table S1. An interactive analysis was conducted to explore genes related to ZG16B, the PAUF protein expression gene, using mutual exclusivity and network analysis with cBioPortal network analysis. cBioPortal for Cancer Genomics (http:// cbioportal.org) is a web resource for exploring, visualizing, and analyzing multi-dimensional cancer genomics data (11). The detailed queries for the analysis are as follows; samples were "cervical squamous cell carcinoma and endocervical adenocarcinoma, " and genomic profile was "mRNA Expression z-Scores (RNA Seq V2 RSEM)". The access date was November 2014.

Patients and Tumor Samples
As previously described, tissue microarrays (TMAs) were constructed using samples from 336 early-stage cervical cancer patients who were primarily treated with radical hysterectomy at the Department of Gynecologic Oncology, Samsung Medical Center, between 2002 and 2009 (12,13). Tissue samples were collected from patients who had signed an informed consent form and the study was approved by the Institutional Review Board at Samsung Medical Center, Seoul, Korea (IRB No: 2015-07-122). No patients had undergone previous radiation or chemotherapy. All patients underwent radical hysterectomy with pelvic lymphadenectomy. Para aortic lymph node sampling was performed on patients who were diagnosed with pelvic lymph node metastasis in frozen biopsy during operation, or patients with para aortic lymph node metastasis in preoperative radiologic findings. Patients with rare histology or limited availability of tissue block specimens were excluded from the TMA cohort. Following radical hysterectomy, adjuvant radiotherapy with or without concurrent chemotherapy was given according to the following risk factors: lymph node metastasis, parametrial involvement, positive resection margins, stromal invasion of more than half of the cervix, lymphovascular space invasion (LVSI), or a tumor larger than 4 cm in diameter. During follow-ups conducted every 3 months for 2 years then every 6 months for the next 3 years, Pap smears, tumor markers, and CT scans were performed every 3-12 months. Patients who relapsed within 3 years after adjuvant chemoradiation were classified as resistant to chemoradiation. Disease-free survival (DFS) was assessed from the date of surgery to the date of recurrence or the date of the last follow-up visit. Overall survival (OS) was measured from the date of surgery to the time of death or the date of last contact (for living patients). The data for patients who had not experienced an event as of the date of the final observation were censored.

Quantitative Evaluation of Immunohistochemical Staining
Immunohistochemical staining was quantitatively evaluated using computer-assisted image analysis software (Visiopharm, Hoersholm, Denmark) as described previously (14). In brief, the slides were scanned using a whole-slide scanner (NanoZoomer 2.0, Hamamatsu Photonics, Hamamatsu City, Japan) and imported into Visiopharm software (for Windows 7, version 4.5.1.324) using the TMA workflow. Staining intensity was categorized as 0, 1+, 2+, or 3+ according to the distribution pattern across cores. A brown staining intensity (0-negative, 1-weak, 2-moderate, or 3-strong) was obtained using the predefined algorithm and optimized settings. The overall immunohistochemical scores (histoscore) for AGR2, BRD7, and POM121 were expressed as the percentage of positive cells multiplied by their staining intensity (possible range, 0-300). For the survival analysis and hierarchical clustering, the expression values were dichotomized (High vs. Low, described in superscript) using cut-off values affording the most discriminative power (histoscore of 39 for AGR2, 46 for BRD7, and 204 for POM121). PAUF protein expression was evaluated as previously described (4) to compare with that of the three proteins identified.

Statistical Analysis
Statistical analysis was performed using R software, version 3.1.3 (R Foundation, Vienna, Austria; http://www.R-project. org). The expression levels of the proteins according to the clinicopathological characteristics were analyzed using Student's t-test or Mann-Whitney U-test. Analysis using Spearman's rho coefficient was used to assess the correlations between the expressions of the proteins. Survival distributions were estimated using the Kaplan-Meier method, and the relationship between survival and each parameter was analyzed with the log-rank test. A Cox proportional hazards model was created to identify independent predictors of survival. Statistical significance was defined as p < 0.05. Hierarchical clustering was performed to determine the clustering of samples according to PAUF and PAUF-related protein expression. In addition, the random survival forest (RSF) method was modified to include both clinical and molecular features to assess the predictive power of integrating the molecular data with clinical variables (15). A detailed description of the analysis has been previously described (12).

Candidate Genes Connected With PAUF (ZG16B)
Analysis of the GSE44001 dataset showed that 21 genes were differentially expressed according to ZG16B gene expression (dichotomized using the median cutoff). These identified genes were validated in the TCGA cervix dataset (Table S1). The cBioPortal network view of these genes showed that anterior gradient protein 2 (AGR2), bromodomain-containing protein 7 (BRD7), and nuclear envelope pore membrane protein (POM121) are interconnected with the ZG16B gene ( Figure 1B). These three proteins were selected to evaluate their prognostic value in cervical carcinoma.

Clinicopathological Characteristics of Patients
The clinicopathological characteristics are summarized in Table 1. The mean age of the patients was 48.9 ± 11.2 years. Of the 336 patients with cervical cancer, 291 (86.6%) were stage IIA or less, and 45 (13.4%) were stage IB2 or IIB. In Korea, many surgeons perform radical hysterectomy on patients who were suspected of having a stage IIB diagnosis in pre-operative radiologic findings or physical examination. In this study, we intended to include possible samples which were obtained from radical hysterectomy even if the sample sizes were small. Of the patients, 256 (76.2%) had squamous cell carcinoma and 80 (23.8%) had adenocarcinoma. Lymph node metastasis was found in 80 (23.8%) patients and parametrial involvement in 31 (9.2%) patients. One hundred and sixty (47.6%) patients were treated with adjuvant radiation with or without concurrent chemotherapy due to the presence of risk factors. Among them, 20 patients were radiation resistant.

AGR2, BRD7, and POM121 Protein Expression
AGR2 and POM121 protein expression was mainly observed in the cytoplasm, and BRD7 expression was observed in both the cytoplasm and nucleus. Representative examples of positive and negative staining are shown in Figure 2A. Subsequently, we examined the correlations between PAUF protein expression and these three proteins using subgroup analysis. Although there were no statistically significant differences, AGR2 and POM121 protein expression showed the tendency of a positive association with PAUF protein expression (r = 0.213; p < 0.001 and r = 0.135; p = 0.013, respectively; Figure 2B). Moreover, AGR2 protein expression was positively associated with POM121 protein expression as shown in Figure 2B (r = 0.330; p < 0.001). In the subgroup analysis, PAUF and POM121 protein expression was positively correlated with AGR2 protein expression (r = 0.305; p = 0.001 and r = 0.456; p < 0.001, respectively) in the radiation-sensitive group (n = 113; Figure S1A). Furthermore, in the radiation-resistant group (n = 20), BRD7 protein expression was positively associated with POM121 protein expression (r = 0.55; p = 0.012) ( Figure S1B). These findings are consistent with the TCGA dataset and suggest that these proteins are correlated in cervical cancer, especially within subgroups according to radiation therapy.
When compared with normal cervix tissues, a significant decrease in BRD7 (mean histoscore 63 vs. 266, p < 0.001) and POM121 (mean histoscore 144 vs. 228, p < 0.001) expression was detected in cancer tissues, but the expression of AGR2 showed no significant difference ( Table 1). The expression of AGR2 and BRD7 was cell type-dependent; AGR2 was more highly expressed in adeno/adenosquamous carcinoma, and BRD 7 expression was more prominent in squamous cell carcinoma (p < 0.001 and p = 0.004, respectively; Table 1). These results suggest that AGR2 and BRD7 potentially have different roles in cervical cancer according to cell type. In addition, the higher expression of AGR2 was negatively correlated with LVSI and the depth of invasion (p = 0.007 and p = 0.042, respectively), and POM121 expression was decreased (p = 0.001) in radiation-resistant cervical cancer.

Unsupervised Hierarchical Clustering of PAUF-Associated Markers
To examine whether patients can be grouped according to protein expression, a total of 336 cervical cancer cases were analyzed using hierarchical clustering with histoscore. As shown in Figure 3A, the patients were divided to three groups. Category 1 (n = 118) consisted exclusively of AGR2 High and POM 121 High expression. Category 2 (n = 194) consisted exclusively of AGR2 Low and PAUF Low expression. Category 3 (n = 24) consisted of PAUF High and POM 121 Low expression. There were significant differences between categories 1, 2, and 3 in terms of cell type, depth of invasion, and radiation therapy resistance (p < 0.001, p = 0.008, and p < 0.001, respectively; Table 2). Notably, category 3 was associated with adenocarcinoma cell type and the radiation-resistant phenotype and showed a poor overall survival (OS) tendency (p = 0.053; Figure 3B).

Improved Prognostic Power of Combined Molecular and Clinical Model
To examine whether a network-based model using expression of PAUF and related genes provided enhanced prognostic power compared to clinical prognostic factors (i.e., FIGO stage, lymph node metastasis, lymphovascular invasion, stromal depth of invasion, parametrial involvement, and resection margin), we

DISCUSSION
In the present study, we identified novel molecular markers interconnected with PAUF using genomic network-based analysis and evaluated their relationships and prognostic potential in a large cohort of cervical cancer patients. The identified proteins (AGR2, BRD7, and POM121) were differentially expressed according to cell type, clinical prognostic factors, and resistance to radiation therapy in cervical cancer patients. Furthermore, we demonstrated for the first time that low expression of AGR2, high expression of BRD7 or POM121 combined with low expression of PAUF predict delayed recurrence in cervical cancer patients. Based on these results, we suggest that a combined clinical/molecular variable model may have more prognostic potential than a model using clinical variables alone to predict disease recurrence. Based on the above findings, results from in silico analysis of publicly available genomic data were confirmed to show clinical significance with regard to protein levels of the new patient cohort in this study (11,16). The pro-oncogenic anterior gradient 2 (AGR2) protein stimulates cancer cell initiation, proliferation, invasion, and metastasis (17), and previous studies have shown that AGR2 is involved in p53 regulation and cell survival control (18). Additionally, elevated AGR2 protein expression is observed in numerous cancer types including gastric adenocarcinoma (19), breast (20), non-small cell lung (21), ovarian (22,23), esophageal, prostate, and pancreatic cancer (24), and overexpression of the AGR2 protein is correlated with poor prognosis in these carcinomas. Monoclonal antibodies against AGR2 markedly reduce tumor growth and metastasis in pancreatic ductal adenocarcinoma, but this reduction has only been observed in preclinical mouse models (25). In gynecological cancer, overexpression of the AGR2 protein was reported in high-grade ovarian carcinoma, which is related with metastasis and poor prognosis (22,23). Furthermore, AGR2 was associated with chemotherapy resistance in an analysis of survival in high-grade serous ovarian carcinoma (22). Consistent with previous studies, AGR2 expression was positively correlated with adenocarcinoma and poor prognosis in cervical cancer in this study.
It is well known that cervical cancers result from a persistent infection with high-risk types of HPV (26). HPV infects the basal cells of the squamous epithelium, causing entry of its DNA into the nucleus of the host cells. The HPV DNA replicates as the basal cells differentiate, and induce a dramatic increase in the expression of the two HPV oncoproteins (E6 and E7). Herfs et al. have suggested that cuboidal cells at the cervical squamocolumnar junction are the origin of cervical cancer and its precursors (27). They also identified the squamocolumnar junction cell-specific biomarkers keratin 7, AGR2, CD63, matrix metalloproteinase 7, and guanine deaminase and documented their expression in HPV-infected Patients with high AGR2 expression (AGR2 High ) and low BRD7 expression (BRD7 Low ) showed worse disease-free survival (log-rank test, p < 0.001 and p = 0.009, respectively) than patients with low AGR2 expression (AGR2 Low ) and high BRD7 expression (BRD7 High ). (C) Patients with low POM121 expression (POM121 Low ) showed a trend of worse disease-free survival. (D-F) Patients with combined PAUF High /AGR2 High , PAUF High /BRD7 Low , or PAUF High /POM121 Low expression showed shorter disease-free survival (log-rank test, p < 0.001, p < 0.001, and p < 0.001, respectively) than patients with the opposite combined expression profile. (G-J) Patients with combined AGR2 High /BRD7 Low , BRD7 Low /POM121 Low , or AGR2 High /BRD7 Low /POM121 Low expression showed shorter disease-free survival (log-rank test, p < 0.001, p = 0.009, and p = 0.044, respectively) than patients with combined AGR2 Low /BRD7 High , BRD7 High /POM121 High , or AGR2 Low /BRD7 High /POM121 High expression.
Frontiers in Oncology | www.frontiersin.org  cell lines, high-grade cervical intraepithelial neoplasms, and cervical carcinoma (27). Recently, Morbini et al. showed that only AGR2 and keratin 7 were associated with active HPV infection in oropharyngeal carcinomas (28). Taken together, these data support AGR2 involvement in cervical cancer tumorigenesis by promoting the proliferation of human cervical cancer cells. Bromodomain-containing protein 7 (BRD7), also known as celtix-1, is a subunit of the SWI/SNF complex specific for polybromo-associated BRG1-associated factor (PBAF), which was first identified as a tumor suppressor in nasopharyngeal carcinoma cells in 2000 (29). Published studies show that BRD7 plays an important role as a transcription co-factor in the expression of the TP53 gene (30) and interferes with the PI3K signaling pathways (31). Down-regulation of BRD7 expression has been observed in several malignancies and is correlated with poor prognosis of various cancer types, including breast (32), colorectal (33), gastric (29), and prostate cancer (34), as well as hepatocellular carcinoma (35), and osteosarcoma (36). Endometrial cancer (37) and ovarian cancer exhibit low expression of BRD7, which is correlated with downregulated beta-catenin accumulation (38). In ovarian cancer, BRD7 expression is decreased in high-grade serous carcinoma relative to normal or low-grade carcinoma (38). In our study, high expression of BRD7 was associated with improved DFS in cervical cancer patients. Cox regression analysis confirmed that BRD7 and its combination with PAUF and other PAUF-related proteins is an important prognostic indicator in cervical cancer.
The nuclear envelope pore membrane protein (POM121) gene is one of the PAX5 fusion genes and is related to the development of B-ALL leukemia through gene rearrangement (39). POM121 is one of the integral membrane components of the nuclear pore complex (NPC) in vertebrate cells, which maintains NPC structure (40). Decreased POM121 expression significantly reduces assembled NPCs on the nuclear envelope and induces clustering of NPCs, which regulate molecular transport through the nuclear envelope (often referred to as the nuclear permeability barrier) (41). In our cervical cancer samples, low expression of POM121 alone was not correlated with decreased DFS. However, in combination with high expression of PAUF or low expression of BRD7, low expression of POM121 predicted significantly reduced DFS when compared with BRD7 or PAUF alone. BRD7 was found predominantly in the nucleus and the distribution ratio of BRD7 in nucleus/cytoplasm has effect on its function of transcriptional regulation (42). Also, in our previous study, PAUF had prognostic value when it was located in cytoplasm (4). Accordingly, we deduced that POM121 may have function of transporting BRD7 and PAUF proteins between nucleus and cytoplasm. However, further research is needed to identify the accurate mechanism of POM121 in cervix carcinoma tumorigenesis.
In this study, we performed integrative analysis using the cBioPortal to identify molecules in the PAUF pathway and found that AGR2, BRD7, and POM121 are possible important molecules of the oncogenic function of PAUF. Furthermore, we demonstrated that AGR2, BRD7, and POM121 protein expression had prognostic value when combined with PAUF expression using an analysis of a large cohort of cervical cancer patients. These findings suggest the usefulness of network-based prognostic markers because the prognostic accuracy of a clinical variable model was enhanced by the addition of molecular markers (i.e., AGR2, BRD7, and POM121). With these efforts to find novel molecular makers in cervical carcinoma prognostication, patients with high risk of recurrence might have opportunities to improve their treatment outcomes with closer follow up.