miRNA array analysis of plasma miRNA alterations in rats exposed to a high altitude hypoxic environment.

In the present study, the microRNA (miRNA) expression profiles of rats exposed to high altitude hypoxia and normal conditions were obtained from miRNA array analysis. Bioinformatics analyses, including the use of the Gene Oncology and Kyoto Encyclopedia of Genes and Genomes databases, were used to identify the genes and pathways, which were specifically associated with high altitude hypoxic environment‑associated miRNAs. A total of 26 miRNAs were differentially expressed in the two groups, comprising six upregulated and 20 downregulated miRNAs. In the present study, a novel pattern of upregulated miRNAs and their associated pathways were constructed, including proteoglycans in cancer, spliceosome, gluamatergic synapse, glycolysis/gluconeogenesis, Foxo, cGMP‑PKG and p53 signaling pathways, which may provide novel targets for diagnosing and understanding the mechanism of high altitude hypoxia‑induced disease.

Introduction

High altitudes create an environment with lower atmospheric pressure, compared with that at sea level. It is reported that there are ~17,000,000 individuals living 3,500 m above sea level worldwide and encounter hypoxia, which causes deficiency of red blood cells or hemoglobin, and can reduce the ability of the body to transfer oxygen to tissues (1–3). Furthermore, exposure to high altitude can cause digestive system disease, resulting in severe damage to the intestinal tract (4).

MicroRNAs (miRNAs), a set of endogenous small non-coding RNAs with a length of 19–24 nt, exert their functions through translational inhibition or the degradation of target mRNAs (5). It has been demonstrated that miRNAs are key in regulating genes, pathways and various biological networks (5). Circulating miRNAs are being intensively investigated for their involvement in various pathogenic processes and may serve as potential diseases biomarkers (6,7). The concentrations of circulating miRNA can be affected by factors, including age and gender, and environmental factors, including living conditions, residence and altitude (8). However, the expression profile and mechanism of miRNAs under hypoxia at high altitudes remain to be fully elucidated. The advent of genome-wide technologies, including gene expression microarrays, has made it possible to obtain a comprehensive view of miRNA alterations associated with high altitude hypoxia, and the use of bioinformatics enables analysis of the pathways of the altered miRNAs. A previous study (9) found that miRNAs were significantly altered in the Tibetan population, compared with the Nanjing Han population, which suggested that miRNAs may function as important regulators in high altitude hypoxic conditions.

The aim of the present study was to use microarray analysis to examine alterations in the plasma miRNAs of rats exposed to hypobaric hypoxia for different durations. In addition, miRNA functions and pathways were analyzed using Gene Ontology terms and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis. The results may assist in identifying novel targets to treat high altitude hypoxia-associated disease and provide biomarkers for diagnosis and prognosis.

Materials and methods

Animals and hypoxic exposure

Adult male Sprague-Dawley rats (n=24) weighing ~180–200 g were purchased from the Academy of Military Medical Sciences (Beijing, China). The animals were maintained in the animal house facility of the Key Laboratory of Disaster and Emergency Rescue Medicine in People's Liberation Army (Tianjin, China) under a standard environment (altitude 0 m, temperature 25±1°C, humidity 50±10%) with a 12 h light/dark cycle, and were provided with food and water for 7 days. To investigate the effects of acute hypobaric hypoxia exposure, the 24 male Sprague-Dawley rats were randomly divided into four groups (I to IV; n=6). Group I served as the normoxia group maintained in a standard environment. Groups II, III and IV served as hypoxia groups, in which the rats were exposed to simulated hypobaric hypoxia for 24, 48 and 72 h, respectively, at 25,000 ft (8,000 m; 7.5 mm Hg; 35.4 kPa) in a specially designed animal decompression chamber in which altitude was maintained by reducing the ambient barometric pressure, and temperature and humidity were precisely controlled. The airflow in the chamber was 2 l/min. The temperature and humidity were maintained at 20±2°C and 60±10%, respectively. The rate of ascent to altitude was maintained at a rate of 300 m/min and it took 25–30 min to reach the desired altitude. A slow rate of ascent in altitude with a gradual decrease in ambient pressure was used, as this is less likely to induce decompression-induced gas bubbles during exposure. Experiments and animal care were approved by the ethics committee of the Affiliated Hospital of the Logistics University of Chinese People's Armed Police Force (Tianjin, China).

Blood sample collection

Blood was collected from the inferior cava vein of rats in the normoxia (control) group, and the rats in groups I–IV immediately following 24, 48 and 72 h exposure to hypoxia. To prepare plasma, anticoagulants (EDTA, heparin or sodium citrate) were added to the blood samples immediately following collection of the blood to prevent clotting. EDTA plasma (10 ml) contained ~1.7 mg potassium EDTA; heparin plasma (5 ml) containing heparin 1 vial; sodium citrate plasma, (10 ml) containing 1 ml 0.118 mol/l (3.2%) citrate solution. The specimens were then centrifuged at 1,500 g for 10 min at 4°C to avoid hemolysis, and were then decanted and transferred into RNAase-free Eppendorf tubes as aliquots. The plasma samples were stored at −80°C until further analysis.

Arterial blood measurements

Arterial blood gases were determined using the i-STAT system with the CG8+ cartridge (Abbott Point of Care, Inc., Princeton, NJ, USA). The arterial blood samples (0.5 ml) were collected and the following blood parameters were determined: pH, partial pressure of arterial carbon dioxide (PCO2), partial pressure of arterial oxygen (PO2), total carbon dioxide (TCO2, base excess (BE), saturation of arterial blood oxygen (SaO2), bicarbonate concentration (HCO3−). The TCO2 was calculated as follows: TCO2 (mmol/l)=HCO3− + 0.03 PCO2.

Lung wet-to-dry weight (W/D) ratio

Upon termination of each experiment, the animals were sacrificed with a lethal dose of i.v. sodium pentobarbital. The lungs were isolated following chest opening, and the right superior lobe of the hemi-lung was excised and blotted on filter paper to remove adherent blood, with extra pulmonary tissue dissected. The wet weight was determined, followed by drying in an oven at 50°C for 72 h when a constant weight was achieved. The water content of the tissue, which was used as an index of pulmonary edema formation, was calculated as wet weight minus dry weight, expressed as mg water per mg dry tissue.

Lung histology

The remainder of right middle lobe lung was immersed in formalin, embedded in paraffin, cut into 6 µm sections and stained with hematoxylin and eosin for histological analysis. The severity of lung injury was inspected using light microscopy (Nikon Eclipse 50i; Nikon Instruments, Inc., Tokyo, Japan) (10).

RNA isolation

TRIzol was used for total RNA extraction according to the manufacturer's protocol (Invitrogen; Thermo Fisher Scientific, Inc., Waltham, MA, USA). Small RNAs (<200 nt) were separated from the total RNA using mirVana miRNA purification columns (Ambion; Thermo Fisher Scientific, Inc.) for microarray analysis and reverse transcription-quantitative polymerase chain reaction (RT-qPCR) analysis according to the manufacturer's protocol. The quality and quantity of each RNA preparation were determined using a Nanodrop ND-1000 spectrophotometer (Agilent Technologies, Inc., Santa Clara, CA, USA).

Expression of miRNAs

Each total RNA sample (700 ng) was labeled and hybridized using a FlashTag™ Biotin HSR RNA labeling kit (manufactured for Affymetrix, Inc., Santa Clara, CA, USA by Genisphere LLC, Hatfield, PA, USA). Total RNA was labeled using poly A polymerase. Biotin-labeled RNAs were hybridized for 16–18 h at 45°C on an Affymetrix miRNA v2.0 array. GeneChips were washed and stained in the Affymetrix Fluidics Station 450, and were then scanned using the Affymetrix GeneChip Scanner 3000 7G. The data were analyzed using the robust multi-array analysis-detection above background (RMA-DABG) normalization method, implemented in Affymetrix® expression console software (version 1.2.0.20; Affymetrix; Thermo Fisher Scientific, Inc.). The normalized and log-transformed intensity values were analyzed using Expression Console (Affymetrix; Thermo Fisher Scientific, Inc.). Fold change filters were set requiring genes to be present in ≥200% of the controls to be considered an upregulated miRNA and <50% of controls to be considered a downregulated miRNA.

RT-qPCR analysis

To validate the fold change results of the miRNA arrays, six miRNAs, comprising three upregulated (miR-25-5p, miR-451-5p and miR-466b-5p) and three downregulated (miR-214-3p, miR-140-3p and let-7a-5p) miRNAs, were randomly selected and examined using RT-qPCR analysis. The RT reaction was performed using mature miRNA-specific primer sets (Applied Biosystems; Thermo Fisher Scientific, Inc.) and an microRNA reverse transcription kit (Applied Biosystems; Thermo Fisher Scientific, Inc.) according to the manufacturer's protocol. PCR was performed using the 7500 Fast Real Time PCR system (Applied Biosystems; Thermo Fisher Scientific, Inc.) and the SYBR Premix Ex Taq™ kit (Takara Biotechnology Co., Ltd., Dalian, China) according to the manufacturer's protocol. The thermocycling conditions were as follows: 95°C for 10 min, followed by 35 cycles of 95°C for 15 sec and 60°C for 1 min. Melting curve analyses were subsequently performed. The fold change for each miRNA was calculated using the 2−ΔΔCq method (11) with U6 small nuclear RNA as the endogenous control. All reactions were performed in triplicate for each sample. The PCR primer sequences were as follows: The universal reverse primer was: CCAGTGCAGGGTCCGAGGT; and the forward primer sequences were as follows: miR-25-5p, TGCGGAGGCGGAGACUUGGG; miR-451-5p, TGCGGAAACCGUUACCAUUA; miR-466b-5p, TGCGGUGAUGUGUGUGUACA; miR-214-3p, TGCGGACAGCAGGCACAGAC; miR-140-3p, TGCGGUACCACAGGGUAGAA; let-7a-5p, TGCGGUGAGGUAGUAGGUUG; and U6, TGCGGGTGCTCGCTTCGGCAGC.

MicroRNA target prediction

The miRNA target sites were predicted by computer-aided algorithms obtained from TargetScan (version 7.1; www.targetscan.org/). The miRNA sequences were downloaded from the miRBase website (http://www.mirbase.org).

Gene ontology (GO) analysis

Based on the GO database (http://www.geneontology.org), the significant GO terms of the upregulated miRNA targeted genes were analyzed using the Database for Annotation, Visualization and Integrated Discovery (http://david.abcc.ncifcrf.gov/home.jsp) (12). Fisher's exact test and the Chi-square test were used to classify the GO categories, and the false discovery rate (FDR) was calculated to correct the P-values using the Benjamini and Hochberg procedure (13). The P-values of each differentially expressed gene in all the GO terms were calculated. P<0.05 was considered to indicate a statistically significant difference.

Pathway analysis

Pathway analysis was performed to determine the significant pathways of the differentially expressed genes, according to the KEGG Orthology Based Annotation System (http://kobas.cbi.pku.edu.cn/), in which pathway enrichment can be analyzed by the KEGG pathway, pathway interaction database, BioCyc, Reactome and Protein ANalysis THrough Evolutionary Relationships. Fisher's exact test and a χ2 test were used to select the significant pathway, and the threshold of significance was defined by the P-value and false discovery rate (13).

Statistical analysis

Numerical data are presented as the mean ± standard deviation. Differences between means were analyzed using Student's t-test. For comparison of multiple groups, one-way analysis of variance followed by Tukey's post-test was performed. All statistical analyses were performed using SPSS 13.0 software (SPSS, Inc., Chicago, IL, USA).

Results

Blood gas analysis

The results of the blood gases analysis confirmed that, compared with control group, PCO2 and lactate were markedly increased in the high altitude hypoxia group, whereas PO2, BE, HCO3−, TCO2 and SaO2 were significantly decreased at 24, 48 and 72 h in a time-dependent manner (P<0.01). The pH level showed marginal decrease with no statistically significant difference (Table I).

Table I.

Measurements of arterial blood gas parameters.

		High altitude hypoxia (h)
Parameter	Control	24	48	72
pH	7.396±0.03	7.369±0.03	7.291±0.04	7.263±0.03
PCO2 (mmHg)	34.4±5.01	37.280±5.23	39.475±2.46[a]	50.5±2.89[a–c]
PO2 (mmHg)	88.7±3.98	61.67±7.12[a]	53.33±2.16[a,b]	41.17±1.47[a–c]
BE (mmol/l)	1.67±0.52	−8±1.22[a]	−8.67±1.51[a]	−9.33±0.82[a]
HCO3− (mmol/l)	30.717±0.86	19.083±1.21[a]	18.825±1.33[a]	16.25±0.49[a]
TCO2 (mmHg)	32.67±1.51	20±1.41[a]	19.17±1.94[a]	17.5±1.58[a]
SaO2 (%)	91.67±1.86	67.5±1.8[a]	57±2.37[a,b]	47.8±2.48[a–c]
Lac (mmol/l)	0.982±0.36	1.488±0.72	1.537±0.74	1.635±0.81

P<0.01, compared with the control group

P<0.01, compared with the 24 h group

P<0.01, compared with the 48 h group (n=6). PCO2, partial pressure of arterial carbon dioxide; PO2, partial pressure of arterial oxygen; BE, base excess; HCO3−, bicarbonate; TCO2, total carbon dioxide; SaO2, saturation of arterial blood oxygen; Lac, lactate.

Lung W/D ratio and histology

The W/D lung weight ratio was measured to evaluate fluid accumulation in the lung specimens (14). Compared with the normal control groups, the W/D of lung tissues in the hypoxia groups were significantly upregulated at 24, 48 and 72 h (P<0.01). The W/D of the lung tissues at 48 and 72 h were higher, compared with that at 24 h. However, no significant difference in the lung W/D ratio was observed between the 48 and 72 h time points (Fig. 1A). Using light microscopy, a detailed examination of lung pathology was performed. Compared with the control groups, the pulmonary interstitium showed marginal dilations or hyperemia in the 24 h group. In the 48 h group, the pulmonary interstitium exhibited increased dilations and thickening. Blood stasis was also identified in the alveolar space. The tissues in the 72 h group exhibited marked thickening of the pulmonary interstitium and visible pink exudant in the alveolar space. Hyperemia and expansion of the alveolar capillary were also observed (Fig. 1B).

Figure 1.

Lung pathology. (A) Wet-to-dry lung weight ratios were measured in the normal control group and high-altitude hypoxia groups at 24, 48 and 72 h. *P<0.01, vs. 24, 48 and 72 h groups. (B) Light micrograph of hematoxylin and eosin-stained paraffin-embedded sections.

Overview of the miRNAs profiles

From the miRNAs expression profiles, differentially expressed miRNAs were identified between the high-altitude hypoxia environment and normal control samples. The miRNA expression profiles were determined by calculating the log FC in the high-altitude hypoxia group/normal group. Among a set of differentially expressed miRNAs, six miRNAs were significantly upregulated and 20 miRNAs were significantly downregulated in the high-altitude hypoxia group, compared with the control group. miR-451-5p, miR-25-3p, miR-879-3p, miR-486, miR-92a-3p and miR-466b-5p were the identified upregulated miRNAs, whereas miR-6216, miR-2985, miR-23a-3p, miR-6328, miR-22-3p and let-7a-5p were the top six significantly downregulated miRNAs (Fig. 2 and Table II).

Figure 2.

(A) Hierarchical clustering analysis of miRNA expression in rats exposed to high altitude hypoxic and normal conditions. miRNAs are presented in rows and samples are presented in columns. Colors indicate relative signal intensities, the color key value indicates the fold-change. (B) Column graph shows the fold changes of the differently expressed miRNAs. miRNA, microRNA.

Table II.

List of miRNAs with altered expression in high altitude hypoxia.

miRNA	Fold change
Upregulated
rno-miR-25-3p	6.796382
rno-miR-451-5p	6.811991
rno-miR-466b-5p	2.243331
rno-miR-486	3.692567
rno-miR-879-3p	3.723761
rno-miR-92a-3p	2.373798
Downregulated
rno-let-7a-5p	−7.25955
rno-let-7c-5p	−3.42339
rno-miR-125a-3p	−2.86652
rno-miR-140-3p	−5.06208
rno-miR-1896	−3.84228
rno-miR-208a-5p	−2.92098
rno-miR-214-3p	−4.88104
rno-miR-22-3p	−7.2768
rno-miR-223-3p	−6.82643
rno-miR-23a-3p	−7.45505
rno-miR-24-3p	−6.97121
rno-miR-27a-3p	−2.41551
rno-miR-2985	−8.11875
rno-miR-3573-3p	−5.68544
rno-miR-483-3p	−4.86634
rno-miR-494-3p	−2.27509
rno-miR-6216	−8.4358
rno-miR-6328	−7.31127
rno-miR-672-5p	−5.84322
rno-miR-762	−3.6636

miRNA/miR, microRNA.

Validation of the microarray analysis data using RT-qPCR analysis

The relative expression changes of the randomly selected miRNAs were analyzed using RT-qPCR analysis. The results of the miR-25-5p, miR-451-5p and miR-466b-5p upregulated miRNAs, and the miR-214-3p, miR-140-3p and let-7a-5p the downregulated miRNAs were generally consistent with the microarray analysis results, as shown in the histograms in Fig. 3 (P<0.05).

Figure 3.

Summary of the results of the reverse transcription-quantitative polymerase chain reaction analyses of upregulated miRNAs (miR-25-3p, miR-451-5p and miR-466b-5p) and downregulated miRNAs (miR-214-3p, miR-140-3p and let-7a-5p). The six randomly selected miRNAs in the control and hypoxia groups are shown with relative expression levels shown on the y-axis as upregulated and downregualated. *P<0.05, vs. control. miRNA/miR, microRNA.

Microarray-based GO analysis

For GO analysis, the target mRNAs for the six upregulated miRNAs were predicted using TargetScan (data not shown). As shown in Fig. 4A-D, GO analysis showed that the targets were involved in several biological processes, and the top five gene-associated processes were cellular process, single-organism process, metabolic process, biological process and regulation of biological process (Fig. 4A). Several genes were involved in the cellular component, the top five of which were cell, cell part, organelle, membrane and organelle part (Fig. 4B). Several genes were also involved in molecular function, the top five of which were binding, catalytic activity, molecular transducer activity, transporter activity and molecular function regulator (Fig. 4C and D). These results supported the hypothesis that these biological processes, cellular components and molecular functions are important in rats exposed to high altitude hypoxia.

Figure 4.

GO analysis of differentially expressed mRNAs in response to the high altitude hypoxic environment. The most related parts are shown for (A) biological process, (B) cellular component, (C) molecular function and (D) percentage of genes associated with the GO terms. GO, Gene Ontology. The node number indicates the corrected P-value (P<0.001). The node color indicates the Benjamini and Hochberg method corrected P-value for the enrichment of the term; the blue indicates biological process, yellow indicates cellular component, and red indicates molecular function, while the white color indicates little statistical enrichment.

Microarray-based pathway analysis

Pathway-based analysis assists in further understanding the biological functions of genes. In the present study, KEGG pathway enrichment analysis was performed to identify significantly enriched metabolic pathways or signal transduction pathways in the differentially expressed genes. Using pathway analysis, which considered the relative change direction and fold change, and had a the threshold of significance of P<0.05, the top 30 significant pathways were identified (Fig. 5). The pathways targeted by upregulated miRNAs with the highest enrichment included the proteoglycans in cancer, spliceosome, gluamatergic synapse, glycolysis/gluconeogenesis, Foxo, cGMP-PKG and p53 signaling pathways, which suggested that these pathways may be involved in the response to a high altitude hypoxic environment.

Figure 5.

Histogram of signaling pathways found to differ significantly between the high altitude hypoxia and normal environments. The x-axis shows the negative logarithm of the P-value (Log P), the y-axis shows the pathway. The higher the Log P, the lower the P-value.

Discussion

Understanding the clinical relevance of miRNA expression patterns in the environment during exposure to high altitude is necessary to circumvent the therapeutic challenges faced in clinical management. The present study used bioinformatics methods to screen the differently expressed miRNAs, and then analyzed the functions and pathways of the upregulated miRNAs. In the present study, a total of 26 aberrantly expressed miRNAs were identified in the blood of rats exposed to high altitude, compared to the samples from rats in normal conditions. Subsequently, RT-qPCR analysis was used to confirm the miRNA array results. To measure the concentrations of miRNAs, three upregulated and three downregulated miRNAs were selected for RT-qPCR analysis, the results of which were in accordance with the miRNA array assays. Among the differently expressed miRNAs, a subset of the altered miRNAs (miR-24, miR-25, miR-486, miR-451 and miR-92a) has also been associated with the hypoxic response in mammalian cells in other studies (15–19), which suggested that a high altitude hypoxic environment has a marked effect on human plasma miRNA patterns. The functional investigation of miRNAs, which respond during hypoxia, may assist in uncovering the molecular basis of hypoxic acclimatization and elucidate the complexity of hypoxic response pathways in humans. The present study also found that the expression of rno-miR-23a-3p was similar with the that reported in a study by Yan et al (9), which showed that human miR-23a-3p was downregulated in the Tibetan population, compared with the Nanjing Han population. This suggested the species conservation of this miRNA and indicated its importance for further investigation.

GO is widely recognized as a premier tool for molecular organization and functional annotation (20). Using the criteria of P<0.05 to identify significant GO terms, the present study revealed that the predicted targets of the six upregulated miRNAs were associated with biological processes, including cellular process, single-organism process, metabolic process, biological process and regulation of biological process; cellular components, including cell, cell part, organelle, membrane and organelle part; and molecular functions, including binding, catalytic activity, molecular transducer activity, transporter activity and molecular function regulator. The above GO terms have also been well represented in Triplophysa dalaica (21).

Pathway analysis can reveal distinct biological processes and identify the significant pathways, which dysregulated mRNAs are involved in. This enables a comprehensive understanding of the interactions of genes, their functions and the association between upstream and downstream genes, and can identify genes, which may be regulated by miRNAs. The appearance of the pathways in proteoglycans in cancer, spliceosome, gluamatergic synapse, glycolysis/gluconeogenesis, Foxo, cGMP-PKG and p53 signaling pathways confirmed their concordance with GO terms and their critical role in high altitude hypoxia. A significant adaptation was identified as an increased evolutionary rate and positive selection of genes involved in the hypoxic response and energy metabolism, analogous to those observed in other organisms in Tibet (22). It was previously reported that a conversion from oxidative glucose metabolism to glycolysis compensated for insufficient levels of oxygen in hypoxic conditions (23–25). The present study also showed the enhanced role of glycolysis. The FOXO subfamily of Forkhead transcription factors has a role in evolutionary conservation in cellular adaptation to stress stimuli, including hypoxic conditions (26). A previous study by Wang et al (27) showed that FOXO1 may be essential in adaptation to high altitudes. p53, as an upstream mediator of p21, may suppress hypoxic human lung fibroblast proliferation and pulmonary arterial remodeling by interacting with hypoxia-inducible factor-1 under hypoxia (28,29). A previous study showed that p53 gene deficiency with a decreased expression of p21 promoted hypoxia-induced pulmonary hypertension in mice (30). cGMP-dependent protein kinase (PKG) is a critical enzyme involved in the regulation of vascular contractility. Impaired PKG-mediated signaling has been found to be responsible for reduced cGMP-mediated pulmonary hypertension and vasodilatation following acute and chronic hypoxia (31–33). As these pathways have been identified to be involved in various hypoxic environments, together with the results of the present study, these pathways may be critical in high altitude hypoxic conditions.

In conclusion, the results of the present study identified six upregulated miRNAs and 20 downregulated miRNAs from two platforms. As upregulated miRNAs may better serve as biomarkers, the six upregulated miRNAs were used to perform GO and pathway analysis, which identified that Foxo, cGMP-PKG and p53 may be critical in the study model of hypoxia. Based on the integrated analysis of transcriptome features, these results may provide an important contribution to future investigations aimed at characterizing the role of specific miRNAs in the pathogenesis of high altitude hypoxia-induced diseases, and contribute to improving diagnosis and treatment.