Hypoxia plays a role in the development and course of a number of pathological conditions, such as cardiovascular disorders [1] and tumor developments [2]. Hypoxia models staged in vitro allow studying the response to hypoxia at the molecular and cellular levels. One of the traditional models involves chemical agents that activate the hypoxia signaling pathways. CoCl₂ is one of the most common chemical hypoxia induction agents, since this substance causes direct and long-term stabilization of the hypoxia-induced factors 1 and 2 (HIF-1, HIF-2) [3].
MiRNAs are short non-coding RNAs; their main function is to suppress genes post-transcriptionally [4]. Usually, one miRNA has dozens of target genes, while 3'-UTR of a gene may have binding sites for hundreds of miRNAs [5]. It has been shown that interactions between miRNAs and their target genes play an important role in intercellular communication [6] and pathogenesis of many diseases, including various types of tumors [7, 8].
A number of studies aimed to investigate the role and functional activity of cellular miRNAs under the hypoxic stress conditions. It was found that some miRNAs, such as miR-210 or miR-27, are altered by hypoxia in many cells, and the differential expression of miRNAs and their targeting usually depend on the mechanism of hypoxia induction and cell type [9]. A relationship between the patterns of hypoxia-driven changes in miRNA expression and tumors has also been established: most of the miRNAs associated with the tumors can be affected by hypoxia [10].
This study aimed to investigate the effect of hypoxia on the miRNA profile and transcriptome in the HT-29 human colorectal adenocarcinoma line cells and to identify the potential key molecules involved in the response to hypoxia.

METHODS

Cell cultivation and processing

Human colorectal adenocarcinoma cells HT-29 (ATCC; USA) were cultured in McCoy's 5A medium (Thermo Fisher Scientific; USA) containing 10% fetal bovine serum (Thermo Fisher Scientific; USA). Penicillin (100 U/ml) and streptomycin (100 mg/ml) were added to the nutrient medium. The cells were seeded into 6-well plates in the amount of 4 × 10⁵ cells per well, then cultured in a humidified atmosphere at +37 °C and 5% CO₂ for 48 h. To induce hypoxia, a water solution of cobalt chloride (CoCl₂) was prepared, added to the medium to the concentration of 300 μM, and incubated for 24 h. Three biological iterations were performed for both control and treatment group cells.

RNA isolation

The cells were lysed in the Qiazol Lysis Reagent (Qiagen; Germany) for subsequent total RNA extraction with the help of the Qiagen miRNeasy Mini Kit (Qiagen, Hilden; Germany). The amount of isolated RNA was determined using a Nanodrop device (Thermo Fisher Scientific; USA). Agilent High Sensitivity DNA Kit (Agilent Technologies; USA) and a Bioanalyzer 2100 capillary electrophoresis instrument (Agilent Technologies; USA) enabled analysis of quality of the isolated RNA samples. The RIN (RNA integrity number) parameter value for all samples was above 9.0.

Preparation of libraries and sequencing

The libraries for mRNA sequencing were obtained from the total RNA samples using the Illumina Stranded mRNA Library Prep Kit Illumina (Illumina; USA). Each sample was sequenced on an Illumina NextSeq 550 to obtain 75 nucleotide reading frames at one end.
The libraries for miRNA sequencing were prepared from the total RNA samples using the NEBNext Multiplex Small RNA Library Prep Kit for Illumina. Each sample was sequenced in the Illumina NextSeq 550 to obtain 50 unidirectional nucleotide reading frames.
We sequenced mRNA and miRNA for three biological iterations, with each of them having four technical iterations set up.

Processing of the sequencing results

FASTQ file quality was assessed with the help of FastQC v0.11.9 (Babraham Bioinformatics; UK). One of the CoCl₂- treated replicates did not pass quality control at the microRNA sequencing stage. The adapters were cut using Cutadapt v2.10 [11]. The resulting mRNA fragment sequences were mapped to the human genome (GENCODE GRCh38.p13) using STAR v2.7.5b [12]. MiRDeep2 v2.0.1.2 package [13] enabled compilation of the miRNA expression matrix.
The sequencing library depths were normalized with the Trimmed Mean of M-values (TMM) algorithm available in the edgeR v3.30.3 package [14] with default background noise filtering. The same package was used to generate normalized mRNA and miRNA expression matrices in Reads Per Kilobase of transcript per Million mapped reads (RPKM) and Reads Per Million mapped reads (RPM) units, respectively. We took the logarithms of the values obtained with 2 as the base. For further processing, we only used the highly expressed transcripts and cut off the bottom 5% of genes and 50% of miRNAs based on their mean RPKM/RPM values.

Evaluation of the differential expression and overrepresented signaling pathways

We used DESeq2 v1.28.1 [15] to analyze the differential expression and applied the Benjamini-Hochberg procedure to determine the false detection rate (FDR). The differences with FDR below the threshold value of 0.05 were considered significant. DAVID v6.8 online service [16] enabled the analysis of overrepresented signal pathways.

Prediction of the miRNA targets

At the first stage of miRNA target prediction we exported a list of miRNA-gene interactions from TargetScan v7.2 [17]. Then, we selected the negative expression correlation miRNA-gene pairs from The Cancer Genome Atlas Colon Adenocarcinoma (TCGA-COAD) cohort [18]. The initial miRNA/mRNA expression matrices for tumor samples were obtained from the GDC Data Portal (https://portal.gdc.cancer.gov/) and converted into the RPKM/RPM tables using the procedure described above. Next, we calculated the Spearman correlation for each miRNA and the predicted target gene. The thresholds of 0.05 and –0.1 were applied to the FDR and correlation values, respectively.

Building the regulatory network of interactions between transcription factors and miRNA

We took the information on the regulatory interactions of transcription factors and miRNAs from the curated TransmiR v2.0 database [19]. The resulting interaction network was built and visualized in the yED Graph Editor (yWorks GmbH; Germany).  

RESULTS

Effect of Cobalt (II) Chloride on gene expression in HT-29 cells

To chemically induce hypoxia, HT-29 cells were treated with cobalt (II) chloride for 24 h. Analysis of sequencing of the RNA isolated from HT-29 control cells and treated with cobalt chloride showed a statistically significant change in the expression of 2511 genes that encode proteins, this change being 2-fold and greater in response to hypoxia. The search for overrepresented metabolic pathways identified 25 genes associated with the switch of aerobic metabolism to anaerobic glycolysis (KEGG pathway hsa00010 Glycolysis/Gluconeogenesis, FDR = 2.04 × 10–4; fig. 1A), and also demonstrated activation of the HIF-1 signaling pathway (KEGG pathway hsa04066 HIF-1 signaling pathway, FDR = 4.45 × 10–3) that allows cobalt chloride to simulate hypoxia (fig. 1B). In addition to the HIF-1 pathway, several other signaling cascades involved in the response to hypoxia were activated, including NF-κB [20] and AMPK [21] (attachment 1).

The most overrepresented category corresponded to the genes encoding proteins of proteasome complexes (KEGG pathway hsa03050 Proteasome, FDR = 2,02 × 10–16). In particular, 33 genes were significantly activated in response to hypoxia, including 6 out of 6 ATPases, 11 out of 12 subunits of the 26S proteasome lacking ATPase activity, 7 out of 8 α-subunits and 7 out of 11 β-subunits of the 20S proteasome, as well as the proteasome maturation protein POMP and proteasome activator subunit PSME4 (attachment 2). We have registered an increase in the expression of UBB, UBC, UBA52 and RPS27A genes, the increase being 3.1, 8.7, 2.0 and 1.6-fold, respectively. These genes encode ubiquitin, which is needed for proteasome-dependent protein degradation.

We have also detected a change in the expression of genes that encode proteins involved in focal adhesion (integrins and laminins). Accordingly, the expression levels of laminin subunits α3, β3, ϒ1, and ϒ2 were increased 5.5, 4.6, 3.1 and 4.5-fold, respectively. Three of them (α3, β3, ϒ2) can form a heterotrimer and thus generate laminin 332, also known as laminin-5 [22]. The expression of integrin subunits changed in a different direction: subunits αE, αV and β1 showed the growth of 1.8, 2.2, 1.8 times, while for subunits α1, α2, β3 and β8 the registered expression level decrease was 2.3, 2.1, 1.7 and 2.7 times.

Effect of hypoxia on the expression of miRNAs and their target genes

We have registered 16 miRNAs (tab. 1) showing a significant change in the expression as a response to the treatment of cells with cobalt chloride. One of them was hsa-miR-210- 3p, the only miRNA that is cited in all papers published as increasing expression in response to hypoxia [23]. The level of several confirmed hsa-miR-210-3p mRNA targets involved in mitochondrial metabolism and apoptosis induction decreased in response to hypoxia: GPD1L by 2.3 times, CASP8AP2 by 1.7 times, and AIFM3 by 8 times.
The following analysis was performed to assess the overall functional effects of the cellular miRNA content fluctuations induced by hypoxia. The TargetScan resource enabled compilation of the list of potential miRNA targets. Since some interactions between miRNA and target mRNA inhibit translation without affecting the level of mRNA expression, analyzing the samples from 441 colon adenocarcinoma patients (taken from the TCGA-COAD database, attachment 3) we searched for miRNA– target mRNA pairs showing a significant negative correlation. Next, we crossed the resulting list with the list of genes that had the expression changed at least two-fold in the direction opposite that peculiar to the corresponding miRNA. The result were six miRNAs with a statistically significant number of deregulated target genes (hypergeometric test; p < 0.05): hsa-miR-18a-5p, hsa-miR-22-3p, hsa-miR-27a-5p, hsa-miR -182-5p, hsa-miR-215-5p, hsa-miR-425-5p (attachment 4).

Building the network of regulatory interactions between transcription factors and miRNA

To better understand the mechanisms underlying the aberrant expression of microRNAs, we analyzed the possibility of regulation of microRNAs with transcription factors (TF). In particular, we considered TFs that had their mRNA representation significantly changing (two-fold or greater) in response to hypoxia and searched for the miRNAs they regulate in the TransmiR database of regulatory TF-microRNA interactions. As a result, we identified 30 TF–miRNA interactions between 15 TFs and 11 miRNAs. We have also considered the reciprocal miRNA-induced TF silence to build a complete regulatory network at these nodes (fig. 2). As is shown, the four TFs encoded by the EGR1, HIF1A, MYC, and RELA genes simultaneously regulate several miRNAs, while most miRNAs are regulated by the TF ensembles.

DISCUSSION

This study investigated changes in the transcriptome landscape of HT-29 cells in response to hypoxia induced by cobalt (II) chloride. The investigation relied on the mRNA and miRNA integrative sequencing technique. In addition to the activation of HIF-1, the canonical signaling pathway, we have shown that the expression of integrins and laminins changes, too. These proteins play a critical role in cell adhesion and interactions with the extracellular matrix. The latter is of particular importance, since recent data indicate a close relationship between the microenvironment formed by hypoxia and the metastatic progression of tumors, including colon adenocarcinoma [24]. One of the possible mechanisms for the metastatic spread of tumors is associated with laminin 332. By interacting with various receptors on the cell surface (including integrins α6β4 and α3β1, epidermal growth factor receptor and syndecan 1), as well as some other components of the basal membrane, laminin 332 regulates the process of oncogenesis, promotes invasion and survival of tumor cells [25]. The increased expression of the ϒ1-arm of laminins (encoded by the LAMC1 gene) may also play a role in the progression of tumor diseases, as has been shown in uterine carcinoma [26].

The analysis of profile of the small noncoding RNAs revealed several miRNAs expressed differentially in response to hypoxia. Some of these miRNAs were already reported as having their expression changed under hypoxia, including hsa-miR-210-3p [23], hsa-miR-27a-5p [27], hsa-miR-182-5p [28]. Four miRNAs (hsa-miR-30b-5p, hsa-miR-32-5p, hsa-miR-425-5p, hsa-miR- 1260a, and hsa-miR-1260b) were not previously covered in the published papers in connection with the cellular response to hypoxia. This can be explained both by the cellular specificity of the response and by the cobalt chloride's stream effects.

Particular attention should be paid to the miRNAs that show differential patterns of expression specific to a particular miRNA arm (arm-specific differential expression patterns). Namely, this applies to the miRNAs that had only the passenger arms hsa-let-7a, hsa-mir-10a and hsa-mir-27a regulated, while the expression of their guide arms did not change (tab. 2). We have recently reported a similar observation for miR-21- 3p (passenger arm) in mouse lungs; it showed an eight-fold increase in expression in response to the SARS-CoV infection, while the guide arm of the same miRNA increased only threefold [29]. One of the most promising theories that could explain this phenomenon is the regulation of miRNA arms by RNA-binding proteins [30].
The analysis of regulatory interactions between TF and miRNA showed that HIF-1, p65, c-Myc, and EGR1 (encoded by the HIF1A, RELA, MYC and EGR1 genes) are the key factors regulating transcription of differentially expressed miRNAs (see fig. 2). Three of the considered miRNAs demonstrated reciprocal capabilities and suppressed some of the TFs. In particular, HIF1A is a confirmed target for hsa-miR-18a-5p that is suppressed by multiple TFs.

CONCLUSIONS

The performed integrative sequencing of miRNA/mRNA allowed revealing significant changes in the transcriptome and miRNA profile in HT-29 cells in response to the CoCl₂-induced hypoxia. We have shown that differential expression of several of the miRNAs can cause significant changes in the expression of their target mRNAs. The analysis of the regulatory interactions between transcription factors and miRNAs revealed possible mechanisms underlying the observed response to hypoxia.