Culture of Cancer Cells at Physiological Oxygen Levels Affects Gene Expression in a Cell-Type Specific Manner

Standard cell culture is routinely performed at supraphysiological oxygen levels (~18% O2). Conversely, O2 levels in most mammalian tissues range from 1–6% (physioxia). Such hyperoxic conditions in cell culture can alter reactive oxygen species (ROS) production, metabolism, mitochondrial networks, and response to drugs and hormones. The aim of this study was to investigate the transcriptional response to different O2 levels and determine whether it is similar across cell lines, or cell line-specific. Using RNA-seq, we performed differential gene expression and functional enrichment analyses in four human cancer cell lines, LNCaP, Huh-7, PC-3, and SH-SY5Y cultured at either 5% or 18% O2 for 14 days. We found that O2 levels affected transcript abundance of thousands of genes, with the affected genes having little overlap between cell lines. Functional enrichment analysis also revealed different processes and pathways being affected by O2 in each cell line. Interestingly, most of the top differentially expressed genes are involved in cancer biology, which highlights the importance of O2 levels in cancer cell research. Further, we observed several hypoxia-inducible factor (HIF) targets, HIF-2α targets particularly, upregulated at 5% O2, consistent with a role for HIFs in physioxia. O2 levels also differentially induced the transcription of mitochondria-encoded genes in most cell lines. Finally, by comparing our transcriptomic data from LNCaP and PC-3 with datasets from the Prostate Cancer Transcriptome Atlas, a correlation between genes upregulated at 5% O2 in LNCaP cells and the in vivo prostate cancer transcriptome was found. We conclude that the transcriptional response to O2 over the range from 5–18% is robust and highly cell-type specific. This latter finding indicates that the effects of O2 levels are difficult to predict and thus highlights the importance of regulating O2 in cell culture.


Introduction
Mammalian cell culture is often used to study cell physiology in health and disease. For this purpose, environmental parameters such as temperature and pH are regulated to recreate as closely as possible the in vivo conditions. However, while cells in most mammalian tissues are exposed to 1-6% oxygen in vivo [1], cell culture is routinely performed in incubators that regulate CO 2 but not O 2 . At sea level, atmospheric O 2 levels are~21%, and headspace O 2 in conventional incubators thus equilibrates to~18% O 2 due to high humidity and the addition of 5% CO 2 . Despite being referred to as 'normoxia', 18% O 2 is substantially hyperoxic relative to the cellular microenvironment in vivo. Increasing evidence indicates that the hyperoxic conditions of cell culture affect multiple biological processes, including reactive oxygen species (ROS) production [2], redox homeostasis [3], proliferation and differentiation [4], bioenergetics [5], mitochondrial network dynamics [5], and response to drugs [6] and hormones [7]. These effects of non-physiologically high O 2 levels can compromise the ability of cell culture models to recapitulate in vivo disease pathophysiology (reviewed in [8]). Cancer cell lines are widely used in basic research to study cancer pathophysiology; however, they are routinely cultured in 18% O 2 . To gain a better understanding of the effects of this supraphysiological O 2 levels in cell culture we studied the transcriptional response of four cancer cell lines to O 2 levels in standard cell culture (18% O 2 ) versus physioxia (5% O 2 ). While physioxia is a range (typically 1-6% O 2 ) rather than a set value, we chose 5% O 2 for this study as we have accumulated data from other studies using this value [2,5,7,9], which is typical of many mammalian tissues in vivo [1]. Furthermore, one challenge of maintaining physioxia in vitro is avoiding pericellular hypoxia due to the lower gradient for O 2 entry into media. We have previously shown that culture at 5% O 2 in the incubator headspace maintains pericellular O 2 in media at above 3.7% O 2 under conditions used here [10].
We used RNA-seq and bioinformatic approaches to analyze differential gene expression of four human cancer cell lines cultured for 14 days at either 5% or 18% O 2 . The cell lines used in this project were LNCaP (prostate adenocarcinoma), Huh-7 (hepatocellular carcinoma), PC-3 (prostate adenocarcinoma), and SH-SY5Y (neuroblastoma). We chose these cell lines because we have previously measured the effects of O 2 on a wide range of cellular activities in them, including energy metabolism, mitochondrial dynamics, ROS production, and response to drugs such as resveratrol [5,9]. Inclusion of SH-SY5Y cells was done to increase the breadth of cell type. Further, these cell lines are not only widely used in cancer research, but also as surrogates of primary cells to study their physiology and pathology. For example, SH-SY5Y cells are used in models of neurodegenerative disease [11] and ischemia/reperfusion models [12] and Huh-7 are frequently utilized to study hepatic xenobiotic metabolism [13] and fatty liver disease [14]. Using this approach, we asked whether the transcriptomes of all four cell lines were sensitive to O 2 levels in the range from 5% to 18% O 2 , and whether the effects of O 2 were similar amongst cell lines, or cell-line specific. We found that the effects of O 2 on the transcriptomes of these cell lines were substantial and largely cell-line specific. Nonetheless, we identified interesting patterns in gene expression among these cell lines, such as a differential but strong induction of genes encoded by mitochondrial DNA (mtDNA) and upregulation of HIF-1/2 targets in physioxia.

Cell Culture
LNCaP, SH-SY5Y, Huh-7, and PC-3 cell lines were purchased from ATCC (Manassas, VA, USA). Cell passages between 15-19 were used throughout this study. Upon thawing, cells were cultured in 10-cm plates with Plasmax (Ximbio, London, UK) supplemented with 2.5% FBS and 1% penicillin/streptomycin (Sigma-Aldrich; St. Louis, MO, USA) for a week in a humidified 5% CO 2 (~18% O 2 ) to allow for acclimatization to Plasmax (physiological cell culture medium; see [15]) and reduced FBS concentration. Afterwards, cell lines were incubated in a humidified 5% CO 2 incubator at either 5% or 18% O 2 . Three replicates per each cell line were used in each condition. For the experimental groups kept at 5% O 2 , Plasmax media was preincubated overnight in the 5% O 2 incubator to allow for gas equilibration. Sub-culture was performed with 0.25% Trypsin-EDTA (Sigma-Aldrich, St. Louis, MO, USA) every 3 or 4 days, when cells reached~80% confluence. Media was refreshed every 24 h and cells were routinely monitored for mycoplasma contamination. Cell culture at either 5% O 2 or 18% O 2 was performed for 14 days. Cells were seeded at a density of 2 × 10 6 cells/plate prior to RNA extraction.

RNA Isolation
Total RNA was extracted using the RNeasy Plus Mini Kit (QIAGEN, Toronto, ON, Canada) according to the manufacturer's instructions. RNA integrity was assessed using 1.5% agarose gel electrophoresis, while RNA concentration and purity were evaluated as A260/280 ratio using a Thermo Fisher Scientific Nanodrop spectrophotometer. RNA samples were snap frozen in liquid nitrogen and stored at −80 • C until being sent to Novogene (Sacramento, CA, USA) for sequencing and analysis.

Sequencing and Differential Gene Expression Analysis
Quality check (QC), library preparation, sequencing, and differential expression analysis were performed by Novogene. Paired-end at 150 bp (PE150) high throughput Illumina sequencing was performed at a sequencing depth of 40 million reads per sample. Reads were aligned to the Homo sapiens reference genome (GRCh38) using Hisat2 v2.0.5 [16]. Gene expression levels were estimated by calculating FPKM (fragments per kilobase of transcript per million mapped sequence reads), which were further adjusted by edgeR program package [17] through one scaling normalized factor. Differential expression analysis was performed using the edgeR R package. A p-value < 0.05 and |log 2 FC|≥ 1 were set as threshold for significantly differential expression, as done previously [10].

Functional Enrichment Analysis
Prior to functional enrichment analysis, the list of differentially expressed genes (DEGs) was further reduced to genes with a Benjamini adjusted p-value (p adj ) < 0.1 and a FPKM ≥ 1 in at least one of the experimental groups in order to produce a concise list of enrichment terms which reflect the most strongly affected genes. Functional enrichment analysis was performed using the Database for Annotation, Visualization, and Integrated Discovery (DAVID) [18]. Enriched Gene Ontology (GO) terms, Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways, and Reactome pathways were selected for functional annotation. A raw p-value and Benjamini adjusted p-value (p adj ) of 0.05 were applied for identifying the most statistically significant enriched annotation terms.

Correlation Analysis of PC-3 and LNCaP Data with Gene Expression Data from the Prostate Cancer Transcriptome Atlas (PCTA)
Data from the Prostate Cancer Transcriptome Atlas (PCTA) was downloaded from [19], containing a dataset of expression levels of 18,390 genes from 2115 human prostate tumor samples. Upon removing genes with FPKM = 0 in any of the experimental groups, our lists of DEGs upregulated at 5% O 2 and 18% O 2 in both LNCaP and PC-3 were independently cross-referenced with the PCTA dataset to find overlapping genes. We then analyzed the correlation between the change in expression levels (log 2 FC) in our samples and the PCTA dataset (mean log 2 FC from the 2115 samples) by performing a Pearson correlation test with the software GraphPad Prism 8.

Oxygen Levels in Culture Strongly Modulated Transcript Abundance Cell-Line Specifically
The abundance of over a thousand transcripts were affected by O 2 in each cell line. In general, more differentially expressed genes (DEGs) showed higher expression at 5% O 2 than at 18% in all cell lines. In addition, there was substantial variation between the four cell lines in their sensitivities to O 2 ( Figure 1A). For example, 2126 DEGs were identified in LNCaP cells, including 433 upregulated at 18% O 2 and 1693 at 5%. In contrast, SH-SY5Y was shown to be the least sensitive to O 2 among the cell lines, with only 386 transcripts upregulated at 18% O 2 and 848 at 5%. The full lists of DEGs for LNCaP, Huh-7, PC-3, and SH-SY5Y cells are available in Tables S1-S4, respectively.
A remarkable result was the extremely limited overlap between cell lines in terms of the identities of the DEGs ( Figure 1B). Only four genes were identified as being O 2 -sensitive in all four cell lines. Of these, BOLA2B was the only protein-coding gene. The BolA protein family has important roles in Fe-S cluster biogenesis, iron and Fe-S cluster trafficking and storage, and iron sensing and regulation [20]. Differential expression of BOLA2B may thus have an impact in redox sensing and signaling. Interestingly, BOLA2B was found to be upregulated at 5% O 2 in all cell lines, except Huh-7, where it was upregulated at 18% O 2 . Even amongst the two prostate cancer cell lines, LNCaP and PC-3, where 2126 and A remarkable result was the extremely limited overlap between cell lines in terms of the identities of the DEGs ( Figure 1B). Only four genes were identified as being O2-sensitive in all four cell lines. Of these, BOLA2B was the only protein-coding gene. The BolA protein family has important roles in Fe-S cluster biogenesis, iron and Fe-S cluster trafficking and storage, and iron sensing and regulation [20]. Differential expression of BOLA2B may thus have an impact in redox sensing and signaling. Interestingly, BOLA2B was found to be upregulated at 5% O2 in all cell lines, except Huh-7, where it was upregulated at 18% O2. Even amongst the two prostate cancer cell lines, LNCaP and PC-3, where 2,126 and 1,461 genes were differentially expressed, respectively, only 192 were shared between both cell lines. Similarly, of the 2,099 transcripts affected by O2 in Huh-7 cells, 1638 (78%) were exclusively affected in this cell line. This indicates that O2 effects on gene expression are highly specific to a given cell line. This in turn makes it difficult to predict how the non-physiological O2 levels of standard cell culture are affecting cell biology in general terms.
Functional enrichment analysis revealed that different biological processes and pathways were enriched by O2 level in the four cell lines ( Figure 2). For example, the most significantly affected pathway in LNCaP cells was TGF- signaling (hsa04350; padj < 0.05), which was found to be enriched at 18% O2. In Huh-7 cells, pathways such as extracellular matrix (ECM) organization (R-HSA-1474244; padj < 0.005) and drug metabolism by the cytochrome P450 (CYP450) enzymes (hsa00982; padj < 0.05) were strongly enriched at 5% O2, while oxidative phosphorylation (hsa00190; padj < 0.005) and oxidative stress-induced senescence (R-HSA-2559580; padj < 0.05) were enriched at 18% O2. Interestingly, in contrast to Huh-7 cells, both PC-3 and SH-SY5Y showed enrichment of annotation terms related to mitochondrial respiration and oxidative phosphorylation at 5% O2 (see Figure 2). Signaling by interleukins (R-HSA-449147; padj < 0.05) and neurogenesis (GO:0022008; padj < 0.005) were among the processes enriched at 18% O2 in SH-SY5Y cells. The full lists of functional annotation terms enriched by O2 level in all cell lines are available in Tables S5-S12. Functional enrichment analysis revealed that different biological processes and pathways were enriched by O 2 level in the four cell lines ( Figure 2). For example, the most significantly affected pathway in LNCaP cells was TGF-β signaling (hsa04350; p adj < 0.05), which was found to be enriched at 18% O 2 . In Huh-7 cells, pathways such as extracellular matrix (ECM) organization (R-HSA-1474244; p adj < 0.005) and drug metabolism by the cytochrome P450 (CYP450) enzymes (hsa00982; p adj < 0.05) were strongly enriched at 5% O 2 , while oxidative phosphorylation (hsa00190; p adj < 0.005) and oxidative stress-induced senescence (R-HSA-2559580; p adj < 0.05) were enriched at 18% O 2 . Interestingly, in contrast to Huh-7 cells, both PC-3 and SH-SY5Y showed enrichment of annotation terms related to mitochondrial respiration and oxidative phosphorylation at 5% O 2 (see Figure 2). Signaling by interleukins (R-HSA-449147; p adj < 0.05) and neurogenesis (GO:0022008; p adj < 0.005) were among the processes enriched at 18% O 2 in SH-SY5Y cells. The full lists of functional annotation terms enriched by O 2 level in all cell lines are available in Tables S5-S12.

The Top Differentially Expressed Genes Have Key Roles in Cancer Cell Biology
By sorting the DEGs according to their adjusted p-value, we found that most of the genes highly affected by O 2 are implicated in cancer cell biology, including several with roles in cancer cell proliferation, tumor progression, metastasis, invasion, and chemosensitivity to anticancer therapy. A selection of these genes is shown in Table 1. The complete list of top 10 protein-coding DEGs in all cell lines at both O 2 conditions and their corresponding log 2 FC values are shown in Figure 3.

The Top Differentially Expressed Genes Have Key Roles in Cancer Cell Biology
By sorting the DEGs according to their adjusted p-value, we found that most of the genes highly affected by O2 are implicated in cancer cell biology, including several with roles in cancer cell proliferation, tumor progression, metastasis, invasion, and chemosensitivity to anticancer therapy. A selection of these genes is shown in Table 1. The complete list of top 10 protein-coding DEGs in all cell lines at both O2 conditions and their corresponding log2FC values are shown in Figure 3. Huh-7 S100A9 ** S100 calcium binding protein A9 TLR4 and RAGE ligand, promotes HCC progression through MAPK and NF-B pathways. +3.37 [23] SLC3A2 ** (GLUT3) Solute carrier family 2 member 3 Selective glucose uniporter. Expression is correlated with HCC growth/invasion. −5.74 [24] PC-3  CHAC1 ** ChaC glutathione specific gamma-glutamylcyclotransferase 1 Degrades glutathione. Involved in ferroptosis; associated with increased chemosensitivity. −2.81 [22] Huh-7 S100A9 ** S100 calcium binding protein A9 TLR4 and RAGE ligand, promotes HCC progression through MAPK and NF-κB pathways.

Oxygen Levels Affected mtDNA-Encoded Transcript Abundances in Most Cell Lines
Oxygen induced differential expression of mtDNA-encoded genes in most cell lines ( Table 2). In Huh-7 cells, 11 mtDNA-encoded gene transcripts were affected by O2, of which six are subunits of respiratory complexes I, IV, and V, while the rest are mitochondrial tRNAs. Interestingly, all of these were upregulated at 18% O2. Eleven mtDNA-encoded genes were affected by O2 in PC-3 cells and 10 in SH-SY5Y, however, all were upregulated at 5% O2, in striking contrast with the observation in Huh-7 cells. Again, these

Oxygen Levels Affected mtDNA-Encoded Transcript Abundances in Most Cell Lines
Oxygen induced differential expression of mtDNA-encoded genes in most cell lines ( Table 2). In Huh-7 cells, 11 mtDNA-encoded gene transcripts were affected by O 2 , of which six are subunits of respiratory complexes I, IV, and V, while the rest are mitochondrial tRNAs. Interestingly, all of these were upregulated at 18% O 2 . Eleven mtDNA-encoded genes were affected by O 2 in PC-3 cells and 10 in SH-SY5Y, however, all were upregulated at 5% O 2 , in striking contrast with the observation in Huh-7 cells. Again, these DEGs encoded subunits of the respiratory chain and tRNAs. In contrast, only two mtDNA-encoded genes were affected in LNCaP. These results suggest that O 2 levels in cell culture affect the expression of mtDNA-encoded genes, but in a highly cell-type specific manner.

HIF Targets Were Upregulated at 5% O 2 in All Cell Lines
Five percent O 2 is not hypoxic, and we have previously shown that, under the conditions used here, pericellular O 2 levels do not fall below 3.7%, even during extensive static incubation periods in culture [10]. A handful of reports have shown that HIF-1/2 activity is detectable at 2-5% O 2 levels (i.e., physioxia) [6,29,30]. Here we identified several HIF-1/2 gene targets upregulated at 5% O 2 in all cell lines, a selection of which is shown in Table 3. For example, in LNCaP, transcripts related to angiogenesis and vasodilation, such as VEGFA and ADM, were enriched. Similarly, genes that encode enzymes involved in the metabolic reprograming of cells towards a glycolytic phenotype were upregulated in Huh-7 cells grown at 5% O 2 . These genes include the glucose transporter SLC2A3 (GLUT3), the glycolytic enzyme ENO2 (enolase), and the gluconeogenic enzyme PCK1 (phosphoenolpyruvate carboxykinase 1). A lower number of HIF-1/2 targets were detected in PC-3 and SH-SY5Y cells, consistent with our initial observation that these two cell lines were less sensitive to O 2 than Huh-7 and LNCaP cells.   We next investigated whether the HIF targets upregulated at 5% O 2 in our cells were regulated by HIF-1α, HIF-2α, or both. To do this, we compared our DEG datasets with the dataset from the study by Downes et al. who sequenced the transcriptional outputs of stabilized forms of HIF-1α and HIF-2α [75]. We obtained a more exhaustive list of HIF-1/2 targets upregulated at 5% O 2 in our cell lines (Table S14). A total of 103, 145, 40, and 29 HIF-1/2 targets were identified in LNCaP, Huh-7, PC-3, and SH-SY5Y, respectively. When comparing the proportion of unique and shared DEGs regulated by HIF-1α and HIF-2α, we found that the number of unique HIF-2α targets was greater in all cell lines, although several gene targets of both were also observed ( Figure 4). We next investigated whether the HIF targets upregulated at 5% O2 in our cells were regulated by HIF-1, HIF-2, or both. To do this, we compared our DEG datasets with the dataset from the study by Downes et al. who sequenced the transcriptional outputs of stabilized forms of HIF-1 and HIF-2 [75]. We obtained a more exhaustive list of HIF-1/2 targets upregulated at 5% O2 in our cell lines (Table S14). A total of 103, 145, 40, and 29 HIF-1/2 targets were identified in LNCaP, Huh-7, PC-3, and SH-SY5Y, respectively. When comparing the proportion of unique and shared DEGs regulated by HIF-1 and HIF-2, we found that the number of unique HIF-2 targets was greater in all cell lines, although several gene targets of both were also observed ( Figure 4).

Expression Patterns of Genes Upregulated at 5% O 2 in LNCaP Cells, but Not PC-3, Better Correlate with the In Vivo Prostate Cancer Transcriptome
To determine if physiological O 2 levels in cell culture produce a transcriptional signature that better resembles a prostate cancer transcriptome in vivo, we compared our transcriptomic data from LNCaP and PC-3 cells (two prostate adenocarcinoma cell lines) with the expression data from the Prostate Cancer Transcriptome Atlas (PCTA), which contains a dataset of expression levels of 18,390 genes from 2,115 human prostate tumor samples [19]. We independently matched the DEGs upregulated at 5% O 2 and at 18% O 2 from both cell lines with the PCTA dataset. Upon obtaining matched gene lists from our data and the PCTA dataset, we performed a correlation analysis of the change in gene expression levels (log 2 FC) in both datasets. In LNCaP cells, gene expression patterns from the DEGs upregulated at 5% O 2 showed a highly statistically significant (p < 0.0001) albeit weak (r = 0.1837) correlation with the PCTA expression data ( Figure 5A). On the other hand, expression levels of DEGs upregulated at 18% O 2 showed no statistically significant correlation (p = 0.1331) with the PCTA expression data ( Figure 5B). Finally, the correlation between the expression levels of DEGs upregulated at either 5% O 2 or 18% O 2 in PC-3 cells, and the PCTA dataset was non-significant (p > 0.05).

Expression Patterns of Genes Upregulated at 5% O2 in LNCaP Cells, but not PC-3, Better Correlate with the In Vivo Prostate Cancer Transcriptome
To determine if physiological O2 levels in cell culture produce a transcriptional signature that better resembles a prostate cancer transcriptome in vivo, we compared our transcriptomic data from LNCaP and PC-3 cells (two prostate adenocarcinoma cell lines) with the expression data from the Prostate Cancer Transcriptome Atlas (PCTA), which contains a dataset of expression levels of 18,390 genes from 2,115 human prostate tumor samples [19]. We independently matched the DEGs upregulated at 5% O2 and at 18% O2 from both cell lines with the PCTA dataset. Upon obtaining matched gene lists from our data and the PCTA dataset, we performed a correlation analysis of the change in gene expression levels (log2FC) in both datasets. In LNCaP cells, gene expression patterns from the DEGs upregulated at 5% O2 showed a highly statistically significant (p < 0.0001) albeit weak (r = 0.1837) correlation with the PCTA expression data ( Figure 5A). On the other hand, expression levels of DEGs upregulated at 18% O2 showed no statistically significant correlation (p = 0.1331) with the PCTA expression data ( Figure 5B). Finally, the correlation between the expression levels of DEGs upregulated at either 5% O2 or 18% O2 in PC-3 cells, and the PCTA dataset was non-significant (p > 0.05).

Discussion
The main goal of this study was to determine how O 2 tension in the standard cell culture environment (18% O 2 ) impacts the cancer cell transcriptome compared to a more representative in vivo environment (5% O 2 ). We were particularly interested in the extent to which any effects of O 2 were shared amongst cell lines versus cell-line specific. Our results indicate broad transcriptional effects of O 2 between 5% and 18% that are highly cell-type specific. Even in LNCaP and PC-3, both prostate cancer cell lines, the overlap in O 2dependent DEGs was only~5%. These results are consistent with a previous study showing little overlap in the proteome of three diffuse large B-cell lymphoma cell lines cultured in the same two O 2 conditions [76]. Increased oxidative stress associated with higher O 2 levels affects a variety of pathways, including the p53 pathway and mitogen-activated protein kinase (MAPK) pathways (reviewed in [77]). Given their different origins and genetic backgrounds, cancer cells have distinct mutations, including gene copy number differences, of genes related to these pathways. Therefore, differential transcriptional response of these cell lines to O 2 is perhaps not surprising. In any case, these observations highlight the need for considering oxygenation status as an important factor in experimental design, since the effects of growing cells at 18% O 2 are broad and may not be easy to predict.
Functional enrichment analysis revealed that biological processes and pathways relevant to the disease etiology of each cell type were altered by O 2 level. For example, prostate cancer commonly metastasizes to bone, forming primarily osteoblastic lesions. TGF-β and bone morphogenetic proteins (BMPs), released by prostate cancer cells, induce osteoblast differentiation, which in turn releases growth factors that stimulate the proliferation of prostate cancer cells [78]. The TGF-β signaling pathway (hsa04350) was the annotation term most significantly enriched (p adj < 0.05) in LNCaP cells at 18% O 2 . Signaling by BMP (R-HSA-201451) and osteoblast differentiation (GO:0001649) were also enriched (p < 0.05). Another interesting observation was the upregulation of the androgen receptor (AR gene) in PC-3 cells grown at 18% O 2 (see Table S3), which is contrasting with fact that PC-3 is derived from an androgen-insensitive tumor and do not express AR. These interactions of key prostate cancer genes and signaling pathways with O 2 levels attest to the potential issues associated with a hyperoxic cell culture environment. Functional studies are necessary to determine how O 2 levels in cell culture affect metastasis and other functional characteristics of prostate cancer cells in vitro One of the main functions of hepatic cells is detoxification of xenobiotics through phase I (CYP450) and phase II enzymes. In Huh-7 cells, drug metabolism by the CYP450 (hsa00982) was among the annotation terms enriched at 5% O 2 (p adj < 0.05), in accordance with a previous observation that CYP1A1, CYP1A2, and CYP2E1, along with a number of phase II enzymes, were upregulated in HepG2 cells cultured in physioxia, compared with cells at atmospheric O 2 [79]. Differential expression of phase I and II enzymes at different O 2 tensions leads to altered hepatic metabolism of drugs and toxins, which in turn results in altered biological responses to these compounds when tested in vitro. Indeed, DiProspero et al. showed that the toxicity and potency of acetaminophen, cyclophosphamide, and aflatoxin B1 were dependent on O 2 tension in HepG2 cells cultured at~18% O 2 , 8% O 2 , and 3% O 2 [79]. Notably, these pharmacodynamic parameters obtained from cells grown at physioxic conditions better matched in vivo primary human hepatocyte data than cells cultured under standard conditions. Therefore, O 2 tension should be considered as an important factor in the design of experiments aimed to study the effects and mechanisms of bioactive molecules in vitro.
We observed several HIF-1/2 targets upregulated at 5% O 2 in LNCaP and Huh-7 cells. Although traditionally known as transcription factors that mediate the response to hypoxia, some studies have shown HIF-1/2 expression and activity in the physiological O 2 range [6,29,30]. A thorough discussion of the roles of HIF and other pathways in physioxia versus normoxia (18% O 2 ) is provided in [8]. The fact that fewer HIF-1/2 targets upregulated at 5% O 2 were observed in PC-3 and SH-SY5Y cells supports the notion that these cells are less sensitive to O 2 than LNCaP and Huh-7. Moreover, because cells were grown at either 5% O 2 or 18% O 2 for 14 days, the induction of HIF targets shown in this study suggests a role for HIF activity in physioxia, rather that it being a mere consequence of an acute reduction in O 2 levels while performing the experiments. Similar observations have been made in other studies [30,80], supporting this notion. Remarkably, the proportion of unique HIF-2α targets, compared to HIF-1α targets, was higher in all the cell lines (see Figure 4). These results are in line with previous observations where HIF-1α is described as the most active isoform in acute and severe hypoxia, whereas HIF-2α has been shown to be predominantly active during chronic and "physiological" hypoxia [81,82]. Downes et al. showed through functional enrichment analysis that the genes induced by HIF-1α are associated with biological processes like glycolysis and NADH regeneration, while HIF-2α-enriched processes include angiogenesis, extracellular matrix organization, pattern specification process, and negative regulation of cell adhesion. Indeed, GO terms and KEGG pathways related to all the above processes were found to be enriched by the HIF-regulated genes upregulated at 5% O 2 in our cell lines (see Tables S15-S18). In addition, given the variety of mechanisms in tumorigenesis regulated by HIF (e.g., metabolism, migration, invasion, survival), complete loss of HIF-1/2 activity at 18% O 2 may compromise experiments focused on cancer biology and chemotherapeutic strategies.
Enrichment of functional annotation terms related to the mitochondria and/or mitochondrial processes (e.g., respiration) was observed in three of the four cell lines (Huh-7, PC-3 and SH-SY5Y), although the directionality of the effects were inconsistent. While mitochondrial terms were enriched at 18% O 2 in Huh-7 cells, they were enriched at 5% O 2 in PC-3 and SH-SY5Y. The same general trend was apparent for mtDNA-encoded genes. Decreased expression of mtDNA-encoded genes in Huh-7 cells at 5% O 2 may reflect decreased mitochondrial abundance. Indeed, mitochondrial biogenesis (R-HSA-1592230) was one of the Reactome pathways enriched at 18% O 2 in Huh-7 cells (p adj < 0.05; see Table S8). In agreement, Moradi et al. observed decreased mitochondrial footprint in Huh-7 cells grown at 5% O 2 compared to 18% O 2 in the same conditions [5]. Further, it has been shown that HEY1, a HIF-1 target, decreases mitochondrial biogenesis by repressing the expression of PTEN-induced kinase 1 (PINK1) in human hepatocellular carcinoma cells through a HIF-1-dependent mechanism. In accordance, our data shows HEY1 expression increased 2.05-fold at 5% O 2 (p < 0.05; see Table 3). Additionally, Zhang et al. reported that HIF-1 inhibits mitochondrial biogenesis in renal carcinoma cells by repressing PGC-1β [83]. Interestingly, expression of PPARGC1B was~3 times lower in Huh-7 cells at 5% O 2 compared to 18% (p adj < 0.005; see Table S2). On the other hand, expression of all mtDNA-encoded genes affected by O 2 was higher in physioxia in both PC-3 and SH-SY5Y cells. Expression of genes related to mitochondrial biogenesis induced by HIF-1 has been observed in the neuroblastoma cell line SK-N-AS when exposed to hypoxia, along with increased mtDNA copy numbers [84]. Moreover, mitochondrial abundance was found to be higher in primary neurons grown at 2% and 5% O 2 compared to atmospheric O 2 [85]. No direct link between O 2 tension and regulation of mitochondrial biogenesis has been reported in prostate cancer. Future research should be directed to investigating the O 2 -dependent mechanisms of mitochondrial gene expression and biogenesis in different cell types.
Another aspect of cancer biology that has been shown to be affected by O 2 levels is cancer cell stemness [86]. As happens with non-transformed stem cells, the stemness of cancer cells is regulated by transcription factors like the octamer-binding transcription factor 4 (Oct4), Nanog, and SRY-box 2 (Sox2) [87]. Crosstalk between the HIF-1/2 pathway and the Oct4/Nanog/Sox2 axis is well documented [86]. For instance, Li et al. demonstrated that HIF-1/2 regulates the tumorigenic capacity of glioma stem cells and showed that HIF-2α colocalizes with cancer stem cell markers [88]. Using data from the study by Sharov et al. [89], we analyzed the effects of O 2 levels in our experimental groups on the expression of Oct4/Nanog/Sox2 related genes. While the expression of Oct4/Nanog/Sox2 was not affected by O 2 levels in most of the cell lines used here-except SH-SY5Y where POU5F1 (Oct4) was upregulated at 5% O 2 -we did observe differential expression of their gene targets at 5% and 18% O 2 (see Table S13). A similar observation was made by Westfall et al.
where they observed constant expression of Oct4/Nanog/Sox2, but differential expression of some of their targets in human embryonic stem cells grown in 4% O 2 versus 18% O 2 [90]. Altogether, these results suggest that culturing cancer cells at non-physiological O 2 levels may alter the transcriptional programs controlling stemness and renewal, and may in turn affect the response of cancer cells to chemotherapeutic drugs in vitro.
Finally, we analyzed the correlation between the transcriptomic data from the two prostate cancer cell lines used here and in vivo prostate tumor transcriptional profiles, by using the dataset from the PCTA [19,91]. In LNCaP cells, we found that the list of DEGs upregulated in physioxia (5% O 2 ) shows a highly significant correlation (p < 0.0001) with the PCTA gene expression data, whereas the DEGs upregulated at 18% O 2 has no significant correlation with the PCTA data. The overall weak correlation degree (Pearson r = 0.1837) of DEGs upregulated in physioxia and PCTA data could be explained by the fact that LNCaP is derived from a metastatic tumor from a single patient whereas the PCTA is a collection of data from a total of 2,115 patients. This shows the high heterogeneity of cancers within different patients. No significant correlation with the PCTA data was found from the transcriptomic profiles of PC-3 cells at either O 2 level, which supports the notion that PC-3 cells are not as sensitive to O 2 tension (between 5-18% O 2 ) as LNCaP and Huh-7 cells. Additional research should be conducted using multi-omics approaches to investigate how culturing cells in physioxia produces more similar genomic, transcriptomic, proteomic, and epigenomic profiles to the ones observed in tissues in vivo.
This study has some limitations. We recognize that transcriptional effects do not necessarily reflect changes in the cell proteome or phenotype, and this may be due to complex epigenetic, post-transcriptional, and post-translational mechanisms governing gene expression and protein function [92,93]. Investigation of the effects of O 2 levels on the epigenome, proteome, metabolome, and lipidome of different cell types is certainly warranted. In addition, data obtained from functional enrichment analysis needs to be interpreted carefully, as their significance in terms of cell behavior is not always clear. As an example, "regulation of apoptotic process" (GO:0042981) is one of the biological processes (GO) found to be enriched in PC-3 cells at 18% O 2 , with 11 DEGs associated with this process being upregulated at this O 2 tension. Out of these eleven genes, five have an antiapoptotic role, five are proapoptotic, and one has dual roles in the regulation of apoptosis. Thus, the functional effect of 18% O 2 in the regulation of apoptosis in PC-3 cells on apoptosis cannot be predicted by our bioinformatic analysis alone, and an experiment challenging PC-3 cells to apoptotic stimuli at both 5% O 2 and 18% O 2 would be necessary. Nonetheless, our results provide insight into the cellular processes and pathways that may be affected by O 2 in these cancer cells, which can direct subsequent research questions.
In conclusion, our results show that supraphysiological O 2 levels in cell culture significantly alter the global transcriptomes of cancer cell lines in highly cell-line specific ways. This makes it difficult to establish general rules regarding how non-physiological O 2 levels might affect experiments. Our results, together with the increasing amount of functional data regarding the effects of physioxia versus standard cell culture hyperoxia [8], should encourage cell culturists to implement the regulation of O 2 levels in their experiments. This would certainly be expected to increase the likelihood that results will translate to in vivo.