Transcriptome sequencing revealed differences in the response of renal cancer cells to hypoxia and CoCl₂ treatment

Introduction

Hypoxia is characterized by reduced oxygen supply and appears in multiple pathological conditions including tumours. However, hypoxia can also have a functional role during normal mammalian development and embryogenesis ¹ . Cells respond to hypoxic conditions both on biochemical and gene expression levels by switching from aerobic metabolism to anaerobic glycolysis and by expression of stress-related genes involved in regulation of cell death, erythropoiesis, angiogenesis and survival ^2–4 . The activation of many O₂-regulated genes is mediated by hypoxia-inducible factor (Hif1a). Under normoxia, Hif1a is hydroxylated by specific prolyl hydroxylases (PHD1, PHD2 and PHD3). This reaction requires oxygen, 2-oxoglutarate and ascorbate ^5,6 . When Hif1a is hydroxylated, it interacts with the von Hippel-Lindau tumor suppressor protein (pVHL). pVHL forms the substrate-recognition module of an E3 ubiquitin ligase complex, which directs Hif1a poly-ubiquitylation and proteasomal degradation ^7,8 . Under hypoxia (less than 5% O₂), PHD activity is inhibited by cytoplasmic reactive-oxygen species (ROS) which alter the oxidation state of Fe²⁺ (a cofactor for PHD activity) to Fe³⁺. This alteration inhibits PHD activity and Hif1a hydroxylation, thus Hif1a cannot interact with pVHL and promotes HIf1a stabilization ^9,10 . This anaerobic condition and stabilization of Hif1a are characteristic of many tumors. The most common molecular abnormality in renal cell carcinoma is the loss of VHL, which is found in about 50–70% of sporadic cases. Consequently, renal carcinomas with mutations in VHL have high steady-state levels of Hif1a expression and are hypoxic ¹¹ . Some divalent cations such as cobalt (Co²⁺), nickel (Ni²⁺), and the iron-chelator deferoxamine (DFX), have been applied to mimic hypoxic conditions in cultured cells as they activate hypoxic signals by stabilizing HIF1a ¹² . Transition metal Co²⁺ could induce hypoxic response by inhibiting PHD activity via iron replacement. Therefore, treatment of a cell culture with cobalt chloride (CoCl₂) is a common model of hypoxia ¹³ . The second classical setup to study hypoxia is hypoxia induction in a CO₂ incubator with a regulated level of oxygen (less than 1% O₂). In this work, we performed RNA sequencing of Caki-1 clear cell renal cancer cell lines treated with hypoxia and with CoCl₂ to understand how adequate CoCl₂ treatment was as a hypoxia model. We propose that CoCl₂ is not a completely correct model for hypoxia, as it aberrantly induces various hydroxylases not involved in hypoxia pathways and fails to induce downstream biochemical pathways normally induced by hypoxia.

Methods

Results

To compare transcriptional effects caused by hypoxia and by CoCl₂ exposure, we performed transcriptome sequencing (RNA-seq) of Caki-1 clear cell renal carcinoma cell line in three conditions: untreated, treated with CoCl₂, and exposed to hypoxic conditions (1% O₂). We searched for differentially expressed genes in two comparisons: CoCl₂-treated against untreated cells and cells under hypoxia against untreated cells.

We first asked how transcriptomic changes after both treatments fit into the established hypoxia gene signature ^1,18,19 . Figure 1A shows relative expression of the genes from hypoxia signature in all sequenced samples. We observed that genes in hypoxia signature were upregulated both in hypoxia and after CoCl₂ treatment, though the effect was much stronger in CoCl₂-treated samples. To check if a higher magnitude of transcriptomic changes in CoCl₂ compared to hypoxic conditions was characteristic to other differentially expressed genes, we compared the distributions of logarithmic fold-changes (absolute value) for gene expression in these two experiments (Figure 1B). Surprisingly, for the majority of genes, hypoxia-induced changes were higher than the ones induced by CoCl₂. In other words, hypoxia resulted in broader transcriptome response than CoCl₂ treatment even though specific hypoxia-related genes were more affected after CoCl₂ treatment.

Figure 1.

(A) CoCl₂ causes much more pronounced expression changes in the expression of key hypoxia regulators compared to real hypoxia treatment. Heatmap represents relative gene expression for key genes involved in hypoxia regulation. (B) Hypoxia results in broader transcriptome response compared to CoCl₂ treatment, i.e., more genes are changing expression under hypoxia. The figure shows absolute log fold change values for gene expression between hypoxia (or CoCl₂) groups relative to control group. Genes are sorted according to absolute log fold change values. P(wilcoxon) < 2.2 × 10^-16.

These results suggested that CoCl₂ treatment could be an incomplete model of hypoxia capturing only upstream signalling events in hypoxia pathways and not reflecting broader downstream effects. To further investigate the differences between CoCl₂ model and real hypoxia, we compared the sets of differentially up- and down-regulated genes in CoCl₂-treatment (against untreated) and in hypoxia-treated (against untreated) cells. Figure 2A summarizes the overlap between those gene sets. To understand what regulatory and biochemical pathways were affected in each treatment, we performed gene category enrichment analysis over KEGG ¹⁷ pathways with DAVID web service (version 6.7) ¹⁶ for genes affected in both treatments or exclusively in one treatment.

Figure 2.

(A) Summary of KEGG pathways, enriched by the genes, up- and down-regulated in CoCl₂ and hypoxia treatment. No enriched pathways were discovered for the genes downregulated in both treatments. (B) Overall expression change in glycolysis/gluconeogenesis (KEGG hsa00260) genes in control, hypoxia conditions and after CoCl₂ treatment. Glycolysis/Gluconeogenesis (KEGG hsa00010) is activated in hypoxia but not after CoCl₂ treatment. (C) Expression change for glycine, serine and threonine metabolism genes (KEGG hsa00010). Glycine, serine and threonine metabolism (KEGG hsa00260) is activated only after CoCl₂ treatment but not under hypoxia.

The genes which were significantly upregulated in hypoxic conditions (1% O₂) but not after treatment with CoCl₂ were significantly enriched in the glycolysis/gluconeogenesis pathway which was known to be related to hypoxia ^1,15,18 . Unexpectedly, we detected no enrichment in glycolysis/gluconeogenesis pathway with the genes differentially expressed after CoCl₂ treatment which was confirmed by a heatmap for glycolysis/gluconeogenesis-related genes (Figure 2B).

The genes upregulated in both hypoxia and under CoCl₂ treatment were enriched in the MAPK pathway. As the MAPK pathway was known to activate hypoxic response, MAPK activation was expected in hypoxia. Even though CoCl₂ affected directly HIF1a pathway, MAPK was activated in CoCl₂-treated samples as well as in real hypoxia. This result supported a previous observation of MAPK-dependent activation of hypoxia response under CoCl₂ treatment ¹² . Surprisingly, we observed systemic lupus erythematosus-related pathway activation in both treatments. The set of genes upregulated in this pathway (histone proteins H2A, H2B, H3, H4 and MHCII antigen-presenting genes) could be unrelated to lupus erythematosus, but rather could indicate increased proliferation and inflammation.

We also explored the genes specific to CoCl₂ treatment but not affected by hypoxia. Surprisingly, the genes upregulated after CoCl₂ treatment but not changed in hypoxia were enriched in the glycine/serine/threonine biosynthesis pathway (Figure 2A and 2B). We hypothesized that Co²⁺ ion could substitute metal cofactors of several oxidoreductases in the pathway and subsequently impair their activity. This, in turn, could require greater amounts of enzyme to be synthesized.

Hypoxic conditions within a tumour have been shown to predict worse clinical outcome ²⁰ . We used our data on whole-transcriptome profiling of kidney cancer cell lines in normal and hypoxic conditions to extract a wide hypoxia signature and validate it with TCGA data on transcriptomes of kidney tumours ²¹ . Major variation in TCGA transcriptomes (Figure 3A) is generated by the difference between tumours and adjacent normal samples. We projected our sequenced samples to the principal components derived from TCGA samples and observed that the direction of transcriptome changes between untreated and hypoxic cell lines is slightly similar to the difference between normal and tumour samples, though hypoxia-related changes couldn’t explain normal-tumour difference. We constructed a transcriptome-based hypoxia signature as described in the Methods section. To test if this hypoxia signature predicted clinical outcome for kidney cancer patients, we explored the distribution of our hypoxia scores between disease free patients and patients in which the disease had recurred or progressed. Hypoxia scores were significantly higher for recurred/progressed patients (Wilcoxon test p-value=0.0009, Figure 3B).

Figure 3.

(A) PCA plot for TCGA samples of kidney tumours, adjacent normal tissue samples, untreated Caki-1 cells, Caki-1 cells treated with CoCl₂ and the cells in hypoxic environment. (B) Hypoxia signature derived from RNA-seq results predicted significantly higher hypoxia scores for recurred or progressed TCGA tumours.

Discussion

The current study for the first time provides RNA-seq data revealing hypoxia-induced transcriptomic changes, which allows broader understanding of the processes related to hypoxia. In our analysis, we explored the limits of applicability of CoCl₂ treatment as the model of hypoxia. Briefly, we observed that CoCl₂ strongly alters expression of few genes important for hypoxia signalling, but fails to influence the essential downstream consequences of hypoxia, particularly the glycolysis/gluconeogenesis pathway. This might suggest the existence of alternative regulation mechanisms which trigger the downstream events in hypoxia together with main VHL/HIF1a pathway. CoCl₂ treatment also abberantly induced pathways which did not respond to hypoxia. These included the glycine, serine and threonine metabolism pathways. We hypothesized that its aberrant activation might be caused by Co²⁺ ion binding to the enzymes (other than PHD proteins) involved in Glycine, Serine and Threonine biosynthesis.

Data availability

F1000Research: Dataset 1. Raw per-gene expression counts for individual genes (see Methods), 10.5256/f1000research.7571.d109570 ²³

RNA-seq data was deposited to NCBI SRA under SRP066934 study accession code. The study contained experiments under the following accession codes: untreated (SRX1459966, SRX1459967, SRX1459969), treated with CoCl₂ (SRX1459974, SRX1459977, SRX1459978), exposed to hypoxia (SRX1459979, SRX1459981, SRX1459984).