KCNJ3 is a new independent prognostic marker for estrogen receptor positive breast cancer patients

Numerous studies showed abnormal expression of ion channels in different cancer types. Amongst these, the potassium channel gene KCNJ3 (encoding for GIRK1 proteins) has been reported to be upregulated in tumors of patients with breast cancer and to correlate with positive lymph node status. We aimed to study KCNJ3 levels in different breast cancer subtypes using gene expression data from the TCGA, to validate our findings using RNA in situ hybridization in a validation cohort (GEO ID GSE17705), and to study the prognostic value of KCNJ3 using survival analysis. In a total of > 1000 breast cancer patients of two independent data sets we showed a) that KCNJ3 expression is upregulated in tumor tissue compared to corresponding normal tissue (p < 0.001), b) that KCNJ3 expression is associated with estrogen receptor (ER) positive tumors (p < 0.001), but that KCNJ3 expression is variable within this group, and c) that ER positive patients with high KCNJ3 levels have worse overall (p < 0.05) and disease free survival probabilities (p < 0.01), whereby KCNJ3 is an independent prognostic factor (p <0.05). In conclusion, our data suggest that patients with ER positive breast cancer might be stratified into high risk and low risk groups based on the KCNJ3 levels in the tumor.


INTRODUCTION
Human G-protein activated inward rectifier potassium channel subunits (GIRKs) are encoded by four genes (KCNJ3; KCNJ5; KCNJ6; KCNJ9). GIRK1-4 proteins form homo-or hetero-tetrameric ion channels, function as G-protein effectors in the plasma membrane and thereby regulate cellular excitability and activity via neurotransmitters and hormones [1]. GIRKs are involved in the regulation of functions as diverse as heartbeat, reward mechanisms, learning and memory functions, insulin secretion, blood platelet aggregation and lipid metabolism [1][2][3][4][5][6]. Increasing evidence suggests an involvement of genes encoding for GIRKs in tumorigenesis and tumor growth. Benign adenomas of adrenal cortex cells, which lead to aldosteronism and severe hypertension, Research Paper www.impactjournals.com/oncotarget have been linked to somatic mutations in the KCNJ9 gene encoding GIRK4 [7,8]. Upregulation of KCNJ3 gene products (i.e. GIRK1 mRNA and protein; synonyms: KGA, Kir3.1) was reported for non-small cell lung cancer [9] and pancreatic adenocarcinomas [10]. Correlation of increased KCNJ3 expression levels and breast cancer progression has been shown by several studies: Stringer et al. [11] reported increased levels of KCNJ3 mRNA in primary invasive breast carcinomas when compared to corresponding normal breast tissue and found a positive correlation between KCNJ3 mRNA expression levels in the tumor and the number of metastatic lymph nodes. Brevet et al. [12] confirmed on protein level that GIRK1 expression is higher in breast tumors than in normal breast tissue. Functional roles of KCNJ3 expression in breast cancer were investigated by Rezania et al. [13], who demonstrated that stable overexpression of KCNJ3 in MCF-7 breast cancer cells results in increased motility, invasiveness and angiogenesis compared to controls. Based on these results, we intended to study and validate KCNJ3 expression in invasive breast carcinoma samples as potential new prognostic biomarker. Consequently, the aim of the current study was a) to compare KCNJ3 expression levels between breast tumors and surrounding normal breast tissue, b) to screen the large patient cohort of The Cancer Genome Atlas (TCGA) for differential expression of KCNJ3 in clinically relevant subsets of breast cancer patients, c) to perform overall and disease free survival analysis to retrieve any possible prognostic value of KCNJ3 for breast cancer patients, d) to validate TCGA data by RNA in situ hybridization on formalin-fixed, paraffin-embedded (FFPE) breast cancer tissue samples of a previously characterized cohort (GEO ID GSE17705, [14]), and e) to get insight into the effects of KCNJ3 upregulation by performing mammosphere formation assays with MCF-7 breast cancer cell lines overexpressing KCNJ3. Our results suggest that KCNJ3 upregulation is an independent prognostic factor for estrogen receptor positive breast cancer.

KCNJ3 expression is upregulated in breast tumors compared to normal breast tissue
First, we investigated whether KCNJ3 mRNA expression is higher in breast tumors when compared to corresponding normal tissue. Analysis of TCGA gene expression data of 105 tumor samples with corresponding normal breast samples showed significantly higher KCNJ3 mRNA levels in the tumors when compared to normal breast tissue (median 14.6 vs. 6.6 normalized counts; p < 0.001; Figure 1A).
Upregulation of mRNAs in tumor cells might be caused by different mechanisms, including gene locus amplification on DNA level [15][16][17]. To explore gene amplification as a possible cause for the observed increase in KCNJ3 mRNA expression levels in breast tumors, we studied the gene copy numbers of 890 TCGA patient samples. Figure 1B shows that only two patients had a gene locus amplification (KCNJ3 copy number > 2) that resulted in increased mRNA expression. Overall, and despite a broad range of KCNJ3 mRNA expression levels, there was no significant increase or decrease in gene copy number (Spearman rank correlation coefficient r S : −0.02; p = 0.536), suggesting that high mRNA levels are generated by other mechanisms than copy number variations.
In order to validate KCNJ3 expression levels in tumor and normal cells, RNA in situ hybridization (ISH) was performed on breast cancer tissue samples to locate KCNJ3 expression in cancerous as well as in surrounding tumor stroma and normal breast epithelial cells. The results proved that KCNJ3 expression is present in tumor cells, but not in non-neoplastic cells including normal mammary ducts ( Figure 1C and 1D).

KCNJ3 expression is associated with estrogen receptor positive tumors
Next, we studied KCNJ3 mRNA expression levels of the TCGA data set to detect possible associations with clinically relevant breast cancer subsets. We observed increased KCNJ3 expression in tumors of patients who presented with positive lymph nodes when compared to those without lymph node metastasis at diagnosis (median 6.9 vs. 25.3 normalized counts; p < 0.001; Figure 2A). KCNJ3 mRNA levels did not differ significantly when patients were grouped based on their tumor size, tumor grade, age, menopausal status, histological subtype or Her2 expression status ( Figure 2B-2F). Then, patients were categorized according to their PAM50 classification (luminal A, luminal B, Her2-enriched, basal-like, normallike): patients of both the luminal A and B subtype had significantly higher KCNJ3 expression levels than patients of the basal or Her2-enriched type (p < 0.001; Figure  3A). No relevant information could be retrieved from the normal-like subtype due to the low patient number (n = 7). Since the luminal A and B subtypes include the hormone receptor positive tumors, we analyzed the patients regarding their estrogen and progesterone receptor (ER and PR) status. The two ER positive groups (ER+/ PR+ and ER+/PR−) had significantly higher KCNJ3 levels than those that were negative for ER and PR (p < 0.001; Figure 3B). Hence we conclude that high KCNJ3 expression levels are associated with positive ER, but not PR status. Finally, comparison of patients grouped solely based on ER status revealed substantially higher KCNJ3 levels in ER positive patients (median 0.4 vs. 48.9 normalized counts; p < 0.001; Figure 3C).
Based on these findings, we investigated whether KCNJ3 expression would correlate with the gene expression levels of the three known estrogen receptors ESR1 (ERα,), ESR2 (ERβ) and GPER (G-protein coupled ER 1; synonym: GPR30). Spearman rank correlation analysis showed that KCNJ3 expression correlates with ESR1 (r S : 0.521; p < 0.001; Figure 3D) but not with the other two estrogen receptors ESR2 (r S : −0.191) and GPER (r S : 0.253). In light of this finding, it was of interest to examine the correlation between KCNJ3 and ESR1 expression in our validation data set (GEO ID GSE17705). A positive correlation between KCNJ3 and ESR1 expression was observed (r S : 0.351; p < 0.001; Figure 3E), but not between KCNJ3 and ESR2 (r S : 0.165) or GPER (r S: 0.09).
Furthermore, we performed hierarchical cluster analysis with KCNJ3, ESR1, ESR2 and GPER ( Figure 4). Based on the expression levels of those four genes, patient samples clustered in three major groups designated A, B and C in Figure 4. Group A, which showed low KCNJ3 expression levels, comprised patients of the basal and the Her2-enriched PAM50 type, being generally negative for pathological ER and PR status and displaying low ESR1 expression levels. Groups B and C were mainly composed of patients of the luminal A and B types, being positive for ER and PR status and displaying high ESR1 expression. The major difference between group B and C was the KCNJ3 expression (low in group B and high in group C; Figure 4). Further, it became evident that, again, KCNJ3 expression clustered with ESR1, but not with ESR2 and GPER.

KCNJ3 is an independent prognostic factor in ER positive patients
Based on the results demonstrated above, we further focused on ER positive patients (groups B and C in Figure 4). Overall survival analysis of ER positive patients of the TCGA data set revealed that those with high KCNJ3 levels in the tumor (group C) had shorter overall survival times than patients of group B with low KCNJ3 levels (n = 647; p < 0.05; HR = 1.77 (1.04-3.02); Figure 5A). In addition, we used a Cox-proportional hazard approach for univariate and multivariate survival analysis of estrogen receptor positive patients ( Table 1). The univariate Cox model showed that age at diagnosis as well as lymph node, metastasis, PAM50 and KCNJ3 status had a significant influence on survival times. Importantly, the multivariate Cox model − with tumor size, lymph node status, metastasis status, histology, Her2 status, menopause status, age at diagnosis, PAM50 classification and KCNJ3 status as co-variates − showed that only the KCNJ3 expression status had a significant impact on survival time (n = 228; p < 0.05; HR = 5.2 (1.3-21.8); Table 1).
To validate our findings and to extend our analysis beyond the TCGA data set, we performed KCNJ3 RNA in situ hybridization (ISH) on 66 breast cancer patient samples that were available with clinical data and followup times from the estrogen receptor positive cohort GSE17705. Patient characteristics of this cohort compared to the TCGA cohort are shown in Table 2. RNA ISH was the method of choice because we previously showed superiority of ISH over classical immunohistochemistry regarding sensitivity and specificity [18]. Twelve of the samples analyzed did not meet the quality control criteria for ISH (see methods section) and were therefore excluded from further analysis (remaining n = 54). Representative images of patient samples with either low KCNJ3 levels and long survival times or high KCNJ3 levels and short survival times are given in Figure 5B and 5C. Kaplan Meier analyses showed significantly shorter overall survival (p < 0.05; HR = 2.39 (1.19-4.82); Figure 5D) and disease free survival (p < 0.01; HR = 3.1 (1.26-7.63); Figure 5E) probabilities for patients with high KCNJ3 expression.

MCF-7 cells overexpressing KCNJ3 display higher self-renewal capacity
In order to gain additional insight into the cellular mechanisms that would lead to worse patient outcomes due to high levels of KCNJ3 in tumor cells, we performed a mammosphere formation assay with the ER positive breast cancer cell line MCF-7. This was of particular interest, as pathological upregulation of different     potassium channels has been shown to affect cancer stem cell properties [19][20][21], and sphere formation is linked to self-renewal capacity and stemness of cancer cells. MCF-7 cells stably overexpressing KCNJ3 (MCF-7 GIRK1a ; [13]) formed significantly higher numbers of mammospheres compared to controls (p < 0.01; Figure 6), indicating that KCNJ3 upregulation might be involved in conferring selfrenewal capacity to cancer cells and thus contributing to higher tumor aggressiveness.

DISCUSSION
We have analyzed > 1000 breast cancer patient samples of two independent data sets (TCGA: n = 905, GSE17705: n = 298) regarding their KCNJ3 expression in order to evaluate a potential prognostic role of this ion channel gene for breast cancer.
Generally, ion channels have gained increased attention as players in cancer development and metastasis, since aberrant expression of as well as mutations in several genes encoding ion channels have been found to influence the hallmarks of cancer towards higher malignancy [22][23][24][25][26]. Potassium channels have been most comprehensively studied regarding their oncogenic potential by promoting proliferation and apoptosis [27][28][29][30][31][32]. Several studies suggested evidence for a role of GIRK1, the G-protein coupled inward rectifier K+ channel encoded by KCNJ3, in breast cancer [10][11][12]. Our results underscore these seminal findings using a substantially larger number of patient samples: we observed a significant increase of KCNJ3 mRNA expression levels in tumors of lymph node positive patients when compared to lymph node negative ones in our evaluation of 867 patients, and KCNJ3 mRNA levels were significantly higher in breast cancer samples than in corresponding normal breast tissue in a set of 105 patient samples.
The positive correlation between KCNJ3 levels and the ER status of breast cancer samples shown here had not been detected previously by the studies of Brevet et al. and Stringer et al. [11,12], supposedly because of the sizes of their patient cohorts (n = 56 and n = 31, respectively). In contrast, Ko et al. [33] performed ion channel profiling in breast cancer and mentioned a decrease of KCNJ3 expression in p53 mutant breast tumors, which are more likely to be ER negative [34], and an increase in KCNJ3 mRNA expression in ER positive tumors when compared to ER negative samples. However, no further details or survival data regarding KCNJ3 expression were given or discussed, as this study aimed to develop an ion channel gene signature (termed IC30) as prognostic tool in breast cancer, and KCNJ3 is not comprised in this final gene panel [33]. In line with this, our data clearly demonstrate that KCNJ3 expression is associated with ER positive breast cancer and that KCNJ3 levels correlate with ESR1 mRNA expression levels but not with expression levels of other estrogen receptors. However, KCNJ3 levels were variable within the ER positive patient cohorts. Most importantly, we could show in two independent data sets and by two different methods (gene expression data and RNA in situ hybridization) that ER positive patients with high KCNJ3 expression levels had worse overall and disease free survival probabilities than ER positive patients with low KCNJ3 levels. This is further corroborated by multivariate Cox proportional hazard analysis, showing that KCNJ3 is an independent prognostic marker for ER positive breast cancer, being also independent from the PAM50 subtype of the patient. Therefore, ER positive patients might be stratified into high risk and low risk groups based on the KCNJ3 levels in the tumor.
Here, we used a highly sensitive and specific RNA in situ hybridization technique for the validation of our findings, that showed positive signals in tumor cells but not in peritumoral tissue, and that allowed to clearly discriminate between patients with low and high KCNJ3 expression. It might be argued that immunohistochemistry would be the method of choice when investigating novel tumor markers. However, and several tested anti-GIRK1 antibodies did not meet the high quality standards regarding sensitivity and specificity that are required for optimal results in immunohistochemistry. Based on a comparison of the two methods, we came to the conclusion that RNA ISH is the superior technique for studying KCNJ3 expression in tissue [18].
To date, the mechanisms leading to KCNJ3 upregulation in breast carcinomas are not understood. There is no evidence for estrogen response elements (ERE) in the promoter region of KCNJ3 [35,36] making a direct activation of KCNJ3 gene transcription via ER unlikely to occur. Our results also argue against amplification of the KCNJ3 gene locus as the underlying mechanism for KCNJ3 upregulation. Further, it is unknown how the GIRK1 protein might interact with the estrogen receptor. GIRK1 is a non-functional channel when expressed as a homomer [37], located intracellularly as it accumulates in the endoplasmatic reticulum upon overexpression [38]. Thus, it remains to be determined whether and how GIRK1 and ER act together in same signaling pathways and further studies are needed to elucidate the mechanism of action of KCNJ3 upregulation in breast cancer. On a more functional level, our research group provided evidence on downstream effects of KCNJ3 upregulation in breast cancer: Rezania et al. showed that MCF-7 cells display higher wound healing capacity, increased invasion towards chemoattractants and higher motility and velocity than controls upon stable overexpression of KCNJ3 [13]. We could further corroborate these findings by showing that MCF-7 cells stably overexpressing KCNJ3 display significantly higher self-renewal capacity. Taken together, we propose a role of KCNJ3 in conferring tumor aggressiveness via invasion, metastasis and increased self-renewal capacity.
In conclusion, the present study a) confirms in two independent data sets, that KCNJ3 is upregulated in breast carcinomas when compared to normal breast tissue, b) shows that increased KCNJ3 expression is significantly associated with estrogen receptor positive breast cancer subtypes, c) highlights that increased KCNJ3 is an independent prognostic marker conferring worse overall and disease free survival probabilities to estrogen receptor positive breast cancer patients, and d) demonstrates that KCNJ3 upregulation might be involved in conferring higher self-renewal capacity to cancer cells. Future studies are required to elucidate the mechanisms that lead to KCNJ3 upregulation in ER positive breast cancer, to unveil its involvement in invasion and metastasis and to evaluate its potential as drugable target.

Gene expression data from the cancer genome atlas
The gene expression levels (RNAseq V2 level 3 data) of 950 invasive breast carcinoma samples and of 105 samples from corresponding healthy tissues were downloaded from the The Cancer Genome Atlas (TCGA) data portal (https://tcga-data.nci.nih.gov) and the upper quartile normalized counts from the RSEM pipeline were used. The corresponding clinical data of each patient were downloaded from the University of California Santa Cruz (UCSC) Cancer Genomics Browser (https://genomecancer.ucsc.edu). Male patients (n = 9) and patients without gene expression data for KCNJ3 (n = 36) were excluded for further analysis (remaining n = 905). Patient characteristics are summarized in Table 2.

Gene expression data from GEO ID GSE17705
The normalized gene expression levels of the estrogen receptor positive cohort GSE17705 (n = 298) were downloaded from the Gene Expression Omnibus (GEO). The data were generated and processed using the Affymetrix Human Genome U133A Array platform as described [14]. Briefly, tissue samples were processed by two different laboratories (MD Anderson Cancer Center, Texas, USA and Jules Bordet Institute, Brussels, Belgium; inter-laboratory reliability was assessed) and probe-level intensities were generated with Microarray Suite (MAS) version 5.0, normalized and log2 transformed.

FFPE patient samples
Detailed corresponding clinical and follow-up data as well as formalin-fixed, paraffin-embedded (FFPE) tissue samples were available for 66 patients of the ER positive cohort GSE17705. FFPE samples were collected from the Medical University of Graz, Austria and the Medical University of Innsbruck, Austria (see Table 2 for patient characteristics) [14]. Differences between the two data sets are given by a) the estrogen receptor status (all positive in the validation cohort, mixed in the TCGA)

RNA in-situ hybridization
FFPE tissue sections (thickness 4 µm) were mounted on Superfrost Plus coated slides (Thermo Scientific, Waltham, MA, USA; #10135642) and peritumoral tissue devoid of cancer was trimmed. Slides were treated according to manufacturer's instructions for the RNAscope ® 2.0 High Definition -BROWN kit (ACD, Hayward, CA, USA; #310035). Briefly, the slides were incubated at 60 °C for 1 hour followed by pretreatment (pretreat 1: 10 min, pretreat 2: 15 min, pretreat 3: 30 min). Three sections of each sample were stained with different probes: the KCNJ3 probe (#Hs-KCNJ3-tv1tv2), the negative control probe bacterial dihydrodipicolinate reductase (DapB; #310043), and the positive control probe DNA-directed RNA polymerase II subunit RPB1 (POLR2A; #310451). Probes were incubated for 2 hours at 40°C using the HybEZ oven from ACD and six signal amplification steps, followed by signal detection with DAB were performed according to manufacturer's protocol. Slides were counterstained with 50% hematoxylin, dehydrated, and cover slipped.

Image analysis
Slides were first assessed by microscopic inspection. Since the staining pattern of RNA in situ hybridization was homogenous across different regions of the same sample, a representative region was selected for each tumor, and z-stacks comprising 10 images were captured at 40× magnification using a Zeiss Observer.Z1 inverted microscope (Zeiss, Jena, Germany). Multiple adjacent single images (3x3 tiles) were acquired and aligned using the MosaiX module of the AxioVision software (Zeiss). Each assembled image covered an area of approx. 0.62 × 0.45 mm (3981x2980 pixels). Image sequences were stacked using the minimal intensity projection type setting of the ImageJ software (http://imagej.nih.gov/ ij/). The SpotStudio software from ACD was used for detection of single cells, detection of spots and clusters and calculation of the estimated number of spots per cell. DapB and POLR2A probes served as technical quality controls that needed to fulfill the below given cut-off criteria in order to a) ensure technical specificity of the probes (negative control) and b) to detect samples with highly degraded RNA (positive control). The maximum threshold for negative controls was set at 0.5 spots/cell and the minimum threshold for positive controls at 2.5 spots/ cell. Samples not fulfilling these cut-off requirements were excluded from further analysis.

Statistical analysis
Statistical analyses were performed using the SigmaPlot/SigmaStat v12.5 software (Systat Software Inc., San Jose, CA, USA) or GraphPad Prism 7.02 for Windows (GraphPad Software, La Jolla, CA, USA) for comparison of groups, Spearman rank correlation analysis, generation of boxplots, bargraphs and scatter plots. All performed tests were two-sided and rank based. For comparisons of paired data, a Wilcoxon signed rank test, for comparisons of two groups a Wilcoxon rank sum test (Mann-Whitney U test), and for comparison for more than two groups a Kruskal-Wallis test followed by Dunn's posthoc tests for pairwise comparisons were used. Different patient characteristics between the TCGA cohort and the validation cohort (subset of GEO ID GSE17705) were compared using a chi-square test and a Wilcoxon rank sum test regarding age at diagnosis and follow-up times. Genesis 1.7.6 (Graz University of Technology, Graz, Austria) was used for log2-transformation and mean centering of gene expression values, calculation of Euclidean distance, hierarchical cluster analysis and heatmap visualization. For overall and disease free survival analysis, the statistical software environment R (www.r-project.org) including the package survival and an adopted code for the auto-cut-off option from Györffy et al. [40] was used. For the construction of survival curves, a Kaplan-Meier estimator was used and survival curves of patient groups with high expression versus low expression were compared by a logrank test. A time-independent Cox-proportional hazard approach was applied for univariate and multivariate survival analysis. Results were considered statistically significant when p < 0.05.