KRT19 is a Promising Prognostic Biomarker and Associates with Immune Infiltrates in Serous Ovarian Cystadenocarcinoma

Introduction

Ovarian cancer is one of the most lethal malignancies in women worldwide, and in which 90% are of an epithelial cell type and comprise multiple histologic types, and the remaining 10% are non-epithelial ovarian cancers, which include mainly germ cell tumors, sex cord-stromal tumors, as well as small cell carcinomas.¹ High-grade serous ovarian carcinoma (HGSOC) is the most common epithelial subtype with a staggeringly low five-year survival rate of 30% to 40% and is responsible for 70% to 80% of ovarian cancer-related deaths.^2,3 Even with bulk surgery and adjuvant chemotherapy, as BRCA1/2 carriers with epithelial ovarian cancer respond better than non-carriers to platinum-based chemotherapies,⁴ however, the mortality rate for late stage patients still reaches a high level.⁵ As the development of biotechnologies, such as transcriptomics, proteomics and genomics,^6,7 these new perspectives can uncover promising therapeutic choices, and contribute to explore prognostic biomarkers and therapeutic targets for serous ovarian cystadenocarcinoma.

Keratin type I cytoskeleton 19 (KRT19) is an intermediate filament protein of 40 kDa which plays an essential role in the maintenance of epithelial cells, and it was believed to be a possible specific marker for epidermal stem cells.⁸ Upregulation of KRT19 has been reported in various tumor tissues, and KRT19 is associated with tumor proliferation and metastasis, such as lung adenocarcinoma (LUAD), breast cancer (BRCA), hypopharyngeal squamous cell cancer and thyroid cancer (THCA).^9–11 However, the expression patterns of KRT19 and its clinical significance in OV remain unclear.

The tumor microenvironment plays a crucial role in cancer progression; different signals from the microenvironment contribute to promote or suppress tumors.¹² For instance, immune infiltrates of solid tumors may be a hallmark of disease progression while the tumor microenvironment may reprogram tumor-infiltrating immune cells to acquire tumor-promoting functions.¹³ Revythis et al¹⁴ also reported that PARP and immuno-checkpoint inhibitors can enhance ovarian cell apoptosis by suppressing T cell activation. Due to the heterogeneity of the microenvironment in OV cancer, clinical treatment for these patients remains challenging.¹⁵ Therefore, exploring promising prognostic markers to explain tumor immune microenvironment heterogeneity can provide significant ways for ovarian cancer therapy.

In this study, TCGA database, GTEx database and GEO database were used to evaluate the correlation between KRT19 expression and prognosis of OV patients. The function of KRT19 was explored by single cell analysis of cancerSEA, and the expression of KRT19 in different tissues was analyzed by using GEO database for further analysis. In addition, our study demonstrated the relationship between KRT19 mRNA levels and tumor-infiltrating immune cells, suggesting that KRT19 may accelerate tumor processing by regulating immune cell infiltration in OV.

Materials and Methods

Gene Expression and Clinicopathological Character Analysis

We investigated the KRT19 expression level in various types of normal tissues and tumor in the integrated datasets combining TCGA¹⁶ (https://portal.gdc.cancer.gov) and the GTEx (Genotype-Tissue Expression) database¹⁷ (https://www.gtexportal.org/ home/-index.html). We obtained the profiles of RNA expression (RNA-Seq2 level 3 data; format: TPM; platform: Illumina HiSeq 2000) and clinical sample of OV patients from the TCGA database. TCGA included 379 OV samples, which contain general information, clinicopathological details, and prognostic information; GTEx included 88 normal samples. The gene expression profiling data sets ([GEO: GSE12470], [GEO: GSE6008]) were obtained from the GEO database (https://www.ncbi.nlm.nih.gov/gds). We analyzed the cancer OMICS data by using The University of Alabama Cancer database (UALCAN) (http://ualcan.path.uab.edu/) which provides protein expression analysis data from the Clinical Proteomic Tumor Analysis Consortium.^18,19 The Human Protein Atlas²⁰ (https://www.proteinatlas.org/) provides comprehensive data on human protein transcripts derived from unique human samples, which include cells, tissues, and organs. In this study, we obtained KRT19 protein immunohistochemistry results from both ovarian cancer tissues and normal ovarian tissues.

Survival Analysis

The Kaplan–Meier survival analysis was performed using the R package survival (v 3.2–10) and survminer (v 0.4.9). Patients in TCGA OV were separated into two groups based on high and low levels of KRT19 expression. We constructed a prognostic classifier to determine if KRT19 expression level affects clinical outcomes of OV patients and assessed three outcome measures: Overall Survival (OS), Progress Free Interval (PFI), Disease Specific Survival (DSS).

Analyses of Univariate and Multivariate Cox Regression

Univariate and multivariate analyses were conducted with Cox proportional hazards regression models. The survival package (v 3.2–10) implemented in R was used. Using univariate Cox regression analysis, we compared the OS of two cohorts of patients with OV and the level of KRT19 expression. Furthermore, we used multivariate analysis to determine whether KRT19 is an independent prognostic marker for OS of OV patients. There is statistical significance for KRT19 in Cox regression analysis, when p < 0.05.

Gene Set Enrichment Analysis

KRT19 mRNA expression was analyzed by R (v 3.6.3), followed by GSEA analysis using the limma, clusterprofiler and ggplot2 package. Adjusted p < 0.05, |ES|>1, and FDR < 0.25 were considered to be statistically significant.

Co-Expression Network Establishment of KRT19

The co-expression genes of KRT19 were screened by R limma package. Spearman correlation coefficient and Z-test from TCGA database were used to examine the correlation between KRT19 expression level and co-expressed genes (related parameter settings: |Spearman correlation coefficient| >0.2 and p-value <0.001).

Functional Enrichment Analysis

GO–KEGG enrichment analysis was performed by Bioconductor package clusterProfiler to better understand KRT19 carcinogenesis. The ggplot2 package was employed for visualization. Adjusted p < 0.05 is considered to be statistically significant in the enrichment results.

Single-Cell Analysis

CancerSEA (http://biocc.hrbmu.edu.cn/CancerSEA/home.jsp) aims to comprehensively decode the functional states of cancer cells at a single-cell resolution.²¹ There was 224 cells used to explore the function of KRT19 expression. They were from the dataset “Guillaumet-Adkins A. Genome Biol. 2017 (PDX_OV)”, and the dataset contains a cell group, which denoted patient-derived xenograft Ovarian carcinoma (OV) sample.

KRT19 Biological Function Analysis

We used two microarray data set from the GEO database: GSE69428, GSE40595. GSE69428 contains 10 HGSOC samples, 10 normal oviduct samples, 2 FTSC samples, 2 FTSC_IM samples, 2 FTSC_T samples and three Xenograft_FTSC_T samples. GSE40595 contains 8 microdissected normal ovarian stroma samples, 6 ovarian surface epithelium (HOSE), 31 laser microdissected cancer associated stroma samples, and 32 epithelial tumor samples from HGSOC patients. Boxplots were constructed in R using the packages ggstatsplot (v 0.9.1). Density plots were produced with the function “geom_desity” present in the ggplot2 (v 3.3.6) package of the R statistical environment.

Immunity Analysis

The tumor-infiltrating immune cells were assessed by ssGSEA (single-sample gene set enrichment analysis). R package GSVA and immune datasets were used to determine the infiltration level of immune cells in OV. The TISIDB database (http://cis.hku.hk/TISIDB), an integrated repository tool, were performed to analyze the interaction between tumor and immune system.²² In addition, we will further explore the immune potency of KRT19 in cancer. First, correlations between KRT19 expression and Immunoinhibitors in multiple types of cancer were determined via the “Immunomodulator” module. Second, statistical assessment was performed using a Spearman correlation to screen the appropriate candidates.

Statistical Analysis

We conducted all statistical analyses by R (v 3.6.3). We generated ROC curves using the pROC package. We used the chi-square test to examine the association between KRT19 mRNA expression and clinical characteristics. The prognostic value of KRT19 mRNA expression was analyzed by multivariate Cox analysis and Kaplan–Meier analysis. Statistical significance was represented by p < 0.05 (*), p < 0.01 (**), and p < 0.001 (***).

Results

OV Patients’ Clinical Characteristics

We transformed all transcript data to Transcripts Per Million (TPM) from GTEx and TCGA databases which are available for public to allow for comparisons across TCGA and GTEx. The data set contains 379 tumor samples from TCGA, 88 normal samples from GTEx. We downloaded clinical data including patients’ FIGO stage, Primary therapy outcome, Race and Age (Table 1).

Table 1 Clinical Characteristics of OV Patients

Expression of the KRT19 in Different Tumors of OV Patients

We first analyzed KRT19 mRNA expression in various tumor and normal tissue types using the Cancer Genome Atlas (TCGA) and Genotype-Tissue Expression (GTEx) databases. As shown in Figure 1A, KRT19 was significantly highly expressed in the vast majority of tumors compared to normal tissues, such as CESC, CHOL, COAD, KIRP, LUAD, LUSC, OV, PAAD, PCPG, READ, STAD. In detail, for OV tumor, a significant upregulation of KRT19 mRNA was found in OV tissues (p < 0.001) (Figure 1B). With UALCAN database, we also assessed the protein levels of KRT19, the results showed that KRT19 protein levels were upregulated in OV cancer compared with that in normal ovarian tissues (p< 0.001) (Figure 1C). Moreover, the immunohistochemical stain of KRT19 validated the higher expression in OV (Figure 1D).

Figure 1 Expression patterns of KRT19 in tumors and normal tissues. (A) KRT19 mRNA levels in distinct tumors and normal controls using TCGA and GTEx databases. (B) KRT19 expression in OV and normal ovarian tissues based on TCGA, GTEx and GEO databases. (C) KRT19 protein levels in OV and normal ovarian tissues analyzed by the UALCAN database. (D) Representative results of immunohistochemically stained KRT19 proteins in OV and normal ovarian tissues from Human Protein Atlas. ns, no significant difference; *p < 0.05; **p < 0.01; ***p < 0.001.

KRT19 Expression is Independently Associated with a Poorer Outcome in OV Patients

The survival curve comparison indicated that the high mRNA level of KRT19 gene was associated with poor OS (p = 0.024), PFI (p = 0.01) and DSS (p = 0.023) (Figure 2A–C). Further subgroup analysis found that KRT19 expression significantly affected the PFI in OV cases of FIGO stages I and II (p = 0.033), FIGO stages III and IV (p = 0.033), PD and SD (p = 0.01), PR and CR (p = 0.001), respectively (Figure 2D–G). The univariate analysis revealed that higher KRT19 mRNA expression in stage IV, SD and CR of OV was correlated with OS (Table 2). Thus, the KRT19 gene may be an independent risk factor for OV cancer (Table 2 and Figure 3).

Table 2 Correlations Between OS and mRNA Expression of KRT19 Analyzed by Univariate and Multivariate Cox Regression

Figure 2 Diagnostic and prognostic value of KRT19 in OV. (A–C) Kaplan–Meier curves for OS, PFI and DSS in OV patients. Effects of KRT19 expression on PFI in pathologic stages of OV for FIGO stages I and II (D), FIGO stages III and IV (E), PD and SD (F), PR and CR (G). (H) ROC curve of KRT19 mRNA expression in OV cohort. Abbreviations: OS, overall survival; PFI, progress free interval; DSS, disease-specific survival; FIGO stage, International Federation of Gynecology and Obstetrics stage; PD, progressive disease; SD, stable disease; PR, partial remission; CR, complete remission; HR, hazard ratio.

Figure 3 Univariate and multivariate regression analysis of KRT19 and clinicopathologic factors with OS in OV patients. Abbreviations: FIGO stage, International Federation of Gynecology and Obstetrics stage; PD, progressive disease; SD, stable disease; PR, partial remission; CR, complete remission; HR, hazard ratio; CI, confidence interval.

KRT19 Expression and Diagnostic Value in OV

The area under the ROC curve (AUC) (Value data: 0.5–1.0) is usually used to evaluate accuracy of detection method or hypothetical model. The AUC data is closer to 1.0, the higher the authenticity of the analytic system is. Our result indicated that the AUC data of KRT19 was 0.976, suggesting high accuracy of our model’s predictions (Figure 2H).

GSEA Identifies KRT19-Related Signaling Pathways

To identify the critical pathways associated with OV progression, we performed a GSEA analysis, and the results showed that the marked KRT19-related signaling pathways were enriched in Apical junction, Inflammatory response, TNFA signaling via NFKB, and epithelial mesenchymal transition (Figure 4).

Figure 4 Enrichment hallmarks from the gene set enrichment analysis (GSEA). (A) Apical junction. (B) Inflammatory response. (C) TNFA signaling via NFKB. (D) Epithelial mesenchymal transition. Abbreviations: NES, normalize enrichment score; FDR, false discovery rate adjustment.

Network Establishment for KRT19-Correlated Genes in OV

To further study the gene networks associated with KRT19 in OV, top 10 genes were positively correlated with KRT19, such as KRT19P1, KRT15, KRT7, S100A14, ELF3, PERP (Figure 5A, Figure S1), and top 10 genes were negatively correlated with KRT19, such as LINC01276, AC060809.1, RFX5-AS1, TIGD1, CA14, NRAP in OV cancer (Figure 5B, Figure S2). Next, GO enrichment and KEGG pathway analysis were performed based on these KRT19-correlated gene network. GO–KEGG enrichment analysis showed that the KRT19 associated genes were mainly involved in the epidermis development, skin development, cell–cell junction, intermediate filament cytoskeleton, intermediate filament, cornification and other signaling pathways (Figure 5C).

Figure 5 Gene expression and enrichment of KRT19 associated genes in OV. (A) Heatmap of the top 50 genes positively correlated with KRT19 in OV. (B) Heatmap of the top 50 genes negatively correlated with KRT19 in OV. (C) Enrichment of biological processes with KRT19 co-expression genes in OV by GO–KEGG analysis.

Biological Functions of the KRT19 in OV

To investigate the biological functions of KRT19 in OV, single-cell analysis was performed using CancerSEA. The results indicated that KRT19 positively regulated the metastasis of OV cells, while negatively regulated the invasion of OV cells (Figure 6A and B). The high expressions of KRT19 gene in the OV cells was also observed by the t-SNE method (Figure 6C).

Figure 6 Biological processes regulated by the KRT19 gene in OV from CancerSEA. (A) The positive or negative biological function associated with KRT19 in OV by single cell analysis. (B) The association of KRT19 gene with metastasis and invasion of OV cells. (C) Expression distribution with t-SNE plot of KRT19 in OV. ***p < 0.001.

The serous subtype of ovarian adenocarcinoma has been hypothesized to arise from at least two possible classes of progenitor cells: the ovarian surface epithelia (OSE) and/or an as yet undefined class of progenitor cells residing in the distal end of the fallopian tube.²³ In our findings, the expression of KRT19 was very similar in HGSOC, the fallopian tube stem cells (includes FTSC, immortalized FTSC, transformed FTSC, xenografting tumors) and normal oviduct cells calculated by using GSE69428 (Figure 7A). Moreover, for the analysis using another dataset, GSE40595, KRT19 expression was prominent down-regulated in OV cancer stroma and epithelial component of HGSOC when compared with the normal OV surface epithelium (Figure 7C). To examine the gene expression in different tissues, we compared the expression of all genes between the different groups in each GEO Series. The geom density plot for GSE69428 showed that the gene expression of FTSC and HGSOC was similar to each other (Figure 7B), and for GSE40595, the gene expression of human ovarian surface epithelium was approximately close to ovarian cancer stroma (Figure 7D).

Figure 7 KRT19 gene expression in various ovarian cells. (A) KRT19 expression in different groups of GSE69428; (B) Density plot in different groups of GSE69428; (C) KRT19 expression in different groups of GSE40595; (D) Density plot in different groups of GSE40595.

KRT19 Expression Correlates with Immune Infiltration Levels

To investigate whether KRT19 expression associate with tumor immune microenvironment, we analyzed the relationship between KRT19 expression and immune infiltration of OV cancer cells. The results showed that the immune infiltration levels of the most immune cells were positively correlated with the KRT19 expression in OV (Figure 8A), such as Neutrophils (r = 0.237, p < 0.001), NKCD56dim (r = 0.150, p = 0.003), Mast cells (r = 0.148, p = 0.004), DC (r = 0.142, p = 0.006), Th17 cells (r = 0.137, p = 0.007), iDC (r = 0.132, p = 0.010), Macrophages (r = 0.127, p = 0.014), Th1 cells (r = 0.126, p = 0.014), Tem cells (r = 0.124, p = 0.016) (Figure 8B–J).

Figure 8 Association analysis of KRT19 gene expression with immune infiltration. (A) Correlation between KRT19 expression and 24 tumor-infiltrating lymphocytes. (B–G) Correlations of KRT19 level with immune infiltration level of (B) Neutrophils, (C) NK CD56dim cells, (D) Mast cells, (E) DC cells, (F) Th17 cells, (G) iDC cells, (H) Macrophages, (I) Th1 cells, (J) Tem cells.

With advances in immunotherapy, Immune checkpoint inhibitors (ICIs) have been used in some novel strategies which have yielded striking effects against various cancers.²⁴ We assessed whether KRT19 expression was associated with immunoinhibitors by the TISIDB database. Each Spearman correlation between KRT19 and a distinct immunoinhibitor in cancer was integrated into the heatmap (Figure 9A). Specifically, in OV, multiple types of immunoinhibitors especially VTCN1 (r = 0.416, p = 1.59e−14), IL10RB, LGALS9, PVRL2, KDR, IL10 were significantly correlated with KRT19 expression (Figure S3). Meanwhile, KRT19 was positively correlated with the expression of several immunostimulators (Figure 9B), such as C10orf54 (r = 0.373, p = 2.04e−11), TNFSF9 (r = 0.297, p = 1.35e−07) (Figure S3). These results suggest that KRT19 may play a role in regulating tumor immunity.

Figure 9 Correlation analysis between KRT19 expression and immunoinhibitors and immunostimulators in pan-cancers. (A) Correlation heatmap between KRT19 and immunoinhibitors in tumors. (B) Correlation heatmap between KRT19 and immunostimulators in tumors.

Chemokines and chemokine receptors mediate the entry of immune cells into the tumor microenvironment.²⁵ We used TISIDB to analyze the correlation between the KRT19 expression and immune cell chemokines and chemokine receptors in OV (Figure 10). The results showed that KRT19 expression was positively correlated with several types of chemokines and receptors especially CXCL17 (r = 0.311, p = 2.95e−08), CX3CL1 (r = 0.266, p = 2.4e−06) and CX3CR1 (r = 0.156, p < 0.00616) (Figure S3).

Figure 10 Correlation analysis between KRT19 expression and chemokines and receptors in pan-cancers. (A) Correlation heatmap between KRT19 and chemokines in tumors. (B) Correlation heatmap between KRT19 and receptors in tumors.

Discussion

In this study, we conducted a systematic bioinformation analysis to investigate its potential regulatory mechanisms and biological functions in OV. According to our findings, KRT19 was significantly highly expressed in OV cancer tissues and may play a vital role in OV development. Our results are consistent with previous research regarding the KRT19 role in other cancers. The AUC of the ROC curve indicated that KRT19 expression had certain diagnostic value in OV. The survival analysis revealed that high KRT19 expression had relatively shorter survival times, revealing that KRT19 is associated with the poor prognosis of OV cancer. Our results strongly indicated that KRT19 could be used as a prognostic and diagnostic biomarker for OV.

KRT19, as a type I CK, has been reported to serve as a marker for metastatic tumors.¹¹ Traditionally, the surface epithelium of the ovary (OSE) has long been considered as the origin of HGSOC. Since the late 1990s, oviduct epithelium has emerged as a potential primary source of HGSOC.²³ However, some ovarian cancers may arise in the fallopian tube which may detract from compelling evidence that ovarian cancer arised from OSE.²⁶ By GSEA analysis, we revealed that EMT is one of the most significant causality of cancer, and is also considered as a principal biological process affected by gene expression in OV. Decreased expression of keratins during EMT contributed to an initiation of metastasis by loosening cell–cell attachment through disassembly of desmosomes,²⁷ which was consistent with our study that KRT19 was down-regulated in epithelial component of HGSOC patients. Furthermore, we found that KRT19 positively regulates the metastasis of OV cancer cells, and it was similarly expressed both in HGSOC and oviduct epithelium. We speculated that HGSOC may originated from both OSE and oviduct epithelium; the former is related to EMT which occurs due to the low-expression of KRT19 while the latter is related to metastasis which positively regulated by KRT19. In the future, we plan to further validate this model in animal experiment and clinical samples.

Furthermore, KRT19 expression is correlated with diverse immune infiltration levels in OV. According to our findings, KRT19 expression was positively correlated with most immune cells such as Neutrophils, NK CD56dim cells, Mast cells, DC cells, Th17 cells, iDC cells, Macrophages, Th1 cells, Tem cells. Interestingly, even though KRT19 expression is positively correlated with multiple immune cells, OV cancer patients with high KRT19 expression still have a poor prognosis, which aroused an important academic problem needed to be further studied.

We found that KRT19 was positively correlated with several types of immunoinhibitors and immunostimulators, such as VTCN1, C10orf54 and TNFSF9. By reducing immune response, V-set domain containing T cell activation inhibitor 1 (VTCN1) could promote cancer progression indirectly. Meanwhile, a previous study showed that low expression of VTCN1 inhibited the viability and metastasis of OV cancer by reducing the expression of CDK2/4 and CDC25C and phosphorylation of JAK2/STAT.²⁸ V-domain immunoglobulin suppressor of T cell activation (VISTA also named C10orf54) is an immune checkpoint receptor expressed on tumor-infiltrating T lymphocytes (TILs) and myeloid cells.²⁹ Mulati et al reported that C10orf54 expression in tumor cells inhibits T cell proliferation and cytokine production, leading to immune evasion.³⁰ Most cancers evade the immune response by overexpression of inhibitory ligands, which inhibit T cell function and then lead to cancer progression.³¹ Hence, the inhibition of T cells by VTCN1 and C10orf54 may provide possible mechanisms which explains why KRT19 expression correlates with immune infiltration and poor prognosis.

Otherwise, the immune infiltration is tightly regulated by the various chemokines, which modulate tumor immunity and the biological phenotype.³² Our study shows that KRT19 expression was positively correlated with several types of chemokines and receptors, especially CXCL17, CX3CL1 and CX3CR1. CXCL17 is a potent chemoattractant for monocytes, dendritic cells and macrophages.³³ Furthermore, chemokine CX3CL1 could stimulate the production of CXCL17 in keratinocyte and monocyte-derived macrophages.^34,35 CX3CR1 is a leukocyte surface protein that is most abundant on monocytes and macrophages in various organs.³⁶ All of the above may explain how KRT19 regulates some of the immune infiltration cells in OV. Overall, we speculated that the negative association of high KRT19 expression with poor prognosis is caused by inhibiting the function of immune cell rather than suppressing the migration of immune cells into the TME.

This study presented significant findings regarding KRT19 in OV, but some limitations still exist. First, despite the utilization of retrospective data from multiple databases, prospective studies are needed to determine the clinical significance of KRT19 in OV. Additionally, in order to fully understand the biological function and the potential molecular mechanisms of KRT19 in OV, further fundamental in vitro and in vivo experiments are required.

In summary, our study found that the expression of KRT19 is significantly upregulated in OV tissues and is closely related to poor prognosis of OV patients. KRT19 expression in different tissues may explain the origin of OV cancer especially HGSOC. Furthermore, given that the biological function and regulatory mechanism in immune cell infiltration, KRT19 could be applied as a novel promising prognosis biomarker for OV cancer therapy.