Skip to main content
Download PDF

Comprehensive analysis of coding and non-coding RNA transcriptomes related to hypoxic adaptation in Tibetan chickens

Journal of Animal Science and Biotechnology volume 12, Article number: 60 (2021) Cite this article

Abstract

Background

Tibetan chickens, a unique native breed in the Qinghai-Tibet Plateau of China, possess a suite of adaptive features that enable them to tolerate the high-altitude hypoxic environment. Increasing evidence suggests that long non-coding RNAs (lncRNAs) and microRNAs (miRNAs) play roles in the hypoxic adaptation of high-altitude animals, although their exact involvement remains unclear.

Results

This study aimed to elucidate the global landscape of mRNAs, lncRNAs, and miRNAs using transcriptome sequencing to construct a regulatory network of competing endogenous RNAs (ceRNAs) and thus provide insights into the hypoxic adaptation of Tibetan chicken embryos. In total, 354 differentially expressed genes (DE genes), 389 differentially expressed lncRNAs (DE lncRNAs), and 73 differentially expressed miRNAs (DE miRNAs) were identified between Tibetan chickens (TC) and control Chahua chickens (CH). GO and KEGG enrichment analysis revealed that several important DE miRNAs and their target DE lncRNAs and DE genes are involved in angiogenesis (including blood vessel development and blood circulation) and energy metabolism (including glucose, carbohydrate, and lipid metabolism). The ceRNA network was then constructed with the predicted DE gene-DE miRNA-DE lncRNA interactions, which further revealed the regulatory roles of these differentially expressed RNAs during hypoxic adaptation of Tibetan chickens.

Conclusions

Analysis of transcriptomic data revealed several key candidate ceRNAs that may play high-priority roles in the hypoxic adaptation of Tibetan chickens by regulating angiogenesis and energy metabolism. These results provide insights into the molecular mechanisms of hypoxic adaptation regulatory networks from the perspective of coding and non-coding RNAs.

Background

MicroRNAs (miRNAs) are small non-coding RNA molecules that inhibit gene expression by binding to specific mRNAs. LncRNAs are another type of non-coding RNA, with a length of more than 200 nucleotides (nt), that regulate gene expression through a number of mechanisms, including epigenetic regulation (genetic imprinting and chromatin remodeling), transcriptional regulation (transcription interference), post-transcriptional regulation (splicing), and so on [1]. Both of these non-coding RNAs play different roles in various aspects of cellular function, including the regulation of hypoxia-related genes. Through their activity on target transcripts, some critical miRNAs regulate angiogenesis, the cell cycle, as well as lactate transport under hypoxic conditions [2,3,4]. HIF-1 and VEGF are regulated by miR-17-92, miR-20b, and miR-199a, which play roles in tumor development, cancer cell proliferation, and cardiac myocytes [5,6,7]. Various studies have shown that lncRNAs can directly bind to HIF-1 and activate its expression, in turn modulating hypoxic responses similar to the suppression of cancer cell adaptation to hypoxia, which provides important insights into how tumor cells sense and adapt to hypoxic stress [8,9,10]. Although there have been various studies on genes involved in high-altitude adaptation in humans and animals, the regulatory mechanism of non-coding RNAs involved in hypoxic adaptation remains largely unknown.

Tibetan chickens, which live in Tibetan highlands and have undergone long-term natural selection, exhibit stable genetic adaptation to high-altitude environments, characterized by high hatchability under hypoxic incubation when compared to lowland chickens [11, 12]. This indicates that Tibetan chickens have special physiological and genetic mechanisms enabling them to adapt to extreme high-elevation environments, and ensuring their normal hatching, reproduction, and growth in plateau environments [13,14,15]. Previous studies using genomic, transcriptomic, and proteomic profiles have identified several functional genes related to hypoxic adaptation in Tibetan chickens [16,17,18]. However, the mechanisms involved in the regulation of hypoxic adaptation by non-coding RNA are unclear.

Chorioallantoic membrane (CAM), which functions as both a respiratory and circulatory organ for the chicken embryo, contains numerous blood vessels. Tibetan chicken embryos have blunted responses to hypoxia on CAM angiogenesis, which might benefit blood flow and oxygen transportation during hypoxic incubation [18]. The characteristics of hypoxic adaptation may involve regulation of expressed mRNAs and non-coding RNAs; however, this remains unclear. In this study, full transcriptome sequencing (RNA-seq) of chicken embryos was performed to explore the regulation of non-coding RNAs and candidate genes responsible for hypoxic adaptation traits in Tibetan chickens as compared to the Chahua chickens. Integration of three types of RNAs (mRNAs, lncRNAs and miRNAs) revealed competing endogenous RNA (ceRNA) interactive network and potential regulatory mechanism to explain hypoxic adaptation in Tibetan chickens. These findings will help us further understand the molecular mechanism of adaptability to hypoxic conditions in Tibetan chickens from a non-coding RNA perspective.

Methods

Sample collection and preparation

The experiments and animal care protocol were approved by the animal welfare committee of the State Key Laboratory for Agro-Biotechnology of the China Agricultural University (approval number, XK257), and all methods were performed in accordance with the relevant guidelines and regulations.

We collected eggs from Tibetan chickens (TCs) and Chahua chickens (CHs) from the Experimental Station for Poultry Genetic Resources and Breeding, China Agricultural University (CAU). The two groups of eggs were incubated in a hypoxic incubator (13% ± 0.2% O2) and the temperature and humidity of the incubator were 37.8 °C and 60%, respectively. Eighteen CAM samples from Tibetan chickens (n = 9) and Chahua chickens (n = 9) were collected from female embryos at day 11 of incubation under hypoxic condition, immediately frozen in liquid nitrogen, and stored at − 80 °C for RNA extraction and sequencing.

RNA extraction, library construction, and sequencing

Total RNA was extracted from the CAM samples using an RNA pure Tissue Kit (Tiangen Biotech Co. Ltd., Beijing, China) according to the manufacturer’s protocol. Six RNA-seq libraries (TC1, TC2, TC3, CH1, CH2, and CH3) were constructed, three biological duplicates were prepared, and each sample consisted of a mix from three individuals of the same breed.

For mRNA and lncRNA sequencing, 1.5 μg of total RNA was used as input material for rRNA removal using the Ribo-Zero rRNA Removal Kit (Epicentre, Madison, WI, USA). Sequencing libraries were generated using the NEBNextR Ultra™ Directional RNA Library Prep Kit for IlluminaR (NEB, USA) following the manufacturer’s instructions, and 150 paired-end reads were generated. For small RNA sequencing, 1.5 μg of total RNA per sample was used to ligate the 3′ SR and 5′ SR Adaptor. Reverse transcription synthetic first chain and PCR amplification were carried out, PAGE gel was used to screen the electrophoresis fragment, then rubber cutting and recycling. The PCR products were purified (AMPure XP system) and library quality was assessed. Clustering of the index-coded samples was performed using TruSeq PE Cluster Kit v4-cBot-HS (Illumina) according to the manufacturer’s instructions, and 50-nt single-end reads were generated.

Read mapping and transcriptome assembly

For mRNA and lncRNA data, using in-house Perl scripts, raw 150 bp paired-end reads were filtered by removing adapters and reads containing poly-N and low-quality sequences to obtain clean reads. Clean reads were aligned to the chicken reference genome Gallus-gallus-6.0 (GRCg6a-96.Gallus-gallus. GRCg6a-96.genome.fa) using the HISAT2 software [19]. The transcriptome was assembled using StringTie and annotated using the gffcompare program. For miRNA data, the 3′ linker sequences were removed from the raw 50-nt single-end reads, which were then filtered by removing reads with ≥10% base N and shorter than 18 or longer than 30 nucleotides to obtain clean reads. The clean reads were aligned with the Silva, GtRNAdb, Rfam, and Repbase databases to filter other non-coding RNAs, such as ribosomal RNA (rRNA), transfer RNA (tRNA), small nuclear RNA (snRNA), small nucleolar RNA (snoRNA), and repeating sequences to obtain unannotated reads containing miRNAs. Unannotated reads were mapped to the reference genome (Gallus-gallus. GRCg6a-96) using the Bowtie software [20]. The Q20, Q30, and GC content of the clean data were calculated for all RNAs.

Identification of lncRNAs and prediction of their target genes

The unknown transcripts were screened as putative lncRNAs or protein-coding RNAs. The screening methods were used according to previous reports [21, 22], which briefly included length ≥ 200 bp, exon number ≥ 2, and FPKM ≥0.1, using tools of CPC2 (coding potential calculator), CNCI (coding-non-coding index), CPAT (coding potential assessment tool), and pfam protein domain analysis tools [23,24,25,26]. The transcripts were aligned with lncRNA databases (lncRNAdb and NONCODE) to identify known lncRNAs. Integration of following two methods to predict target genes of the lncRNAs: (1) lncRNA regulation of the expression of its neighboring genes, which is based on the distance (within 100 kb) between lncRNA and potential target genes in the genome, and (2) correlations between lncRNA and potential target gene expression in all samples that were significantly positive.

Identification of miRNAs and prediction of their target genes

Known miRNAs were identified by mapping the matched mature miRNAs in miRBase (v22) to the reference genome. Novel miRNAs were predicted using miRDeep2 from the unannotated reads [27]. Target genes of miRNAs were predicted by miRanda and TargetScan based on the miRNA and gene sequence information [28, 29].

RNA quantification and differential expression analysis

Expression level of mRNAs and lncRNAs were quantified based on fragments per kilobase of transcript per million fragments mapped (FPKM) using the StringTie software [30]. Expression levels of miRNAs were quantified and normalized using the TPM algorithm (TPM = Readcount × 1,000,000/mapped reads) [31]. The DESeq2 software package was used to detect differentially expressed lncRNAs (DE lncRNAs), differentially expressed genes (DE genes), and differentially expressed miRNAs (DE miRNAs) with criteria of log2| fold change| > 1 and P < 0.05 [32].

GO and KEGG enrichment analysis of differentially expressed RNA

Gene Ontology (GO) enrichment analysis of the DE genes was performed using the GOseq R package based on Wallenius’ non-central hypergeometric distribution. Kyoto Encyclopedia of Genes and Genomes (KEGG) is a database resource for understanding high-level gene functions and utilities of a biological system, such as the cell, organism, and ecosystem, using molecular-level information, especially large-scale molecular datasets generated by genome sequencing and other high-throughput experimental technologies (http://www.genome.jp/kegg/). We used the KOBAS software to test the statistical significance of DE gene enrichment in KEGG pathways.

Construction of the ceRNA network

Competitive endogenous RNAs (ceRNAs) are novel transcriptional regulators that can modulate the expression of transcripts through competitive binding. The interactions of miRNA-mRNA and miRNA-lncRNA was predicted using TargetScan and miRanda, and the competitive pairs of mRNA-lncRNA were identified using hypergeometric test and Pearson correlation coefficient. Based on the sharing of miRNA binding sites and competitive pairs between lncRNAs and mRNAs, ceRNA network was constructed as follows: P < 0.01 and FDR < 0.01 in hypergeometric test [33,34,35]. Visualization of ceRNA using Cytoscape software (http://cytoscape.org/).

Expression analysis of DE lncRNAs, DE miRNAs, and DE genes using qRT-PCR

To confirm the differentially expressed RNAs identified through RNA-seq, quantitative real-time PCR (qRT-PCR) was carried out for six DE lncRNAs, six DE miRNAs, and seven DE genes. We used the PrimeScript RT reagent kit with gDNA (Tiangen Biotech Co. Ltd., Beijing, China) to convert total RNA to cDNA, using random hexamers for mRNA and lncRNAs, and stem-loop RT primers for miRNAs. The qPCR was performed using the SYBR Green PCR Kit (Tiangen Biotech Co. Ltd., Beijing, China) according to the manufacturer’s instructions. The internal controls included GAPDH (glyceraldehyde-3-phosphate dehydrogenase) for lncRNAs, chicken 5S rRNA for miRNAs, and HPRT (hypoxanthine guanine phosphoribosyl transferase) for mRNAs. All primers used in qRT-PCR are shown in Additional file 1: Table S1. The BioRad CFX96 (Bio-Rad, CA, USA) was used to perform qRT-PCR with the SYBR Green PCR Master Mix (Tiangen Biotech Co. Ltd.). Each qRT-PCR experiment was performed in triplicate, relative RNA expression levels were calculated using the 2−ΔΔCT method, and ANOVA testing via SPSS 25.0 was applied. The results are expressed as mean ± SE, and P ≤ 0.05 was considered statistically significant. Software and associated parameters used are shown in Additional file 2: Table S2.

Results

Overview of RNA transcriptomic profiles of chicken CAM tissues

After low-quality reads were removed, 110.3–114.1 M paired-end reads (2 × 150 bp in length) for the mRNAs and lncRNAs, and 23.8–29.2 M single-end reads (50 bp in length) for the miRNAs, were obtained from the CAM of Tibetan and Chahua chicken samples (Additional file 3: Table S3). The majority of the lncRNAs (60.7%) were intergenic long-stranded non-coding RNA (lincRNA), followed by intronic-lncRNA (17.7%), antisense-lncRNA (12.4%), and 9.2% sense-lncRNA (Fig. 1a). The length of mRNAs ranged from 51 to 30,0721 bp, while that of lncRNAs ranged from 202 to 223,626 bp; the majority of mRNAs and lncRNAs were 201–400 bp and > 3000 bp long, respectively (Additional file 4: Fig. S1). In total, 20,912 expressed genes, 15,523 lncRNAs (3194 known), and 1439 miRNAs were identified in the CAM tissues of Tibetan and Chahua chickens. The 1439 miRNAs included 471 known and 968 novel miRNAs, and most of the miRNA reads were between 20 and 24 nt in length (Fig. 1b, Additional file 5: Table S4). The values of FPKM of mRNA/lncRNA and TPM of miRNA distribution were similar among the six samples (Fig. 1c, d, e).

Fig. 1
figure 1

Comparison of the features of RNAs in the chorioallantoic membrane of chicken embryos. (a) Statistics of long non-coding RNA, including intergenic non-coding RNAs (lincRNAs), anti-sense lncRNAs, intronic lncRNAs, and sense lncRNAs; (b) Length distribution of 1439 miRNAs; (c) Distribution of mRNA expression values (FPKM) of six samples; (d) Distribution of lncRNA expression values (FPKM) of six samples; (e) Distribution of miRNA expression values (TPM) of six samples. The ordinate shows the log10(FPKM + 1) (c, d) and the log10(TPM + 1) (e) of each RNA from the six sets of RNA-seq data; the middle line in the box represents the median of FPKM/TPM. TC, Tibetan chicken; CH, Chahua chicken

Full size image

Differential expression analysis of RNAs between TC and CH

We found 389 DE lncRNAs, 73 DE miRNAs, and 354 DE genes between the TC and CH groups, among which 194 DE lncRNAs, 46 DE miRNAs, and 180 DE genes were upregulated, whereas 195 DE lncRNAs, 27 DE miRNAs, and 174 DE genes were downregulated in the Tibetan chickens compared to Chahua chickens (Fig. 2, Additional file 6: Table S5). The DE lncRNAs, DE miRNAs, and DE genes were classified into two categories for the six samples, in which the expressed quantity of the DE lncRNAs, DE miRNAs, and DE genes had good repeatability within groups respectively (Additional file 4: Fig. S2).

Fig. 2
figure 2

Identification of differentially expressed (DE) lncRNAs, DE miRNAs, and DE genes between Tibetan (TC) and Chahua chicken (CH). Number of DE lncRNAs, DE miRNAs, and DE genes (a); Volcano plots displaying DE lncRNAs (b), DE miRNAs (c), and DE genes (d) between Tibetan and Chahua chicken. Upregulated and downregulated genes are shown in red and green, respectively. Black dots represent genes with similar expression levels

Full size image

Functional analysis of differentially expressed non-coding RNAs (DE lncRNAs and DE miRNAs)

For 389 DE lncRNAs, we predicted 722 cis- and trans-target genes that were mainly involved in the GO terms of angiogenesis, blood vessel development, regulation of vasoconstriction, regulation of blood pressure, response to oxygen levels, ATP metabolic process, fatty acid metabolic process, glucuronate metabolic process, enzyme binding, and ATP binding (Fig. 3a, Additional file 7: Table S6a). The KEGG analysis revealed the involvement of target genes in carbon metabolism, fatty acid metabolism, vascular smooth muscle contraction, and the Wnt, mTOR, MAPK, VEGF, and calcium signaling pathways (Fig. 3b, Additional file 7: Table S6b).

Fig. 3
figure 3

GO and KEGG enrichment analysis of differentially expressed (DE) lncRNAs. (a) GO enrichment analysis of DE lncRNAs between Tibetan (TCs) and Chahua chickens (CHs). KS represents the statistical significance of enriching the GO term, the smaller the KS value mean the more significant the enrichment; (b) KEGG enrichment analysis of DE lncRNAs between Tibetan and Chahua chickens

Full size image

A total of 1766 target genes were predicted for 73 DE miRNAs. The target genes were mainly enriched in the GO terms of vasculogenesis, positive regulation of fatty acid oxidation, glycogen metabolic process, angiogenesis, vasculature development, response to hypoxia, oxygen metabolic process, vasodilation, ATP binding, and enzyme binding (Fig. 4a, Additional file 8: Table S7a). Notably, the enriched KEGG pathways were cysteine and methionine metabolism, amino sugar and nucleotide sugar metabolism, apoptosis, vascular smooth muscle contraction, fatty acid metabolism, carbon metabolism, vitamin B6 metabolism, and the MAPK, Notch, VEGF, mTOR, and calcium signaling pathways (Fig. 4b, Additional file 8: Table S7b).

Fig. 4
figure 4

GO and KEGG enrichment analysis of differentially expressed (DE) miRNAs. (a) GO enrichment analysis of DE miRNAs between Tibetan (TCs) and Chahua chickens (CHs). The enrichment factor is the ratio of DE miRNA numbers annotated in the GO term to the total gene numbers annotated in the same GO term. (b) KEGG enrichment analysis of DE miRNAs between Tibetan and Chahua chickens

Full size image

Comparative parsing of lncRNAs and mRNAs and functional analysis of DE genes

To compare the exon number, open reading frame (ORF), and expression levels of lncRNA and mRNA, we needed to consider the differences in structure and sequence. The number of exons corresponding to lncRNAs was mainly less than 8, whereas the number of exons in mRNAs was relatively large, ranging from 1 to > 30 (Additional file 4: Fig. S3a, b). The majority of lncRNA ORFs were mainly between 50 and 250 nt long, while the length of the ORFs of mRNAs was mainly between 100 and 1,000 nt (Additional file 4: Fig. S3c, d). Expressed mRNAs and lncRNAs were located on all chromosomes of the chicken genome (Additional file 4: Fig. S4). The average and maximum values of mRNA expression were higher compared to those of lncRNA expression (Fig. 5a).

Fig. 5
figure 5

Comparative parsing of lncRNAs and mRNAs and functional analysis of differentially expressed (DE) genes. (a) A comparison of expression between lncRNAs and mRNAs. The box chart statistics from top to bottom are the maximum value, upper quartile, middle value, lower quartile, and the minimum value; (b) GO and KEGG enrichment analysis of DE genes between Tibetan (TC) and Chahua chicken (CH)

Full size image

The 354 DE genes were mainly enriched in GO terms of defense response, homeostatic process, regulation of lipid metabolic process, glucose homeostasis, blood microparticle, oxygen transport, blood vessel endothelial cell migration, blood circulation, and respiratory system development. Enriched KEGG pathways included fatty acid metabolism, vascular smooth muscle contraction, as well as the calcium, PPAR and MAPK signaling pathways (Fig. 5b, Additional file 9: Table S8). Twelve DE genes (ACTC1, KCNMB4, NCS1, NGFR, ADAM8, CASQ2, IRF4, PTPRZ1, CALML3, ERBB4, ARR3, and NTSR1) were mainly associated with angiogenesis, blood circulation, hematopoiesis, as well as calcium and MAPK signaling pathways. Fifteen DE genes (SSTR5, NR1H4, HTR2C, APOA1, KCNB1, P3H2, CHST8, LYZ, HAO2, ACER1, ACSBG1, ELOVL2, ELOVL3, GBE, and NOX3) were related to glucose and carbohydrate metabolism, fatty acid metabolism, and lactate oxidation (Table 1).

Table 1 Potential key differentially expressed mRNAs (DE genes), their targeted differentially expressed miRNAs (DE miRNAs)/differentially expressed lncRNAs (DE lncRNAs), and their functions related to hypoxic adaption in the Tibetan chicken
Full size table

Construction of the ceRNA network

We constructed a ceRNA network that yielded 529 pairs of candidate ceRNAs (lncRNA-miRNA-mRNA) from 162 DE lncRNAs and 108 DE genes through 25 miRNA target-mediated relationships (Additional file 4: Fig. S5). Considering that the network contains enormous information, and each relationship cannot be displayed in the figure, we constructed a mini-ceRNA (DE lncRNA-DE miRNA-DE mRNA) network of important RNAs. We identified 10 known mRNAs that were related to hypoxic adaptation and involved in angiogenesis and blood circulation (NGFR, ACTC1, CASQ2, ERBB4, KCNMB4, NCS1, and NTSR1) (Fig. 6a), as well as energy metabolism (SSTR5, NR1H4, and GBE) (Fig. 6b). A total of 39 lncRNA-miRNA-mRNA interactions were identified in the mini-ceRNA network constructed with these 10 DE genes and their targeted 37 DE lncRNAs and 9 DE miRNAs (Table 1, Additional file 10: Table S9).

Fig. 6
figure 6

Mini-ceRNA network constructed with 10 differentially expressed (DE) genes and their 37 target DE lncRNAs and 9 target DE miRNAs, involved in angiogenesis (a) and energy metabolism (b). Red color represents the up-regulated and blue color represents down-regulated levels

Full size image

Expression levels of DE lncRNAs, DE miRNAs, and DE genes using qRT-PCR

The expression levels of six DE lncRNAs (MSTRG.118362.24, MSTRG.25780.4, MSTRG.19949.13, MSTRG.29590.15, MSTRG.32082.2, and MSTRG.73129.28), six DE miRNAs (gga-miR-726-3p, gga-miR-6608-3p, gga-miR-6606-5p, novel-miR-676, gga-miR-460b-3p, and gga-miR-7442-5p), and seven DE genes (ARR3, NGFR, NTSR1, ADAM8, ACTC1, CASQ2, and CALML3) were measured using qRT-PCR to validate expression differences identified through RNA-seq. The results indicated that the expression levels of all RNAs were significantly different between TC and CH (P < 0.05), and fold-changes in expression followed the same trend in qRT-PCR and RNA-seq (Fig. 7).

Fig. 7
figure 7

The expression of six differentially expressed (DE) lncRNAs (a), six DE miRNAs (c) and seven DE genes (e) validated with qRT-PCR in chicken embryo chorioallantoic membrane (CAM) tissue. Note: The Y-axis represents the expression value of RNA in chicken embryo CAM, and the X-axis represents the types of RNAs. Error bars represent SE of expression. * on the bars indicate P < 0.05 and ** indicate P < 0.01 between TC and CH breeds. Tibetan chickens (TCs, n = 4–8), and Chahua chickens (CHs, n = 4–8). The fold changes of the six DE lncRNAs (b), six DE miRNAs (d) and seven DE genes (f) showed that qRT-PCR results were consistent with the RNA-seq data

Full size image

Discussion

Through evolution, Tibetan chickens have developed genetic adaptability to survive in low-oxygen environments. Previous reports have shown that Tibetan chickens exhibit unique phenotypic and physiological characteristics, including enhanced blood oxygen binding capacity and increased blood circulation [12, 17]. The ceRNA hypothesis was reported a few years ago and has now been widely accepted, which has improved our understanding of human diseases, as well as led to important advances in animal research [36, 37]. The recently discovered non-coding RNAs represent a new type of master regulator that affects gene expression and modulates a variety of cellular processes. Based on these theories, we constructed a ceRNA (lncRNA-miRNA-mRNA) regulatory network to clarify the molecular mechanism of Tibetan chicken embryo CAM adaptation to hypoxic environments.

The CAM, an important respiratory and circulatory organ for chicken embryo development, contains a large number of blood vessels, and hypoxia can induce increased blood vessel density [18, 38]. From the constructed ceRNA (lncRNA-miRNA-mRNA) network, we identified several miRNAs, including down-regulated novel-miR-589, novel-miR-815, novel-miR-85, novel-miR-669, as well as up-regulated miR-6606-5p, novel-miR-676, novel-miR-589, novel-miR-815, novel-miR-85, and novel-miR-567. Among these, miR-6606-5p was significantly up-regulated in the Tibetan chicken, and its target differentially expressed genes (ACTC1, NCS1, and KCNMB4) were found to be involved in vascular development under hypoxic conditions. ACTC1 and NCS1, up-regulated in Tibetan chickens, were enriched in the GO term of blood circulation and regulated calcium ion transport;Ca2+ is necessary for the activity of HIF-1, which is a major transcriptional regulator of cells and development in response to hypoxia [39, 40]. KCNMB4 has been reported to regulate blood pressure and was enriched in the GO terms of regulation of vasoconstriction, ion transport, and vascular smooth muscle contraction pathway [41]. Consequently, these findings indicate that Tibetan chickens can improve blood circulation and stimulate the activity of HIF-1 through upregulation of ACTC1 and NCS1, and Tibetan chickens promote vasodilation through downregulation of KCNMB4 in CAM to adapt to hypoxic conditions. It is worth mentioning that with the exception of miR-6606-5p, other differentially expressed miRNAs, including miRNA-155, miRNA-302a, miR-302b-3p, miR-302b-5p, miR-302c-3p, miR-460b-3p, and miR-460b-5p, also play important roles in angiogenesis under hypoxia. The miRNA-155 promotes HIF-1α activity during prolonged hypoxia and participates in the PI3K/AKT pathway [42, 43]. Previous studies have shown that the miR-302 family (miR-302a/b/c/d) suppresses the proliferation, migration, and angiogenesis of vascular endothelial cells by targeting VEGFA [44,45,46]. The miR-460b-3p and miR-460b-5p were identified in animal models of hypoxic pulmonary hypertension and were involved in the regulation of HIF-1α [47, 48]. Moreover, other differentially expressed genes (NGFR, ADAM8, CASQ2, IRF4, PTPRZ1, CALML3, ERBB4, ARR3, and NTSR1) were enriched in blood vessel development, angiogenesis, blood circulation, hematopoiesis, response to hypoxia, oxygen transport, vascular smooth muscle contraction, calcium, and MAPK signaling pathways. The regulatory mechanism of these genes in hypoxic adaptation of Tibetan chickens is the same as that of ACTC1, NCS1 and KCNMB4 genes. In general, they regulate angiogenesis through the up or down-regulation of genes in their pathways. For example, NGFR is the receptor of NGF, which can induce chick CAM neovascularization [49]. Studies have shown that ADAM8 is significantly induced by hypoxia and plays a role in the proliferation and migration of endothelial cells during angiogenesis [50, 51]. The explicit regulatory mechanism of these genes in hypoxic adaption of Tibetan chickens needs to be further studied. Therefore, we speculate that the enhanced tolerance of Tibetan chicken embryo CAM under hypoxic conditions is ascribed to miRNA-mediated modulation of the related target mRNAs, further regulating HIF, which enables Tibetan chicken to maintain hypoxic adaptation via the promotion of angiogenesis and blood circulation.

Glucose uptake and carbohydrate metabolism are the basis for the maintenance of normal physiological functions in humans and animals. Under hypoxic conditions, oxygen and carbon dioxide metabolism mainly depend on mitochondrial respiration and make use of adenosine triphosphate (ATP). Recently, a number of studies have shown the key regulatory role of energy metabolism, including glucose, carbohydrate, and lipid metabolic processes, among others, during hypoxic adaptation [52,53,54,55]. In the current research, from the constructed ceRNA network, we identified seven differentially expressed miRNAs (novel-miR-819, novel-miR-676, novel-miR-85, novel-miR-693, novel-miR-775, novel-miR-669, and novel-miR-867) and sixteen differentially expressed genes related to energy metabolism including glycosylation, glucose metabolism process, carbohydrate metabolism process, fatty acid metabolism, and ATP binding. For example, SSTR5 and NR1H4 are essential for glucose homeostasis and play a pivotal role in glucose metabolism in animals [56, 57]. Cortisol regulates the metabolism of mouse adipose cells through the serotonin receptor gene HTR2C, and genetic variation of the APOA1 gene is linked to lipid metabolism and cardiovascular disease risk [58, 59]. It should be considered that more differentially expressed genes (KCNB1, P3H2, CHST8, LYZ, HAO2, ACER1, ACSBG1, ELOVL2, ELOVL3, GBE, and NOX3) are targeted by key miRNAs, and several studies have shown that these mRNAs are also involved in energy metabolism. Therefore, we speculate that Tibetan chickens have enhanced energy metabolism due to the function of these RNAs, allowing for adaptation to hypoxic conditions.

Accumulating evidence indicates that lncRNAs play the role of ceRNAs (or miRNA sponges) in a variety of biological processes, including high-altitude adaptation. Such as polymorphisms of LINC-PINT and LINC00599 are associated with high-altitude pulmonary edema in Chinese populations [60]. Another study reported the expression profiles of lncRNAs in mice with high-altitude hypoxia-induced brain injury and provided new insights into the molecular mechanism of its treatment [61]. In addition, a previous study reported the expression profiles of lncRNAs responsible for fatness and fatty acid composition traits in Tibetan pigs [62]. In the current work, we found that the aforementioned differentially expressed genes involved in angiogenesis and energy metabolism were targeted by 37 differentially expressed lncRNAs in the ceRNA network, suggesting that these lncRNAs may also function as miRNA sponges and may play a role in the hypoxic adaptation of Tibetan chicken embryos with regard to angiogenesis and energy metabolism.

Based on the ceRNA theory and the ceRNA network constructed in this study, we propose a mode of action of differentially expressed lncRNAs, miRNAs, and mRNAs during hypoxic adaptation of Tibetan chicken embryos (Fig. 8). Under hypoxic conditions, miRNAs act as key regulators to modulate the up or down-regulation of important differentially expressed genes. As a consequence, the angiogenesis/blood circulation of chorioallantoic capillaries and energy metabolism, such as glucose/carbohydrate metabolism, are stimulated, leading to enhanced hypoxia adaptability of Tibetan chicken embryos. In this process, some lncRNAs act as ceRNAs to competitively bind the miRNA response element of miRNAs, which may indirectly affect the expression of mRNA.

Fig. 8
figure 8

Proposed model of hypoxic adaptation in Tibetan chicken embryos. Under hypoxic conditions, miRNAs act as key regulators that directly target important differentially expressed (DE) genes. In this process, some lncRNAs can act as ceRNAs to competitively bind the miRNA response element (MRE) of miRNAs, indirectly affecting the expression of mRNAs

Full size image

Conclusions

In conclusion, 389 DE lncRNAs, 73 DE miRNAs, and 354 DE genes were identified between Tibetan and Chahua chickens. The transcriptomic data revealed several key candidate ceRNAs (DE lncRNAs-DE miRNAs-DE genes) that may play high-priority roles in the hypoxic adaptation of Tibetan chickens by regulating angiogenesis and energy metabolism. These results provide insights into the molecular mechanisms of hypoxic adaptation regulatory networks from the perspective of coding and non-coding RNAs.

Availability of data and materials

All RNA-seq data were deposited in the Gene Expression Omnibus under accession number GSE160324.

Abbreviations

lncRNAs:

Long non-coding RNA

miRNA:

MicroRNA

TC:

Tibetan chicken

CH:

Chahua chicken

CAM:

Chorioallantoic membrane

DE:

Differentially expressed gene

DE lncRNA:

Differentially expressed lncRNA

DE miRNA:

Differentially expressed miRNA

ceRNA:

Competing endogenous RNA

GO:

Gene Ontology

KEGG:

Kyoto Encyclopedia of Genes and Genomes

qRT-PCR:

Quantitative real-time PCR

ORF:

Open reading frame

References

  1. Sun W, Yang Y, Xu C, Guo J. Regulatory mechanisms of long noncoding RNAs on gene expression in cancers. Cancer Genet. 2017;216-217:105–10. https://0-doi-org.brum.beds.ac.uk/10.1016/j.cancergen.2017.06.003.

    Article  CAS  PubMed  Google Scholar 

  2. Kochan-Jamrozy K, Kroliczewski J, Moszynska A, Collawn JF, Bartoszewski R. miRNA networks modulate human endothelial cell adaptation to cyclic hypoxia. Cell Signal. 2019;54:150–60. https://0-doi-org.brum.beds.ac.uk/10.1016/j.cellsig.2018.11.020.

    Article  CAS  PubMed  Google Scholar 

  3. Gupta A, Ragumani S, Sharma YK, Ahmad Y, Khurana P. Analysis of Hypoxiamir-gene regulatory network identifies critical miRNAs influencing cell-cycle regulation under hypoxic conditions. Microrna. 2019;8(3):223–36. https://0-doi-org.brum.beds.ac.uk/10.2174/2211536608666190219094204.

    Article  CAS  PubMed  Google Scholar 

  4. Zhao LL, Wu H, Sun JL, Liao L, Cui C, Liu Q, et al. MicroRNA-124 regulates lactate transportation in the muscle of largemouth bass (micropterus salmoides) under hypoxia by targeting MCT1. Aquat Toxicol. 2020;218:105359. https://0-doi-org.brum.beds.ac.uk/10.1016/j.aquatox.2019.105359.

    Article  CAS  PubMed  Google Scholar 

  5. Taguchi A, Yanagisawa K, Tanaka M, Cao K, Matsuyama Y, Goto H, Takahashi T. Identification of hypoxia-inducible factor-1 alpha as a novel target for miR-17-92 microRNA cluster. Cancer Res. 2008;68(14):5540–5. https://0-doi-org.brum.beds.ac.uk/10.1158/0008-5472.CAN-07-6460.

    Article  CAS  PubMed  Google Scholar 

  6. Lei Z, Li B, Yang Z, Fang H, Zhang GM, Feng ZH, Huang B. Regulation of HIF-1alpha and VEGF by miR-20b tunes tumor cells to adapt to the alteration of oxygen concentration. PLoS One. 2009;4(10):e7629. https://0-doi-org.brum.beds.ac.uk/10.1371/journal.pone.0007629.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  7. Rane S, He M, Sayed D, Vashistha H, Malhotra A, Sadoshima J, et al. Downregulation of miR-199a derepresses hypoxia-inducible factor-1alpha and Sirtuin 1 and recapitulates hypoxia preconditioning in cardiac myocytes. Cir Res. 2009;104(7):879–86. https://0-doi-org.brum.beds.ac.uk/10.1161/CIRCRESAHA.108.193102.

    Article  CAS  Google Scholar 

  8. Huang L, Wang W, Hu Z, Guan C, Li W, Jiang X. Hypoxia and lncRNAs in gastrointestinal cancers. Pathol Res Pract. 2019;215(12):152687. https://0-doi-org.brum.beds.ac.uk/10.1016/j.prp.2019.152687.

    Article  CAS  PubMed  Google Scholar 

  9. Shih JW, Chiang WF, Wu ATH, Wu MH, Wang LY, Yu YL, et al. Long noncoding RNA LncHIFCAR/MIR31HG is a HIF-1alpha co-activator driving oral cancer progression. Nat Commun. 2017;8(1):15874. https://0-doi-org.brum.beds.ac.uk/10.1038/ncomms15874.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Wang X, Wang Y, Li L, Xue X, Xie H, Shi H, et al. A lncRNA coordinates with Ezh2 to inhibit HIF-1alpha transcription and suppress cancer cell adaption to hypoxia. Oncogene. 2020;39(9):1860–74. https://0-doi-org.brum.beds.ac.uk/10.1038/s41388-019-1123-9.

    Article  CAS  PubMed  Google Scholar 

  11. Bao HG, Zhao CJ, Li JY, Zhang H, Wu C. A comparison of mitochondrial respiratory function of Tibet chicken and silky chicken embryonic brain. Poult Sci. 2007;86(10):2210–5. https://0-doi-org.brum.beds.ac.uk/10.1093/ps/86.10.2210.

    Article  CAS  PubMed  Google Scholar 

  12. Gou X, Li N, Lian L, Yan D, Zhang H, Wei Z, et al. Hypoxic adaptations of hemoglobin in Tibetan chick embryo: high oxygen-affinity mutation and selective expression. Comp Biochem Phys B. 2007;147(2):147–55. https://0-doi-org.brum.beds.ac.uk/10.1016/j.cbpb.2006.11.031.

    Article  CAS  Google Scholar 

  13. Zhang H, Wu CX, Chamba Y, Ling Y. Blood characteristics for high altitude adaptation in Tibetan chickens. Poult Sci. 2007;86(7):1384–9. https://0-doi-org.brum.beds.ac.uk/10.1093/ps/86.7.1384.

    Article  CAS  PubMed  Google Scholar 

  14. Zhang H, Wang XT, Chamba Y, Ling Y, Wu CX. Influences of hypoxia on hatching performance in chickens with different genetic adaptation to high altitude. Poult Sci. 2008;87(10):2112–6. https://0-doi-org.brum.beds.ac.uk/10.3382/ps.2008-00122.

    Article  CAS  PubMed  Google Scholar 

  15. Li M, Zhao C. Study on Tibetan chicken embryonic adaptability to chronic hypoxia by revealing differential gene expression in heart tissue. Sci China Ser C. 2009;52(3):284–95. https://0-doi-org.brum.beds.ac.uk/10.1007/s11427-009-0005-8.

    Article  CAS  Google Scholar 

  16. Zhang Q, Gou W, Wang X, Zhang Y, Ma J, Zhang H, et al. Genome Resequencing identifies unique adaptations of Tibetan chickens to hypoxia and high-dose ultraviolet radiation in high-altitude environments. Genome Biol Evol. 2016;8(3):765–76. https://0-doi-org.brum.beds.ac.uk/10.1093/gbe/evw032.

    Article  PubMed  PubMed Central  Google Scholar 

  17. Zhang Y, Gou W, Ma J, Zhang H, Zhang Y, Zhang H. Genome methylation and regulatory functions for hypoxic adaptation in Tibetan chicken embryos. PeerJ. 2017;5:e3891. https://0-doi-org.brum.beds.ac.uk/10.7717/peerj.3891.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Zhang Y, Zheng X, Zhang Y, Zhang H, Zhang X, Zhang H. Comparative transcriptomic and proteomic analyses provide insights into functional genes for hypoxic adaptation in embryos of Tibetan chickens. Sci Rep. 2020;10(1):11213. https://0-doi-org.brum.beds.ac.uk/10.1038/s41598-020-68178-w.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Li J, Ma W, Zeng P, Wang J, Geng B, Yang J, Cui Q. LncTar: a tool for predicting the RNA targets of long noncoding RNAs. Brief Bioinform. 2015;16(5):806–12. https://0-doi-org.brum.beds.ac.uk/10.1093/bib/bbu048.

    Article  CAS  PubMed  Google Scholar 

  20. Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10(3):R25. https://0-doi-org.brum.beds.ac.uk/10.1186/gb-2009-10-3-r25.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Lv J, Cui W, Liu H, He H, Xiu Y, Guo J, Liu H, Liu Q, Zeng T, Chen Y, Zhang Y, Wu Q. Identification and characterization of long non-coding RNAs related to mouse embryonic brain development from available transcriptomic data. PLoS One. 2013;8(8):e71152. https://0-doi-org.brum.beds.ac.uk/10.1371/journal.pone.0071152.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  22. Kelley D, Rinn J. Transposable elements reveal a stem cell-specific class of long noncoding RNAs. Genome Biol. 2012;13(11):R107. https://0-doi-org.brum.beds.ac.uk/10.1186/gb-2012-13-11-r107.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Finn RD, Bateman A, Clements J, Coggill P, Eberhardt RY, Eddy SR, Heger A, Hetherington K, Holm L, Mistry J, Sonnhammer ELL, Tate J, Punta M. Pfam: the protein families database. Nucleic Acids Res. 2014;42(Database issue):D222–30. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkt1223.

    Article  CAS  PubMed  Google Scholar 

  24. Kong L, Zhang Y, Ye ZQ, Liu XQ, Zhao SQ, Wei L, et al. CPC: assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007;35(Web Server issue):W345–9.

    Article  PubMed  PubMed Central  Google Scholar 

  25. Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, Liu Y, Chen R, Zhao Y. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013;41(17):e166. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkt646.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  26. Wang L, Park HJ, Dasari S, Wang S, Kocher JP, Li W. CPAT: coding-potential assessment tool using an alignment-free logistic regression model. Nucleic Acids Res. 2013;41(6):e74. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkt006.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  27. Friedlander MR, Mackowiak SD, Li N, Chen W, Rajewsky N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. 2012;40(1):37–52. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkr688.

    Article  CAS  PubMed  Google Scholar 

  28. Betel D, Wilson M, Gabow A, Marks DS, Sander C. The microRNA.org resource: targets and expression. Nucleic Acids Res. 2008;36(Database issue):D149–53. https://0-doi-org.brum.beds.ac.uk/10.1093/nar/gkm995.

    Article  CAS  PubMed  Google Scholar 

  29. Lewis BP, Shih IH, Jones-Rhoades MW, Bartel DP, Burge CB. Prediction of mammalian microRNA targets. Cell. 2003;115(7):787–98. https://0-doi-org.brum.beds.ac.uk/10.1016/S0092-8674(03)01018-3.

    Article  CAS  PubMed  Google Scholar 

  30. Pertea M, Kim D, Pertea GM, Leek JT, Salzberg SL. Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat Protoc. 2016;11(9):1650–67. https://0-doi-org.brum.beds.ac.uk/10.1038/nprot.2016.095.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  31. Li B, Ruotti V, Stewart RM, Thomson JA, Dewey CN. RNA-Seq gene expression estimation with read mapping uncertainty. Bioinformatics. 2010;26(4):493–500. https://0-doi-org.brum.beds.ac.uk/10.1093/bioinformatics/btp692.

    Article  CAS  PubMed  Google Scholar 

  32. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15(12):550. https://0-doi-org.brum.beds.ac.uk/10.1186/s13059-014-0550-8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Wang J, Cao Y, Lu X, Wang X, Kong X, Bo C, Li S, Bai M, Jiao Y, Gao H, Yao X, Ning S, Wang L, Zhang H. Identification of the regulatory role of lncRNA SNHG16 in myasthenia gravis by constructing a competing endogenous RNA network. Mol Ther Nucl Acids. 2020;19:1123–33. https://0-doi-org.brum.beds.ac.uk/10.1016/j.omtn.2020.01.005.

    Article  CAS  Google Scholar 

  34. Liang R, Han B, Li Q, Yuan Y, Li J, Sun D. Using RNA sequencing to identify putative competing endogenous RNAs (ceRNAs) potentially regulating fat metabolism in bovine liver. Sci Rep. 2017;7(1):6396. https://0-doi-org.brum.beds.ac.uk/10.1038/s41598-017-06634-w.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. Zuo J, Wang Y, Zhu B, Luo Y, Wang Q, Gao L. Analysis of the coding and non-coding RNA transcriptomes in response to bell pepper chilling. Int J Mol Sci. 2018;19(7):1–15.  https://0-doi-org.brum.beds.ac.uk/10.3390/ijms19072001.

  36. Salmena L, Poliseno L, Tay Y, Kats L, Pandolfi PP. A ceRNA hypothesis: the Rosetta stone of a hidden RNA language? Cell. 2011;146(3):353–8. https://0-doi-org.brum.beds.ac.uk/10.1016/j.cell.2011.07.014.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Su X, Xing J, Wang Z, Chen L, Cui M, Jiang B. microRNAs and ceRNAs: RNA networks in pathogenesis of cancer. Chinese J Cancer Res. 2013;25(2):235–9.

    Google Scholar 

  38. Dusseau JW, Hutchins PM. Hypoxia-induced angiogenesis in chick chorioallantoic membranes: a role for adenosine. Respir Physiol. 1988;71(1):33–44. https://0-doi-org.brum.beds.ac.uk/10.1016/0034-5687(88)90113-2.

    Article  CAS  PubMed  Google Scholar 

  39. Schaad NC, De Castro E, Nef S, Hegi S, Hinrichsen R, Martone ME, et al. Direct modulation of calmodulin targets by the neuronal calcium sensor NCS-1. P Natl Acad Sci USA. 1996;93(17):9253–8. https://0-doi-org.brum.beds.ac.uk/10.1073/pnas.93.17.9253.

    Article  CAS  Google Scholar 

  40. Hui AS, Bauer AL, Striet JB, Schnell PO, Czyzyk-Krzeska MF. Calcium signaling stimulates translation of HIF-alpha during hypoxia. FASEB J. 2006;20(3):466–75. https://0-doi-org.brum.beds.ac.uk/10.1096/fj.05-5086com.

    Article  CAS  PubMed  Google Scholar 

  41. Feng Y, Liu Y, Cao PX, Sun X, Li KX, Li XY, et al. Estrogen-dependent MicroRNA-504 Expression and Related Baroreflex Afferent Neuroexcitation via Negative Regulation on KCNMB4 and KCa1.1 beta4-subunit Expression. Neuroscience. 2020;442:168–82.

    Article  CAS  PubMed  Google Scholar 

  42. Bruning U, Cerone L, Neufeld Z, Fitzpatrick SF, Cheong A, Scholz CC, et al. MicroRNA-155 promotes resolution of hypoxia-inducible factor 1alpha activity during prolonged hypoxia. Mol Cell Biol. 2011;31(19):4087–96. https://0-doi-org.brum.beds.ac.uk/10.1128/MCB.01276-10.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Wu X, Li J, Yang X, Bai X, Shi J, Gao J, et al. miR-155 inhibits the formation of hypertrophic scar fibroblasts by targeting HIF-1alpha via PI3K/AKT pathway. J Mol Histol. 2018;49(4):377–87. https://0-doi-org.brum.beds.ac.uk/10.1007/s10735-018-9778-z.

    Article  CAS  PubMed  Google Scholar 

  44. Qin C, Zha W, Fan R, Ding H, Xu Y, Wang C. MicroRNA302a inhibits cell proliferation and invasion, and induces cell apoptosis in hepatocellular carcinoma by directly targeting VEGFA. Mol Med Rep. 2017;16(5):6360–7. https://0-doi-org.brum.beds.ac.uk/10.3892/mmr.2017.7312.

    Article  CAS  PubMed  Google Scholar 

  45. Sun Y, Xiong X, Wang X. RELA promotes hypoxia-induced angiogenesis in human umbilical vascular endothelial cells via LINC01693/miR-302d/CXCL12 axis. J Cell Biochem. 2019;120(8):12549–58. https://0-doi-org.brum.beds.ac.uk/10.1002/jcb.28521.

    Article  CAS  PubMed  Google Scholar 

  46. Cao YP, Pan M, Song YL, Zhang HL, Sui HT, Shan BC, et al. MiR-302 a/b/c suppresses tumor angiogenesis in hepatocellular carcinoma by targeting MACC1. Eur Rev Med Pharmaco. 2019;23(18):7863–73.

    Google Scholar 

  47. Fraidenburg DR, Machado RF. A review of Transcriptome analysis in pulmonary vascular diseases. Methods Mol Biol. 2018;1783:259–77.

    Article  PubMed  PubMed Central  Google Scholar 

  48. Liu P, Yang F, Zhuang Y, Xiao Q, Cao H, Zhang C, et al. Dysregulated expression of microRNAs and mRNAs in pulmonary artery remodeling in ascites syndrome in broiler chickens. Oncotarget. 2017;8(2):1993–2007. https://0-doi-org.brum.beds.ac.uk/10.18632/oncotarget.12888.

    Article  PubMed  Google Scholar 

  49. Ribatti D, Nico B, Belloni AS, Vacca A, Roncali L, Nussdorfer GG. Angiogenic activity of leptin in the chick embryo chorioallantoic membrane is in part mediated by endogenous fibroblast growth factor-2. Int J Mol Med. 2001;8(3):265–8. https://0-doi-org.brum.beds.ac.uk/10.3892/ijmm.8.3.265.

    Article  CAS  PubMed  Google Scholar 

  50. Gao Y, Yu X, Zhang F, Dai J. Propofol inhibits pancreatic cancer progress under hypoxia via ADAM8. J Hepato-Bil-Pan Sci. 2019;26(6):219–26.

    Google Scholar 

  51. Mahoney ET, Benton RL, Maddie MA, Whittemore SR, Hagg T. ADAM8 is selectively up-regulated in endothelial cells and is associated with angiogenesis after spinal cord injury in adult mice. J Comp Neurol. 2009;512(2):243–55. https://0-doi-org.brum.beds.ac.uk/10.1002/cne.21902.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Meir JU, York JM, Chua BA, Jardine W, Hawkes LA, Milsom WK. Reduced metabolism supports hypoxic flight in the high-flying bar-headed goose (Anser indicus). Elife. 2019;8. https://0-doi-org.brum.beds.ac.uk/10.7554/eLife.44986.

  53. Yang J, Li WR, Lv FH, He SG, Tian SL, Peng WF, et al. Whole-genome sequencing of native sheep provides insights into rapid adaptations to extreme environments. Mol Biol Evol. 2016;33(10):2576–92. https://0-doi-org.brum.beds.ac.uk/10.1093/molbev/msw129.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Schippers MP, Ramirez O, Arana M, Pinedo-Bernal P, McClelland GB. Increase in carbohydrate utilization in high-altitude Andean mice. Curr Biol. 2012;22(24):2350–4. https://0-doi-org.brum.beds.ac.uk/10.1016/j.cub.2012.10.043.

    Article  CAS  PubMed  Google Scholar 

  55. Gangwar A, Paul S, Ahmad Y, Bhargava K. Intermittent hypoxia modulates redox homeostasis, lipid metabolism associated inflammatory processes and redox post-translational modifications: benefits at high altitude. Sci Rep. 2020;10(1):7899. https://0-doi-org.brum.beds.ac.uk/10.1038/s41598-020-64848-x.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Wang XP, Yang J, Norman MA, Magnusson J, DeMayo FJ, Brunicardi FC. SSTR5 ablation in islet results in alterations in glucose homeostasis in mice. FEBS Lett. 2005;579(14):3107–14. https://0-doi-org.brum.beds.ac.uk/10.1016/j.febslet.2005.04.069.

    Article  CAS  PubMed  Google Scholar 

  57. Ma K, Saha PK, Chan L, Moore DD. Farnesoid X receptor is essential for normal glucose homeostasis. J Clin Invest. 2006;116(4):1102–9. https://0-doi-org.brum.beds.ac.uk/10.1172/JCI25604.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Priyadarshini S, Pradhan B, Griebel P, Aich P. Cortisol regulates immune and metabolic processes in murine adipocytes and macrophages through HTR2c and HTR5a serotonin receptors. Eur J Cell Biol. 2018;97(7):483–92. https://0-doi-org.brum.beds.ac.uk/10.1016/j.ejcb.2018.07.004.

    Article  CAS  PubMed  Google Scholar 

  59. Lai CQ, Parnell LD, Ordovas JM. The APOA1/C3/A4/A5 gene cluster, lipid metabolism and cardiovascular disease risk. Curr Opin Lipidol. 2005;16(2):153–66. https://0-doi-org.brum.beds.ac.uk/10.1097/01.mol.0000162320.54795.68.

    Article  CAS  PubMed  Google Scholar 

  60. He X, Zheng J, He Y, Wang Y, Wang L, Bai M, et al. Long non-coding RNA LINC-PINT and LINC00599 polymorphisms are associated with high-altitude pulmonary edema in Chinese. Arch Bronconeumol. 2020;56(6):360–4. https://0-doi-org.brum.beds.ac.uk/10.1016/j.arbres.2019.09.021.

    Article  PubMed  Google Scholar 

  61. Zhang Y, Liu L, Liang C, Zhou L, Tan L, Zong Y, et al. Expression profiles of long noncoding RNAs in mice with high-altitude hypoxia-induced brain injury treated with Gymnadenia conopsea (L.) R. Br. Neuropsych Dis Treat. 2020;16:1239–48. https://0-doi-org.brum.beds.ac.uk/10.2147/NDT.S246504.

    Article  CAS  Google Scholar 

  62. Shang P, Li W, Liu G, Zhang J, Li M, Wu L, et al. Identification of lncRNAs and genes responsible for fatness and fatty acid composition traits between the Tibetan and Yorkshire pigs. Int J Genomics. 2019;2019:5070975.

    Article  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

The authors would like to thank BioMarker Technology (Beijing, China) for the sequencing services and bioinformatic analysis.

Funding

This work was supported by the National Natural Science Foundation of China (31972532), the China Agricultural Research System (CARS-40-K05), and the Innovation Base Cultivation and Development Project (Z171100002217072).

Author information

Authors and Affiliations

  1. National Engineering Laboratory for Animal Breeding, Plateau Animal Genetic Resources Center, China Agricultural University, No. 2 Yuanmingyuan West Rd, Haidian District, Beijing, 100193, China

    Ying Zhang, Woyu Su, Bo Zhang, Yao Ling & Hao Zhang

  2. Department of Poultry Science, University of Georgia, 303 Poultry Science Building, Athens, GA, 30602, USA

    Woo Kyun Kim

Authors
  1. Ying Zhang
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Woyu Su
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Bo Zhang
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Yao Ling
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Woo Kyun Kim
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Hao Zhang
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

Y.Z., and H.Z. designed the research; Y.Z., W.S., Y.L., and B.Z. collected samples, extracted RNAs, and performed qRT-PCR verification experiments; Y.Z. performed most of the bioinformatics studies; Y.Z., H.Z., and W.K. wrote and revised the manuscript; and all authors reviewed the manuscript. The authors read and approved the final manuscript.

Corresponding authors

Correspondence to Woo Kyun Kim or Hao Zhang.

Ethics declarations

Ethics approval and consent to participate

All animal work was conducted following the guidelines for the care and use of experimental animals. The experiments and animal care protocol were approved by the animal welfare committee of the State Key Laboratory for Agro-Biotechnology of the China Agricultural University (approval number, XK257).

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Supplementary Information

Additional file 1: Table S1.

The qRT-PCR primers of six DE lncRNAs, six DE miRNAs and seven DE genes for identification by RNA-seq.

Additional file 2: Table S2.

List of software and the parameters used in this study.

Additional file 3: Table S3.

Summary of sequencing reads aligned with the Gallus gallus genome.

Additional file 4: Figure S1.

Lengths of mRNAs and lncRNAs distribution. Figure S2. Figure of hierarchical clustering of 389 DE lncRNAs (a), 73 DE miRNAs (b), and 354 DE genes (c) in Tibetan and Chahua chickens. Figure S3. Comparative analysis of exon and ORF lengths of lncRNAs and mRNAs. Figure S4. Distribution of the DE lncRNAs and DE mRNAs on the different chromosomes. Figure S5. Network of ceRNA (lncRNA-miRNA-mRNA) including 162 DE lncRNAs, 25 DE miRNAs, and 108 DE genes.

Additional file 5: Table S4.

The 3194 known lncRNAs, 471 known miRNAs, and 968 novel miRNAs.

Additional file 6: Table S5.

The 389 differentially expressed lncRNAs (DE lncRNAs), 73 DE miRNAs, and 354 DE genes.

Additional file 7: Table S6.

GO and KEGG pathway enrichment of differentially expressed lncRNA (DE lncRNA) target genes.

Additional file 8: Table S7.

GO and KEGG pathway enrichment of differentially expressed miRNA (DE miRNA) target genes.

Additional file 9: Table S8.

GO and KEGG enrichment of differentially expressed mRNAs (DE genes).

Additional file 10: Table S9.

Important difference in circRNA-miRNA-mRNA interactions based on the competing endogenous RNA (ceRNA) network.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Zhang, Y., Su, W., Zhang, B. et al. Comprehensive analysis of coding and non-coding RNA transcriptomes related to hypoxic adaptation in Tibetan chickens. J Animal Sci Biotechnol 12, 60 (2021). https://0-doi-org.brum.beds.ac.uk/10.1186/s40104-021-00582-2

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://0-doi-org.brum.beds.ac.uk/10.1186/s40104-021-00582-2

Keywords

Journal of Animal Science and Biotechnology

ISSN: 2049-1891