Yak miR-2285o-3p attenuates hypoxia-induced apoptosis by targeting caspase-3.

miRNAs are a class of hairpin-derived RNAs, 21-24 nucleotides in length, which are involved in a range of biological processes. The bta-miR-2285 family has over 40 members spanning the entire bovine genome. We previously found that bta-miR-2285o-3p was highly expressed in yak heart and lung when compared with cattle, which prompted us to investigate its potential function in high-altitude adaptation of yaks. In this study, we detected wide-spread high expression of bta-miR-2285o-3p in yak tissues. Further experiments revealed that the protein tyrosine phosphatase receptor type M (PTPRM) gene was the host gene of bta-miR-2285o-3p and that two linked SNPs in bta-mir-2285o precursor affected the biogenesis of mature miRNA (bta-miR-2285o-3p). Functional analysis in vitro indicated that bta-miR-2285o-3p attenuated hypoxia-induced apoptosis by targeting very low-density lipoprotein receptor (VLDLR), phosphatase and tensin homolog (PTEN) and caspase-3. Expression level analysis in vivo revealed the high negative Pearson's correlation between bta-miR-2285o-3p and caspase3 in yak, highlighting the potential important roles of bta-miR-2285o-3p in yak high-altitude adaptation. Our study provides a typical model for deciphering the function of miRNAs in environmental adaptation.

annexin V‐FITC

BCL2 family apoptosis regulator

propidium iodide

phosphatase and tensin homolog

protein tyrosine phosphatase receptor type M

very low‐density lipoprotein receptor

Introduction

Yak (Bos grunniens), an endemic livestock animal of the Qinghai–Tibet Plateau (Ding et al. 2014), has evolved numerous anatomical and physiological characteristics to equip it for the extreme harsh environment of cold, high ultraviolet levels and especially hypoxia at high altitude, including larger lungs and heart, higher energy metabolism levels and strong environmental sensing (Ding et al. 2014; Zhang et al. 2016). Yak and cattle are distinct species that diverged 5 million years ago (Qiu et al. 2012), but their strong genome similarity and extreme differences in living environments contribute to the use of these two species to illuminate the mechanism of high‐altitude adaptation. Extensive studies have been carried out in a variety of species regarding high‐altitude adaptation at the level of morphology (Hoppeler et al. 1990), anatomy (Ahmad et al. 2016), hemodynamics (Ahmad et al. 2016), physiology (Monge & Leónvelarde 1991) and genomics (Qiu et al. 2012; Emilia et al. 2013), but the miRNA regulatory networks of high‐altitude adaptation remain poorly understood.

miRNAs are hairpin‐derived RNAs ranging from 21 to 24 nucleotides, which post‐transcriptionally regulate the expression of target genes by binding to the 3′ UTR of messenger RNA (mRNA; Murchison & Hannon 2004; Yang et al. 2014a). As an important class of gene regulators, miRNAs are involved in a range of biological processes, including development, differentiation, proliferation and aging. A limited number of reports have illuminated the crucial role of miRNAs in hypoxia adaptation. miR‐15a can promote mesenchymal ablation and adaptation to hypoxia by suppressing Bcl‐2 during lung development in chickens (Hao et al. 2014), and hypoxia‐induced miR‐210 controls metabolic adaptation under hypoxia by targeting the iron–sulfur cluster assembly proteins ISCU1/2 in arterial (Chan et al. 2009) and pulmonary vascular endothelium cells (Chan et al. 2009). Although many animal miRNAs show conservation in sequence, function and expression pattern (Ninova et al. 2015), miRNA evolution is relatively dynamic, resulting in rapid appearance and extinction and lineage‐specific expansions (Andrew Grimson, 2008; Christodoulou et al. 2010; Yang et al. 2014a). Previously, the miR‐2285 family was shown to be lineage specific only in ruminants with relatively low levels of expression, including cattle (Sun et al. 2014), yak (Guan et al. 2017) and sheep (Morrison et al. 2015), and composed of 40 members spanning the entire bovine genome (Lawless et al. 2014). Although miR‐2285 family expression has been demonstrated to be widespread in ruminants, few studies have focused on biological functions and regulatory networks owing to their extremely low levels of expression.

In our previous study, we identified a specific member of the bta‐miR‐2285 family, bta‐miR‐2285o‐3p, which exhibited high levels of expression in the yak heart and lung when compared with cattle (Guan et al. 2017), prompting us to investigate the evolutionary mechanism and biological function of bta‐miR‐2285o in yaks. Here, we systematically compared the tissue expression pattern of bta‐miR‐2285o‐3p between yaks and cattle, and explored the molecular mechanism underlying its yak‐biased expression. We found that two linked SNP sites in pre‐miR‐2285o resulted in the disturbance of bta‐miR‐2285o maturation. Gain and loss functional analysis in vitro revealed that bta‐miR‐2285o‐3p was potentially associated with high‐altitude adaptation through attenuating hypoxia‐induced apoptosis by targeting phosphatase and tensin homolog (PTEN) and caspase‐3. Expression level analysis in five tissues of yak revealed a high negative Pearson’s correlation between bta‐miR‐2285o‐3p and caspase‐3.

Materials and methods

Cell culture and hypoxia treatment

We purchased the rat cardiomyoblast cell line H9C2 from the cell bank of the Chinese Academy of Sciences and cultured it in Dulbecco’s Modified Eagle Medium (Hyclone) supplemented with 10% fetal bovine serum. Normoxic cells were cultured in a humidified 37°C incubator with 5% CO2. The cells were exposed to hypoxic conditions for 48 h, at which point they reached approximately 50% confluency. Hypoxic cells were cultured in the hypoxic chamber using AnaeroPack (Mitsubishi Gas Company) to mimic a hypoxic environment in vitro, with other conditions consistent with those of the normoxic cells. The AnaeroPack started to absorb oxygen within 1 min and oxygen tension inside the chamber dropped to 1 mmHg within 1 h (O2 < 1%, CO2 ≈ 5%).

Flow cytometry analysis of cell apoptosis

Flow cytometry was used to detect cell apoptosis rates with an Annexin V apoptosis detection kit (BD Biosciences). Cells were trypsin‐digested and collected by centrifugation at 1000 for 5 min. Each sample was double stained for 10 min with Annexin V‐FITC (AV) and propidium iodide (PI) and then analyzed by flow cytometry. Cells stained with AV and PI or AV only were considered necrotic and apoptotic respectively. Experiments were performed in triplicate.

Genotyping

Genomic DNA was extracted from the 48 yak and 50 cattle heart tissue using the TIANamp Genomic DNA Kit (Tiangen). PCR amplification of a 1027 bp fragment containing the pri‐miR‐2285o was carried out with reagents from the TaKaRa Taq™Polymerase Kit. The PCR reaction was 30 cycles of denaturation at 97°C for 12 s, annealing at 58°C for 30 s and extension at 72°C for 1 min, and the resulting PCR products were sequenced (TSINGKE).

Plasmid construction

The 298 bp fragment containing pri‐bta‐miR‐2285o from either yak or cattle and harboring XhoI and BamHI restriction enzyme cleavage sites was ligated into the pSilencer™ 3.1 H1 neo plasmid (Ambion).

Dual‐luciferase activity assay

The potentially targeted mRNA (VLDLR, PTEN and caspase‐3) containing the specific miRNA binding sites (wild type (Wt) or mutant (Mut)) was synthesized from TSINGKE. The sequences were cleaved using Sac I/Xho I and cloned into the pmirGLO plasmid (Promega) at the 3′‐end of the firefly luciferase reporter (luc2) gene. HeLa cells were cultured in 96‐well plates until the cell density reached about 80% confluence. Recombinant pmirGLO vector with Wt or Mut was cotransfected with miRNA mimics into cells using Lipofectamine3000 (Invitrogen). Cells were collected after 48 h and dual‐luciferase activity was measured using the Dual‐Luciferase Reporter Assay System kit (Promega), according to the manufacturer’s instructions.

RT‐PCR

Total RNA (including miRNA) was extracted using TRIzol Reagent (Invitrogen). mRNA and miRNA were reverse‐transcribed respectively using PrimeScript RT reagent Kit with gDNA Eraser and Mir‐X™ miRNA First Strand Synthesis Kit (Takara), according to the manufacturers’ recommendations. RT‐PCR was performed using an SYBR Premix Ex Taq kit (Takara) and a CFX96 system (Bio‐Rad). The PCR reaction was 35 cycles of 95°C for 10 s, 56°C for 20 s and extension at 72°C for 1 min. All reactions were performed in triplicate. Relative expression levels of mRNAs and miRNAs were calculated using the 2−ΔΔ method. GAPDH and U6 were used as housekeeping genes for normalizing mRNA and miRNA respectively. Sequences of the primers used for RT‐PCR are shown in Table S1.

Statistical analysis

All data are expressed as means ± SD. Statistical significance was calculated by one‐way analysis of variance (ANOVA) with Tukey post‐hoc test for multiple group or Student’s t‐test for two group comparisons, using SPSS 19.0 software (SPSS Inc.). Values of P < 0.05 were considered statistically significant (*P < 0.05; **P < 0.01).

Results

bta‐miR‐2285o‐3p was highly and ubiquitously expressed in yak

In our previous study (Guan et al. 2017), we found that bta‐miR‐2285o‐3p was specifically expressed in yaks but not in cattle. Bta‐miR‐2285o‐3p showed 8810 and 8033 reads per million in the yak heart and lung respectively, but it was almost undetectable in cattle (Fig. 1a). In addition, the expression level of bta‐miR‐2285o‐3p ranked within the top 20 miRNAs in both the heart and lung of yaks (data not shown), indicating the potential importance of this miRNA in yaks. To further confirm the species and tissue specificity of bta‐miR‐2285o‐3p, we determined its expression profile in 12 tissues from yaks and cattle using RT‐PCR. As a result, bta‐miR‐2285o‐3p showed the highest expression levels in yak longissimus dorsi, followed by heart, lung, omentum majus, psoas major, subcutaneous fat, liver, spleen, perirenal fat, duodenum, kidney and ileum, but remained nearly undetectable in cattle tissues (Fig. 1b).

Figure 1

Expression of bta‐miR‐2285o‐3p in yaks and cattle. (a) bta‐miR‐2285o‐3p showed high expression in yak heart and lung. (b) Expression profile of bta‐miR‐2285o‐3p across 12 tissues in the yak.

Mutations in yak pre‐miR‐2285o affects the biogenesis of mature miRNA in vitro

To identify candidate precursor sequence of bta‐miR‐2285o, we mapped mature bta‐miR‐2285o‐3p to the yak genome and examined whether the extended sequence surrounding the mature sequence of bta‐miR‐2285o‐3p could potentially form hairpins. Extended sequences that potentially form hairpins were selected as candidate precursor sequences of bta‐miR‐2285o. As a result, five potential precursor sequences of bta‐miR‐2285o across five genomic locations were identified in the bovine genome, including two intergenic and three intronic locations (Fig. 2a,b). Among the five candidate miRNA precursors, the two intergenic ones showed higher thermodynamic free energy values larger than −30 kcal/mol, implying that the other three intronic precursors which showed lower thermodynamic free energy values, were more likely to be the real origin of bta‐miR‐2285o (Fig. 2c). We therefore calculated the Spearman’s correlation between bta‐miR‐2285o‐3p and respective corresponding host genes and revealed that PTPRM showed a significantly moderate correlation coefficient (R 2 = 0.56, P < 0.01) with bta‐miR‐2285o‐3p, whereas no significant correlation was detected for Flt and RRBP (Fig. 2d). These results indicated that PTPRM was the most likely host gene for the bta‐miR‐2285o precursor (pre‐bta‐miR‐2285o).

Figure 2

Mutation in pre‐miR‐bta‐2285o affects the biogenesis of mature miRNA in vitro. (a) Five pre‐miRNAs across the five genomic locations in bovine genome. (b) Multiple alignment of the five pre‐miR‐2285o sequences in bovine. (c) Secondary structures of five miR‐2285o precursors predicted by mfold. (d) Spearman’s correlation of miR‐2285‐3p and three potential host genes. (e) Secondary structures of the yak and bovine pre‐miR‐2285o‐4 as predicted by mfold. The miR‐2285o‐3p mature sequence is highlighted in red. The mutation is indicated by yellow squares. Free energies calculated using mfold are indicated. (g) The positions and sequences of the primers used for RT‐ PCR to quantify the level of pri‐miRNA and pre‐miRNA. Note that the primer for pre‐miRNA (F2) also recognizes the corresponding pri‐miRNA and therefore measures both forms. F1 and R represent pri‐miRNA, whereas F2 and R represent both pri‐ and pre‐miRNA. (h) The two linked SNPs significantly alter the expression of mature miR‐2285o‐3p in HeLa cells. U6 and miR‐16 were used as reference genes.

We next cloned and sequenced the intronic sequence of PTPRM across 98 individuals, including 48 yaks and 50 cattle (Table S2). We identified two single nucleotide polymorphic sites in the precursor sequence of bta‐miR‐2285o (defined as pre‐miR‐2285o), 67A/G and 70G/A. 67A/G and 70G/A were found to have strong linkage disequilibrium (D′ = 1, R 2 = 0.885), so haplotype analysis was carried out. A total of three haplotypes were identified in yaks and cattle, and all showed significant frequency bias between yaks and cattle (Table 1). We predicted the secondary structures and free energies of pre‐miR‐2285o for yaks and cattle using Mfold and RNAfold, and revealed altered thermodynamic stability and secondary structures, which may affect both the biogenesis and stability of pre‐bta‐miR‐2285o (Fig. 2e).

Table 1

Haplotype analysis of pre‐miR‐2285o polymorphisms

Haplotype	Cattle (frequency)	Yak (frequency)	χ 2	Fisher’s P‐value	Pearson’s P‐value	Odds ratio [95% CI]
A/G*	11.00 (11.00%)	82.00 (85.40%)	108.78	5.00 × 10−15	1.95 × 10−25	0.02 [0.01–0.05]
G/A*	89.00 (89.00%)	8.00 (8.30%)	127.50	1.40 × 10−25	1.56 × 10−29	89.00 [34.17–231.80]
G/G*	0.00 (0.00%)	6.00 (6.30%)	6.45	0.01	0.01	—

To explore the potential roles of these two linked SNPs in miRNA biogenesis, we constructed expression vectors for yak and cattle bta‐miR‐22855o (Fig. 2f). Specific primers were designed to recognize either primary bta‐miR‐2285o (pri‐bta‐miR‐2285o) alone or pri‐bta‐miR‐2285o and pre‐bta‐miR‐2285o together (Fig. 2g). Upon transfection into human HeLa cells which have no endogenous bta‐miR‐2285o expression, we detected the mature form of yak bta‐miR‐2285o‐3p (A/G) by RT‐PCR; however, we failed to detect the mature form of cattle bta‐miR‐2285o‐3p (G/A). Additionally, whereas yak and cattle pri‐miR‐2285o showed similar expression levels, the overall expression level of pri‐ and pre‐bta‐miR‐2285o in cattle was only ~14% of that in yak (P < 0.01; Fig. 2h). Our results suggested that the SNP around the 3′ end of the precursor sequence altered the processing of pri‐bta‐miR‐2285o and affected the biogenesis of pre‐bta‐miR‐2285o, which then contributed to the yak‐specific high expression of mature bta‐miR‐2285o‐3p.

bta‐miR‐2285o‐3p attenuated hypoxia‐induced apoptosis through targeting the Caspase‐3

To investigate the biological function of yak‐specific bta‐miR‐2285o‐3p, we performed target gene prediction and identified 504 candidate target genes using targetscancow7.1 (http://www.targetscan.org/vert_71/). Functional enrichment analysis revealed that these target genes were mainly enriched in hypoxia‐related terms (Fig. 3a), including the p53 signaling pathway (seven genes, P = 6.8 × 10−2), the cellular response to hypoxia (five genes, P = 3.9 × 10−2) and the response to hypoxia (five genes, P = 1×10−2) (Table 2). As a tetrameric transcription factor, p53 mediates the cellular response to various detrimental stresses and induces the transcription of hundreds of target genes, which are involved in regulating apoptosis, the cell cycle and DNA repair. p53 is maintained at low levels by murine double minute 2‐mediated ubiquitination and subsequent proteasomal degradation, whereas it is stabilized by various cellular stresses such as DNA damage, hypoxia and oxidative stress ( Sermeus & Michiels, 2011). Some apoptosis‐related genes were predicted to be the target genes of bta‐miR‐2285o‐3p, including PTEN, caspase‐3 and BOK (a BCL2 family apoptosis regulator). Therefore, we suppose that bta‐miR‐2285o‐3p could target these genes and attenuate the hypoxia‐induced apoptosis.

Figure 3

Overexpression of bta‐miR‐2285o‐3p attenuate hypoxia‐induced apoptosis by targeting of VLDLR, PTEN and Caspase‐3. (a) Functional enrichment analysis of 504 candidate target genes predicted using targetscancow7.1. (b) Relative miRNA expression in H9C2 cells transfected with a bta‐miR‐2285o‐3p mimic and inhibitor. (c–e) Apoptosis rate (c), cell viability (d) and membrane integrity (e) were evaluated by CCK8, flow cytometry analysis and LDH release assay, respectively. (f) Expression levels of apoptosis‐related genes were measured by RT‐PCR. (g) Apoptosis rate was determined through flow cytometry analysis. (h) Sequence comparison of three candidate genes across mutiple species. (i) Luciferase reporter assay data for three candidate genes. (j) Spearman’s correlation between the expression level of bta‐miR‐2285o‐3p and three target genes in vivo. All experiments performed in triplicate and all data are expressed as means ± SD. *P < 0.05, **P < 0.01.

Table 2

The detailed information about functional annotation of target genes

Terms	Gene name	Gene ID
Response to hypoxia	CBFA2T3	CBFA2/RUNX1 translocation partner 3
	SMAD3	SMAD family member 3
	ANG	Angiogenin, ribonuclease, RNase A family, 5
	EDNRA	Endothelin receptor type A
	MECP2	Methyl‐CpG binding protein 2
	NOS1	Nitric oxide synthase 1
	SLC2A8	Solute carrier family 2 member 8
Cellular response to hypoxia	MDM4	MDM4, p53 regulator
	RORA	RAR‐related orphan receptor A
	PTEN	Phosphatase and tensin homolog
	PDGFB	Platelet‐derived growth factor subunit B
	UBQLN1	Ubiquilin 1
P53 signaling pathway	MDM4	MDM4, p53 regulator
	CASP3	Caspase 3, apoptosis‐related cysteine peptidase
	CCND2	Cyclin D2
	GADD45G	Growth arrest and DNA damage inducible gamma
	PTEN	Phosphatase and tensin homolog
Positive regulation of reactive oxygen Species metabolic process	DUOXA1	Dual oxidase maturation factor 1
	PDGFB	Platelet‐derived growth factor subunit B
	TGFBR2	Transforming growth factor beta receptor 2
Intrinsic apoptotic signaling pathway by p53 class mediator	BOK	BOK, BCL2 family apoptosis regulator
	DDX5	DEAD‐box helicase 5
	E2F2	E2F transcription factor 2
Multicellular organismal response to stress	NOS1	Nitric oxide synthase 1
	PTEN	Phosphatase and tensin homolog
	PPP3CA	Protein phosphatase 3 catalytic subunit alpha
Regulation of glucose metabolic process	RORA	RAR‐related orphan receptor A
	IGFBP4	Insulin‐like growth factor binding protein 4
	NCOA2	Nuclear receptor coactivator 2
	PDK3	Pyruvate dehydrogenase kinase 3

To verify our hypothesis that bta‐miR‐2285o‐3p attenuates hypoxia‐induced apoptosis, we investigated the biological function of bta‐miR‐2285o‐3p in H9C2 cells. Following the transfection of a miRNA mimic or inhibitor, we performed gain‐ and loss‐of‐function analysis. The transfection of a miR‐2285o‐3p mimic or inhibitor significantly increased or decreased miR‐2285o‐3p expression respectively (Fig. 3b), confirming the suitability of the miRNA transfection system. After exposing transfected cells to hypoxia for 48 h, we investigated cell viability, membrane integrity and the apoptosis rate. Under hypoxic conditions, the overexpression of miR‐2285o‐3p promoted cell survival (Fig. 3c,g), significantly suppressed the reduction of cell viability (P < 0.05; Fig. 3d) and reduced cell membrane damage (P < 0.05; Fig. 3e). Additionally, three pro‐apoptotic genes, including p53, PTEN and caspase‐3, were significantly inhibited whereas the anti‐apoptotic gene Bcl‐2 was increased (Fig. 3f). Taken together, these results suggest that miR‐2285o‐3p could attenuate hypoxia‐induced apoptosis in H9C2 cells.

Next, the potential apoptosis‐related target genes of miR‐2285o‐3p were screened according to the results of the target gene prediction software. Intriguingly, we found the direct apoptosis‐related gene (PTEN and caspase‐3) and VLDLR. A previous study reported hypoxia‐induced apoptosis through regulating endoplasmic reticulum stress (Yang et al. 2014b), which might elucidate the antiapoptotic effects of miR‐2285o‐3p. Afterwards, we used a luciferase reporter assay to verify whether VLDLR, PTEN and caspase‐3 were the authentic target genes of miR‐2285o‐3p. Sequence analysis revealed that the target sequence of all three genes was highly conserved across multiple species (Fig. 3h–j), including humans, mice and rats. This enabled us to perform functional analysis in HeLa cells. We cloned and constructed a dual luciferase reporter plasmid encoding the yak gene 3′‐UTR containing bta‐miR‐2285o‐3p binding sites. This reporter plasmid (with Wt or Mut yak gene 3′‐UTRs) was co‐transfected with a bta‐miR‐2285o‐3p mimic into HeLa cells, and luciferase activity was detected 48 h after transfection.

In the case of VLDLR, luciferase activity was downregulated by 90% when transfecting a bta‐miR‐2285o‐3p mimic compared with the negative control. Similarly, luciferase activity was reduced by 33.2 and 50% for PTEN and caspase‐3 respectively (Fig. 3h–j). These results demonstrated that VLDLR, PTEN and caspase‐3 were the target genes of bta‐miR‐2285o‐3p. VLDLR has previously been shown to be involved in hypoxia‐induced endothelial endoplasmic reticulum stress and apoptosis (Carracedo et al. 2011), whereas HIF‐1‐mediated VLDLR upregulation influenced intracellular lipid accumulation by regulating LDL and VLDL uptake under hypoxia. The tumor suppressor PTEN has also been implicated in various processes such as cell cycle progression, migration, spreading and growth (Carracedo et al. 2011). Moreover, PTEN activation increased neuronal apoptosis in the developing rat brain after hypoxia–ischemia (Zhao et al. 2013). Finally, as a well‐established initiator of apoptosis, caspase‐3 has been implicated as a key downstream molecule that executes the apoptotic cascade after hypoxia (Sun et al. 2006).

To test the potential roles of bta‐miR‐2285o‐3p in hypoxia adaptation in vivo, we determined the expression level of bta‐miR‐2285o‐3p and three target genes in five tissues. Thereafter, the Pearson’s correlation was calculated between bta‐miR‐2285o‐3p and target genes. As a result, VLDLR and PTEN showed a non‐significantly positive correlation with bta‐miR‐2285o‐3p (VLDLR, R = 0.39, P > 0.05; PTEN, R = 0.04, P > 0.05); in comparison, caspase‐3 showed a significant negative correlation with bta‐miR‐2285o‐3p (R = −0.68, P < 0.01; Fig. 3h–j). These results indicate the potential roles of bta‐2285o‐3p in yak high‐altitude adaptation through targeting of caspase‐3.

Discussion

Yaks, the largest native highland mammals in Tibetan Plateau, provide protein food and labor to native Tibetans. In recent years, whole‐genome sequencing, transcriptome sequencing and epigenetic analysis have been used to reveal the molecular mechanisms underlying yak adaptation to high‐altitude environments. As a result, large numbers of related candidate genes and pathways have been identified, mainly including hypoxia‐related pathways, metabolism‐related pathways such as PI3K–Akt signaling pathway and fatty acid metabolism (Azad et al. 2017). In addition, many of the miRNAs targeting these pathways have been identified as being involved. In a previous study, our research group compared the miRNA profiles of the heart and lung between yaks and cattle, and target genes of differentially expressed miRNAs were associated with nucleotide excision repair, cell cycles and apoptosis‐related pathways (Guan et al. 2017). Interestingly, we accidentally found that bta‐miR‐2285o‐3p exhibited high levels of expression in yaks compared with cattle. In this study, we determined the expression levels of bta‐miR‐2285o‐3p across 12 tissues comparing yaks and cattle. This results revealed a high expression level of bta‐miR‐2285o‐3p in yak muscles, lung and adipose tissues but relatively low expression level in cattle, suggesting a divergence of expression features of miR‐2285o‐3p between yaks and cattle. However, the underlying mechanism is yet to be determined.

Interestingly, five genomic locations (two intergenic and three intronic locations) were identified as potential precursor sequences of bta‐miR‐2285o. It is difficult to distinguish the real precursor sequence among the five candidate locations. Here, we arbitrarily exclude the two intergenic locations according to their higher thermodynamic free energy values and the fact that intron is a ‘hot spot’ for the emergence of novel miRNAs owing to the non‐necessity of the evolution of a new promoter unit (Berezikov & Eugene 2011). Next, Spearman’s correlation analysis between bta‐miR‐2285o‐3p and three candidate host genes in yak revealed a significant moderate correlation coefficient (R 2 = 0.56) for PTPRM, with no significant correlation for Flt and RRBP. The moderate correlation coefficient indicates that PTPRM was the most likely host gene of bta‐miR‐2285o‐3p; however it is unknown whether the expression level of PTPRM could interpret the divergence of expression level between yaks and cattle. We quantified the expression level of PTPRM and revealed no significant change between the two species (data not shown), suggesting that the expression change of PTPRM is not a reason for the dramatic divergence of bta‐miR‐2285o‐3p between yaks and cattle. It was reported that SNPs in miRNA precursors may affect the biogenesis of mature miRNAs; thus we cloned and sequenced intronic sequences of PTPRM across 98 individuals for yaks and cattle. Two linked nucleotide polymorphic sites were identified in pre‐miR‐2285o as affecting the biogenesis of pre‐bta‐miR‐2285o. Overall, we revealed the involvement of two linked nucleotide polymorphic sites in bta‐miR‐2285o‐3p biogenesis. miRNA evolution shows a high rate of appearance and extinction, and lineage‐specific expansion, which are involved in the development of species‐specific phenotypes (Somel et al. 2011). The miR‐2285 family was reported to be specific to the ruminant lineage and to consist of more than 40 members (Marion et al. 2015). However, the functional roles of bta‐miR‐2285o‐3p remain to be elucidated.

Target gene prediction revealed 504 candidate target genes using targetscancow7.1, and functional enrichment analysis revealed the enrichment of hypoxia‐related terms for target genes. Thus, we performed gain‐ and loss‐of‐function analysis to reveal the function of miR‐2285o‐3p. The results showed that the transfection of miR‐2285o‐3p mimics attenuated hypoxia‐induced apoptosis by promoting cell survival, suppressing the reduction of cell viability and reducing cell membrane damage in vitro. We compared the 3′‐UTR sequences of all three candidate target genes across multiple species (human, mice and rats) and revealed high levels of conservation. Therefore, we performed a luciferase reporter assay in HeLa cells. Furthermore, a luciferase reporter assay revealed that VLDLR, PTEN and caspase‐3 were the authentic target genes of miR‐2285o‐3p. Expression level analysis in five yak tissues revealed the high negative Pearson’s correlation between bta‐miR‐2285o‐3p and caspase3, indicating the potentially important roles of miR‐2285o‐3p in vivo.

Conclusion

The bta‐miR‐2285 family encodes more than 40 mature miRNAs in the bovine genome, but the genes targeted by the miRNAs in this family are currently unknown. Here, we report the widespread expression of miR‐2285o‐3p in yaks but not cattle, and reveal the mechanism underlying this differential expression. Functional analysis indicated the potential important roles of miR‐2285o in high‐altitude adaptation in the yak. Our study provides a typical model for investigating miRNAs in adaptation biology.

Funding

This work was supported by grants from the National Natural Science Foundation of China (31872335, 31802044 and U19A2036) the Sichuan Science and Technology Program (2021YFYZ0009).

Conflict of interest

The authors declare no conflict of interest.

Author contributions

Conceptualization, M.L., K.L., A.J. and J.Z.; methodology, J.Z., J.W. and J.Z.; validation, J.Z. and J.W.; formal analysis, Q.T.; data curation, L.J.; writing – original draft preparation, J.Z. and K.L.; writing – review and editing, J.Z., K.L., J.W. and M.L.; visualization, X.L. All authors have read and agreed to the published version of the manuscript.

Ethical approval

All animal experiment used in this study were conducted according to the Institutional Animal Care and Use Committee in the College of Animal Science and Technology, Sichuan Agricultural University, Sichuan, China. All animal samples were purchased from local farmers.

Table S1. Information for primers used in this study.

Click here for additional data file.

Table S2. Genotype and allele frequencies of pre‐miR‐2285o polymorphisms in yaks and cattle.

Click here for additional data file.