key: cord-0026351-uaqfz80e
authors: Ma, Jun; Zhang, Tianliu; Wang, Wenxiang; Chen, Yan; Cai, Wentao; Zhu, Bo; Xu, Lingyang; Gao, Huijiang; Zhang, Lupei; Li, Junya; Gao, Xue
title: Comparative Transcriptome Analysis of Gayal (Bos frontalis), Yak (Bos grunniens), and Cattle (Bos taurus) Reveal the High-Altitude Adaptation
date: 2022-01-11
journal: Front Genet
DOI: 10.3389/fgene.2021.778788
sha: 9b18b96d8a4e6b12ec8d0a05a09b15bbd51a6504
doc_id: 26351
cord_uid: uaqfz80e

Gayal and yak are well adapted to the local high-altitude environments, yet the transcriptional regulation difference of the plateau environment among them remains obscure. Herein, cross-tissue and cross-species comparative transcriptome analysis were performed for the six hypoxia-sensitive tissues from gayal, yak, and cattle. Gene expression profiles for all single-copy orthologous genes showed tissue-specific expression patterns. By differential expression analysis, we identified 3020 and 1995 differentially expressed genes (DEGs) in at least one tissue of gayal vs. cattle and yak vs. cattle, respectively. Notably, we found that the adaptability of the gayal to the alpine canyon environment is highly similar to the yak living in the Qinghai-Tibet Plateau, such as promoting red blood cell development, angiogenesis, reducing blood coagulation, immune system activation, and energy metabolism shifts from fatty acid β-oxidation to glycolysis. By further analyzing the common and unique DEGs in the six tissues, we also found that numerous expression regulatory genes related to these functions are unique in the gayal and yak, which may play important roles in adapting to the corresponding high-altitude environment. Combined with WGCNA analysis, we found UQCRC1, COX5A are the shared differentially expression hub genes related to the energy supply of myocardial contraction in the heart-related modules of gayal and yak, and CAPS is a shared differentially hub gene among the hub genes of the lung-related module, which is related to pulmonary artery smooth muscle contraction. Additionally, EDN3 is the unique differentially expression hub gene related to the tracheal epithelium and pulmonary vasoconstriction in the lung of gayal. CHRM2 is a unique differentially expression hub gene that was identified in the heart of yak, which has an important role in the autonomous regulation of the heart. These results provide a basis for further understanding the complex transcriptome expression pattern and the regulatory mechanism of high-altitude domestication of gayal and yak.

Species living at high altitudes are exposed to strict selection pressures and physiological challenges owning to harsh environmental conditions, such as thin air, cold temperatures, ultraviolet exposure, and low pressure (Miao et al., 2015) . Despite the harsh conditions, numerous species including Tibetans (Ge et al., 2012; Tashi et al., 2014) , Tibetan pigs (Jia et al., 2016; Zhang B. et al., 2017) , Tibetan sheep , Tibetan chickens (Gou et al., 2014) , and yak (Qiu et al., 2012; Qi et al., 2019) have evolved unique physiological characteristics, such as superior blood oxygen transport systems and high metabolic efficiency, to adapt to the harsh living pressure (Lan et al., 2018) . At present, published studies have identified EPAS1, EGLN1, and PPARA genes that play important roles in highaltitude adaptation (Haasl and Payseur, 2016; Heinrich et al., 2019) . Exploring the molecular mechanisms underlying hypoxia adaptation has long garnered attention.

Gayal, also known as Drung cattle (Bos frontalis) (Winter et al., 1986; Uzzaman et al., 2014) , is a unique semi-wild cattle breed, mainly found in the typical alpine valleys and subtropical rain forests in the Drung and Nujiang river basins of Yunnan Province, southwestern China. Studies have shown that the level of hemoglobin in gayal, as well as the number of red blood cells and white blood cells, were equivalent to those of yak, which was the physiological basis for gayal to resist the invasion of bacteria, viruses, and parasites in the alpine environment (Tian et al., 1998) . Yak (Bos grunniens), a native mammal on the Qinghai-Tibet Plateau and its adjacent regions, provides meat and other necessities for Tibetans. Compared with lowland cattle, yaks have larger alveolar area per unit area, thinner alveolar spacing, thinner gas-blood barrier (Wei and Yu, 2008) , and stronger expression of some genes related to oxygen supply as well as defense system under hypoxia pressure (Wang et al., 2016) . Therefore, a thorough understanding of the genetic basis of the physiological characteristics in gayal and yak will provide insight into their adaptation to the high-altitude environment.

With recent rapid progress in next-generation genome sequencing (NGS) technologies, high-throughput RNAsequencing (RNA-seq) technology represents a powerful and cost-efficient approach to explore species domestication, breeding, and genetic molecular mechanisms in organisms (Salleh et al., 2017; Keren et al., 2018) . RNA-seq can detect and quantify gene expression with digital measurements, which are especially sensitive for low-expressed genes (Mortazavi et al., 2008) . In addition, it is also used to improve gene annotation, discover novel genes or transcripts, and survey alternative splicing (AS) events at single-nucleotide resolution (Sultan et al., 2008; Ozsolak and Milos, 2011) . Thus far, transcriptome studies on gayal and yak have been widely conducted. The muscle transcriptome analysis of Indian Mithun showed that hub genes including MTMR3, CUX1, LONRF3, PLXNB2, KMT2C, ZRSR2Y, PRR14, USP9Y, SLMAP, and KANSL2 might contribute to muscle growth (Mukherjee et al., 2020) . Transcriptomic analysis of yaks with different altitudes indicates PI3K-Akt, HIF-1, focal adhesion, and ECM-receptor interaction pathway were significantly enriched, and the EPAS1 expression increases with altitude (Qi et al., 2019) . Comparative transcriptome analysis of four organs (heart, kidney, liver, and lung) between yak and cattle revealed that DEGs associated with the blood supply system, modulation of cardiac contractility, vascular smooth muscle proliferation, and the glutamate receptor system probably play an important role in yak adaptation to hypoxia environments (Wang et al., 2016) . However, gayal is a species that is well adapted to the local alpine valley's environment, transcriptomic studies focusing on gayal adaptation to the local high-altitude environment have thus far not been reported. In addition, studies have found specific adaptive genes and mechanisms are distinct between species and populations (Heinrich et al., 2019; Pamenter et al., 2020) . Therefore, we want to know how changes in altitude affect the transcriptional regulation of gayal, and whether the relevant transcriptional regulation mechanism is the same as that of yak in the same bovine subfamily.

In this study, we characterized gene expression profiles using RNA-seq data from six major tissues including heart, lung, liver, kidney, spleen, and muscle among gayal, yak, and cattle. We performed a comparative transcriptome analysis of six tissues in gayal vs. cattle and yak vs. cattle, respectively, and identified candidate DEGs related to adaptability to the plateau environment, as well as analyzed their functions and expression patterns. In addition, weighted gene co-expression network analysis (WGCNA) was further used to uncover the hub genes that regulate tissue function in high-altitude adaptation. The results could lay a foundation for further explorations of the genetic changes underlying hilly adaptation in subtropical environment mammals.

Six tissues (heart, kidney, liver, lung, skeletal muscle, and spleen) were collected from three adult female gayal living in the semiwild preserved field at an altitude of 2300-3500 m in the Drung and Nujiang river basins. All tissue samples were taken and immediately snap-frozen in liquid nitrogen. In addition, the corresponding tissue transcriptome data of yak and cattle were retrieved from a previous publication (Tang et al., 2017) from the GEO databases under the accession numbers GSE93855 and GSE77020 (Supplementary Table S1 ). In addition, one gayal multi-tissue transcriptome data that we previously published was also included in this study to expand the sample size, which was stored in National Genomics Data Center under the accession code PRJCA002143. An overview of the samples used and the associated statistics are provided in Supplementary Table S1 .

Total RNA was isolated and purified with the TRIzol reagent (Life Technologies, Carlsbad, CA) according to the manufacturer's protocols. RNA purity was qualified using NanoDrop ND 2000 spectrophotometer at 260 and 280 nm (Thermo Fisher Scientific, Wilmington, DE), and RNA integrity was assessed using Agilent 2,100 Bioanalyzer (Agilent Technologies, Palo Alto, CA). The OD260/280 ratios of all samples were greater than 1.8, and the RNA integrity number (RIN) value of >7.0. Then, mRNA was broken down into small fragments using a magnesium RNA cleavage module (NEB, Ipswich, United States). Random primers and related reverse transcriptase were used to synthesize the firststrand complementary DNA (cDNA). Subsequently, the secondstrand cDNA was synthesized. The average insert size for the final cDNA library was 300 ± 50 bp. Finally, we conducted the 2 × 150 bp paired-end sequencing (PE150) on an Illumina Novaseq ™ 6000 (LC-Bio Technology CO., Ltd., Hangzhou, China) following the supplier's recommended protocol.

To identify orthologous genes among gayal, yak, and cattle, we first downloaded protein sequence data of yak and cattle from the Ensemble genome browser 104 (http://asia.ensembl.org/index. html), and combined with the protein sequence data of the gayal reference genome assembled by ourselves (the accession code of National Genomics Data Center is PRJCA004132) for downstream analysis. We used ORTHOFINDER v2.2.6 (Emms and Kelly, 2015) software to delimit orthologous groups, in which all predicted protein sequences were compared using a BLAST all-against-all search (Camacho et al., 2009 ). The single-copy genes, duplicated genes, and species-specific genes were extracted from the ORTHOFINDER output, respectively. Among them, the 1:1:1 homologous genes of the three species were defined as the conserved single-copy homologous genes of each species, which were retained for the following analysis.

To ensure the accuracy of subsequent biological information analysis, the quality of the raw data was first checked using FastQC (http://www.bioinformatics.babraham.ac.uk/projects/ fastqc/). According to the evaluation results of raw reads, the software FASTP (https://github.com/OpenGene/fastp) with default parameters was used to remove adapters, reads with a high proportion of unknown nucleotides, and low-quality sequences (≤Q20) to obtain high-quality clean reads. To avoid biases due to sequence divergence across bovine subfamily genomes, we downloaded yak (BosGru_v2.0, http://ftp. ensembl.org/pub/release-104/fasta/bos_mutus/) and cattle (ARS_UCD1_2, http://ftp.ensem -bl.org/pub/release-104/fasta/ bos_taurus/) reference genomes from the Ensemble database, and retrieved gayal (Drung_v1.2) reference genomes from National Genomics Data Center (NGDC) (https://bigd.big.ac. cn/). Then, the clean reads for each sample were mapped to the corresponding reference genome using HISAT2 software (Pertea et al., 2016) . FeatureCount from the Subread package (Version 2.0.0) (Liao et al., 2014) was applied to calculate the numbers of reads mapped to each gene. Read pairs were counted only if they uniquely mapped and were properly paired. The gene expression level estimates may be biased among species due to factors such as the size of mRNA transcript and the qualities of different genome annotations (Blake et al., 2020) . To evade these issues, we used a custom script to only retain reads that mapped to the 1: 1 orthologous which can be used for each of the three genomes (Blake et al., 2020) . Taking into account the difference in sequencing depth of samples from different sequencing platforms, we first filtered out genes with low expression, and only keep genes with expression level CPM> 1 in at least half of the samples (Blake et al., 2020) . We normalized the raw read counts using the weighted trimmed mean of M-values algorithm (TMM) in the edgeR package (Robinson et al., 2010) . Then scaling factors were used to scale the expression levels of all orthologous. The normalized expression data for the trimmed orthologous were used for subsequent gene expression analysis.

After performing normalization, we obtained the TMMnormalized log 2 -transformed CPM values for each gene. The principal component analysis (PCA) was performed with the prcomp function in R using normalized log 2 -transformed CPM values. The hierarchical clustering of Pearson's correlation coefficients for each pair of samples was performed using the complete-linkage agglomerative method on the correlation distance matrix and generated symmetrical heatmap plot using the R gplots package.

To characterize the high-altitude adaptive gene expression pattern of gayal and yak, we performed differential expressed analyses among heart, kidney, liver, lung, skeletal muscle, and spleen tissues in gayal vs. cattle group and yak vs. cattle group by using the R package edgeR (Robinson et al., 2010) . Genes with Benjamini Hochberg false discovery rate (FDR) values less than 0.05 and fold change more than 2 were identified as DEGs. If the log 2 fold changes more than 1 they were defined as up-regulated DEGs, or less than -1 they were defined as down-regulated. We performed Gene Ontology (GO) function annotation and KEGG pathway enrichment analysis for DEGs by the Database for Annotation, Visualization and Integrated Discovery (DAVID) v6.8 using the Bos taurus genome as the background set (Huang Da et al., 2009) . The significantly enriched GO and KEGG terms using a 0.05 cutoff for the p-value with Benjamini-Hochberg test adjustment. The jvenn and UpSetR R package (version 1.4.0) were used to identify and visualize the DEGs shared by two comparison groups across the multiple tissues and between two comparison groups in the same tissue, respectively (Conway et al., 2017) . In addition, to better understand the hypoxia adaptation of species, we also collected a set of known or potential high-altitude adaptation candidate genes from previously published literature Wang et al., 2016; Peng et al., 2017) .

The adaptation of species to the hypoxia environment is a complex biological process, which involves the coordinated regulation of multiple genes. WGCNA is a feasible method to clarify the regulatory relationship between genes and identify important hub genes based on gene expression profiles to explore the mechanism of hypoxia adaptation (Wang et al., 2021) . Therefore, to further identify the relevant regulatory gene modules related to the functions of each tissue and explore the core driver genes, we used the CPM values of filtered genes to perform WGCNA analysis on three species, respectively (Langfelder and Horvath, 2008) . The hierarchical cluster analysis on all tissues was carried out using the hclust function. The soft-power threshold beta was determined by the function "sft$powerEstimate" . The minimal gene module size was set to 30 to obtain appropriate modules, and the threshold to merge similar modules was set to 0.25. Gene modules were detected based on the TOM matrix. The correlations between modules eigengene and tissues were investigated using the Pearson correlation, and only modules with correlation coefficients >0.7 and p < 0.001 were considered tissue-relevant modules (Touma et al., 2016) . Here, we designate the highest top 5% of genes as the hub genes based on the intramodular connectivity in each functional module, which is calculated by summing the connection strengths with other intramodular genes (Yang et al., 2014) . Functional enrichment analyses for tissue-related module genes were performed using DAVID v6.8 using the Bos taurus genome as the background set, and the p-values were corrected for multiple testing using Benjamini-Hochberg correction (Huang Da et al., 2009) .

The transcriptome data of six tissues (heart, kidney, liver, lung, skeletal muscle, and spleen) from four gayals, three yaks, and three cattle were used in this study (Supplementary Table S1 ). One gayal heart sample (gayal3_heart) was excluded due to low sequencing quality. In total, ∼327.55G clean data has been generated after removing adapter-embedded or low-quality reads. We calculated Q20 (99% base call accuracy), Q30 (99.9% base call accuracy), and GC-content of the clean data.

The average values of Q20 and Q30 in the samples were 97.04% (ranging from 94.60 to 98.47%) and 92.38% (ranging from 89.39 to 95.43%), respectively (Supplementary Table S2 ). The average base GC content was 50.80%, ranging from 47.69 to 55.03%, indicating that the base composition and quality of the sequencing data were qualified. Then the high-quality reads of each sample were mapped to the corresponding reference genome (i.e., Drung_v1.2, BosGru_v2.0, ARS_UCD1_2) using the HISAT2 software (Pertea et al., 2016) . The average reads mapping rates in gayal, yak, and cattle were 86.35, 93.60, and 95.58%, respectively. Detailed statistics and mapping information for each sample were summarized in Supplementary Table S2 .

To reduce the effect of chromosome number difference (Bos frontalis, 2n 58; and Bos grunniens, 2n 60; Bos taurus, 2n 60) and the qualities of different genome annotations among three bovine species, we first used the blast function of ORTHOFINDER2 (Emms and Kelly, 2015) to identify orthologous genes among the gayal, yak, and cattle. A total of 19,174 protein orthogroups were identified. Among them, only 6371 orthogroups were high-confidence single-copy orthogroups of three species (Supplementary Table S3 ). Then the corresponding gene names were extracted by a custom script, and retained reads that mapped to the single-copy homologous genes of three species genome, and the ensemble gene names of cattle were uniformly used for subsequent analysis. In the 6371 single-copy orthologous gene expression datasets, genes with low expression (CPM <1) in more than half of the samples were further excluded, and 4712 genes were obtained (Blake et al., 2020) . Then we normalized 4712 single-copy orthologous genes using the weighted TMM of edgeR (Hao et al., 2019) .

To investigate broad gene expression patterns in multiple tissues across different bovine species, we performed hierarchical cluster analysis based on the gene expression of six tissues among gayal, yak, and cattle. As shown in Figure 1A , the tissue-specific expression pattern of lung, liver, spleen, and kidney shows that the same tissues of different bovine species gather earlier than different tissues of the same bovine species. The result indicated that tissue differentiation may precede species differentiation, which was consistent with previous studies (Necsulea and Kaessmann, 2014; Hao et al., 2019) . However, the samples of heart and muscle from each species were more similar than the same tissue between different species (speciesspecific expression pattern), which reflects striated muscular tissues may have similar expression profiles (Merkin et al., 2012; Tang et al., 2017; Feng et al., 2020 ) ( Figure 1A ). In addition, our expression profiles suggested that the differences in global gene expression among tissues were more significant than that among altitudes or the same tissue of different species, which was consistent with the results of previous studies (Necsulea and Kaessmann, 2014; Tang et al., 2017; Hao et al., 2019) . Furthermore, principal component analysis (PCA) was performed on the gene expression of 59 samples. The score plot of PCA analysis for six main organs shows that 45.4% of the variance can be explained by the first two principal components, accounting for 26.4 and 19% of the variance, respectively ( Figure 1B) . Basically, 59 samples were clustered based on similar tissue types, which confirmed the tissue differentiation may precede species differentiation.

To detect expression changes of gayal and yak in response to the high-altitude stress, we used edgeR (Robinson et al., 2010) to identify DEGs within each tissue. According to the tissue and species we considered, 1218 to 1585 interspecies DEGs were Frontiers in Genetics | www.frontiersin.org January 2022 | Volume 12 | Article 778788 4 identified in the comparison of six types of tissues in gayal and cattle ( Table 1) . Among the six tissues, the heart has the largest number of DEGs (618 up-regulated and 967 down-regulated genes). By upset analysis of all DEGs in each tissue, 446 DEGs shared in six tissues were identified, and each gene showed consistent up-regulation or down-regulation in all tissues ( Figure 2A ). In the yak vs. cattle group, we identified 1995 DEGs among six types of tissues. The number of DEGs among tissues ranges from 543 in liver tissue (202 upregulated and 341 down-regulated genes) to 919 in lung tissue (387 up-regulated and 532 down-regulated genes) ( Table 1) . There were 154 DEGs with consistent regulation direction shared among six tissues, of which 64 genes were also present in all tissues of the gayal vs. cattle group ( Figure 2B) .

To further investigate the potential biological function of DEGs for environmental adaptation in gayal and yak, we performed functional annotation analysis based on the DAVID database. In the gayal vs. cattle group, DEGs were annotated in energy metabolism, cardiovascular system, blood coagulation system, and immune system among tissues (Supplementary Table S4 ). Such as the DEGs of the heart were enriched in the oxidation-reduction process (GO: 0055114), positive regulation of angiogenesis (GO:0045766), and positive regulation of blood coagulation (GO:0030194). Muscle DEGs were enriched in blood coagulation (GO: 0007596), patterning of blood vessels (GO:0001569), and skeletal muscle cell differentiation (GO:0035914) ( Figure 3) . We also observed DEGs in the gayal vs. cattle group were significantly enriched in hypoxia response, metabolism, and immune system pathway. For instance, hematopoietic cell lineage (bta04640) was enriched in four tissues (heart, spleen, lung, and kidney) and galactose metabolism (bta00052) was enriched in four tissues (heart, liver, spleen, and lung) (Supplementary Table S4 ). In the yak vs. cattle group, DEGs identified across six tissues pairwise comparisons were also enriched in the cardiovascular system, energy metabolism, blood coagulation system, such as response to oxidative stress (GO:0006979) and oxidation-reduction process (GO:0055114) in the heart, and patterning of blood vessels (GO:0001569), and positive regulation of canonical wnt signaling pathway (GO: 0090263) in the muscle (Figure 4) . Similar results were found in the KEGG analysis, pathways related to energy metabolism, immunity system were significantly enriched (Supplementary Table S4 ).

To explore the similarities and differences of expression regulation genes related to environmental adaptability among gayal and yak, we further analyzed the shared and unique DEGs between the gayal vs. cattle group and the yak vs. cattle group. Among the six tissues of gayal and yak, the heart and lung shared the largest number of DEGs, with 397 (101 up-regulated and 296 down-regulated genes) and 396 genes (146 upregulated and 250 down-regulated genes), respectively ( Figure 5 ). In addition, we have identified 822-1331 and 298-500 unique DEGs in the six tissues of gayal and yak Table S4 ). Frontiers in Genetics | www.frontiersin.org January 2022 | Volume 12 | Article 778788 8

As for the unique DEGs of the gayal vs. cattle group, we found five unique DEGs in the heart related to apoptosis (bta04210), including CASP9, ENDOG, AIFM1, CASP6, TNFSF10, RIPK1, FASLG, and TNF. Seven unique genes (BLOC1S6, SERPINE2, TFPI2, HPS4, SERPING1, GP1BA, TFPI) are associated with blood coagulation (GO:0007596) in the lung. Spleen-unique DEGs were enriched in lymphocyte chemotaxis (GO:0048247), inflammatory response (GO:0006954). In addition, we found that positive regulation of angiogenesis (GO:0045766) was significantly enriched in specific genes in the lung, kidney, and spleen; hematopoietic cell lineage (bta04640) was also significantly enriched in the unique genes of the heart, spleen, lung, and kidney. In addition, we found that the unique DEG in the heart, liver, and spleen was also related to galactose metabolism (bta00052). These unique DEGs may play an important role in the adaptation of gayal to the alpine and valley environment (Supplementary Table S4 ). Among the unique DEGs of the yak vs. cattle group, we found that four genes (CD40LG, SRF, F2R, F5) related to platelet activation (GO: 0030168) in the heart. Blood coagulation (GO:0007596) related genes (PROCR, THBD, F10, PROS1, F3, KNG1) were also significantly enriched in the lung. The unique genes in the spleen were enriched in immune-related items, such as NODlike receptor signaling pathway (bta04621), and TNF signaling pathway (bta04668). In addition, cell adhesion molecules (CAMs) (bta04514) were uniquely enriched in the heart and lung of yak, which was reported as the intensity of hypoxia increases lead to marked reduced cell adhesion (Kaiser et al., 2012) (Supplementary Table S4 ).

To further establish the hypoxia adaptation of gayal and yak, we collected a candidate gene set with known or potential functions in hypoxia responses from previously published literature Peng et al., 2017) . Of those, 295 and 187 genes were DEGs in at least one tissue of gayal and yak, respectively (Supplementary Table S5 ). Among which the heart has the largest number of high-altitude related genes in both comparison groups DEGs (149 in gayal vs. cattle group and 85 genes in yak vs. cattle group). Notably, in the altitude-related Frontiers in Genetics | www.frontiersin.org January 2022 | Volume 12 | Article 778788 9 gene set, we found that several positively selected genes related to hypoxia response in Tibetans were differentially expressed in two comparison groups. For example, EGLN1, also known as PHD2, as a known target gene of the hypoxia-sensing pathway (Haasl and Payseur, 2016; Heinrich et al., 2019; Storz, 2021) , showed a consistent down-regulation across all tissues of yak and gayal compared with lowland cattle, but the expression level and pattern vary in different tissues. HFE (homeostatic iron regulator), a gene involved in indirectly regulating iron balance in animals via reducing transferrin-mediated iron uptake (Salter-Cid et al., 1999) , was uniquely down-regulation in all tissues of gayal except for kidney compared with the lowland cattle. Additionally, HOXB6, a gene associated with the regulation of the erythropoietic system, was particularly up-regulated in the spleen and kidney of the yak (Zimmermann and Rich, 1997) .

To explore genes and co-expression modules with similar expression profiles related to high-altitude adaptation in various tissues, WGCNA analysis was independently performed on three species based on the CPM values of 4712 filtered genes. The sample clustering dendrogram of the gayal, yak, and cattle are shown in Figure 6A , Figure 7A , and Figure 8A , respectively. In both species, all samples were included in the clusters. The soft-power threshold of 10, 18, and 16 (scale free R 2 of 0.85) were selected for further analysis of the gayal, yak, and cattle, respectively ( Figure 6B , Figure 7B , and Figure 8B ). Next, the gene modules were detected based on the topological overlap matrix (TOM), and 16 modules were detected under the gayal, while 10 modules were detected under the yak ( Figure 6C , Figure 7C , and Figure 8C ). To identify tissuespecific modules, we calculated the correlation coefficients between the module and each tissue ( Figure 6D , Figure 7D , and Figure 8D ). Further, GO and KEGG analyses were conducted to investigate the biological function of the genes in each tissue-related module.

Indeed, the relevant modules enriched in 6 tissues are highly related to the tissue-specific functions. For example, in the tissuerelated modules of gayal, the MEcyan module was significantly correlated with heart (r 0.88, p 3 × 10 −8 ), which was enriched in functional annotations such as negative regulation of cardiac muscle cell apoptotic process (GO:0010667), and pathway in MEcyan module were involved in cardiac muscle contraction (bta04260). Lung-related module (MEgreenyellow, r 0.99, p 4 × 10 −18 ) was enriched in HIF-1 signaling pathway (bta04066) (Supplementary Table S6 ). Similar results were found in the yak, MEgreenyellow was the most significantly correlated with the heart of yak (r 0.98, p 3 × 10 −12 ), which was enriched in functional annotations such as positive regulation of cardiac muscle cell proliferation (GO:0060045), positive regulation of Table S6 ). For low-altitude cattle, heart-related module (MElightyellow, r 0.99, p 1 × 10 −15 ) was enriched in cardiac right ventricle morphogenesis (GO:0003215), negative regulation of cardiac muscle cell apoptotic process (GO:0010667), and heart development (GO: 0007507). Lung-related module (MEdarkorange, r 0.99, p 7 × 10 −16 ) was enriched in lung development (GO:0030324). Additionally, in the tissue-related module of each species, the terms related to bile secretion were found in the liver-specific module, and functions related to ion transport and absorption are enriched in the kidney-specific module, and immune and skeletal muscle cell differentiation related functions were enriched in the spleen-specific and muscle-specific module, respectively (Supplementary Table S6) . Furthermore, genes with the highest intramodular connectivity were selected as hub genes in each tissue-related module in both gayal and yak. We identified 4-51, 6-44, and 2-46 hub genes in 6 tissue-related modules of gayal, yak, and cattle, respectively (Supplementary Table S7 ). Among the heart-related modules of the three species, there are 14 common hub genes between gayal and yak, of which 5 genes are involved in cardiac contraction (bta04260), including UQCRC1, UQCRQ, UQCRFS1, UQCR11, and COX5A. Among them, UQCRC1, UQCRQ, UQCRFS1, and UQCR11 are important subunits of mitochondrial complex III (Fernández-Vizarra and Zeviani, 2015) . Mitochondria play a key role in cardioprotection and heart contraction energy needed for pumping blood to oxygenate the body organs (Lemieux and Hoppel, 2009) . UQCRC1 directly affect mitochondria function related to cardioprotection (Yi et al., 2020) . COX5A is the component of the cytochrome c oxidase, which is central to oxidative phosphorylation and ATP generation (Huang et al., 2019) . However, no common hub gene was found between gayal and low-altitude cattle, and only four hub genes are shared by yak and cattle, including KCNJ8, MTUS2, P4HTM, and BTG3. KCNJ8 encodes a subunit of an ATP-sensitive potassium channel and regulates potassium inward rectifier, which is related to heart development (Erginel-Unaltuna et al., 1998) . There are 1 (COX6A2) and 2 (CYC1 and UQCR10) non-common hub genes in the heart-related modules of gayal and yak, respectively, which are involved in cardiac contraction. In addition, 4 hub genes (HAND2, CHRM2, EDNRA, and ZFPM2) associated with heart development and regulation were uniquely identified in the greenyellow module of the yak heart. Among these, CHRM2 was an up-regulated DEG in the yak heart (log 2 fold changes of 1.14), which plays a fundamental role in autonomic regulation of the heart, such as heart rate recovery after maximal exercise (Hautala et al., 2006) . HAND2 is essential for cardiac morphogenesis, vascular development, and the regulation of angiogenesis (Anderson et al., 2016) . EDNRA encodes the endothelin 1 (EDN1) receptor, a gene involved in the vasoconstriction mechanism (Sharma et al., 2014) , as well as closely related to pulmonary hypertension and HIF activity, which is positively selected in Tibetan (Simonson et al., 2010) . ZFPM2, regulates the activity of GATA family proteins, plays a vital function in cardiac morphogenesis and coronary vascular development (Luo et al., 2020) . Among the hub genes of the lung-related module, 6 common hub genes were identified among gayal and yak, including CAPS, TSNAXIP1, GPRIN2, SPEF1, BMP3, and EFHC1. Gayal and yak have 2 (BMP3 and SEC14L3) and 3 (BMP3, DNAAF1, and N4BP3) common hub genes with low-altitude cattle, respectively. Among them, BMP3, a hub gene shared by three bovine subfamily species, has been reported to play a regulatory role in the morphogenesis and/or function of the human lung (Vukicevic et al., 1994) . CAPS encodes a calcium-binding protein and plays a role in the regulation of Ca 2+ transport (Pan et al., 2019) , which may be correlated with pulmonary artery smooth muscle contraction (Waypa and Schumacker, 2002) . SEC14L3 is an up-regulated DEG in the gayal lung and plays an important role in maintaining the homeostasis of airway epithelial cells (Shan et al., 2009) . DNAAF1 gene is associated with lung development (Duquesnoy et al., 2009) . As for the noncommon hub genes in the lung-related module of gayal, we discovered several genes (EDN3, GJA5, and CCDC39) related to the tracheal epithelium and pulmonary vasoconstriction. Such as, EDN3, an endothelium-derived vasoconstrictor peptide, which was up-regulated in the lungs of gayal and may play an important role in enhancing pulmonary vasoconstrictor response (Masaki, 2000) . GJA5 is the predominant gap junction protein present in vascular endothelium, plays an important role in coupling between cells in the vascular wall (Ebong et al., 2006) . CCDC39 was related to lung development (Blanchon et al., 2012) . Therefore, these genes related to promoting cardiac blood circulation and pulmonary vasoconstriction may play crucial roles in the environmental adaptability of gayal and yak.

In this study, we explored the high-altitude adaptability of Bovidae using the transcriptome data of 6 tissues of gayal and published data of yak and lowland cattle (Tang et al., 2017) . Compared with numerous previous studies that only used the same reference genomes (yak or cattle genome) for different bovine subfamily species (Tang et al., 2017; Lan et al., 2018; Xin et al., 2019; Wang et al., 2020) , here we optimized the experimental design, using the corresponding reference genomes for 3 species to reduce the gene sequence differences of different species due to phylogenetic evolution and only retain reads that mapped to the 1:1 orthologous gene which had highconfidence alignments across the 3 genomes for subsequent analysis. To reduce the difference in the sequencing depth of transcripts from different sequencing platforms, we used the TMM standardized method of edgeR for correction with reference to the method of Blake et al. (2020) . Consistent with previous studies (Necsulea and Kaessmann, 2014; Tang et al., 2017; Hao et al., 2019) , we found that the differences in global gene expression among tissues were more significant than that among altitudes or the same tissue of different species. Among them, 4 tissues (lung, liver, spleen, kidney) showed tissue-specific expression patterns, while bovine species-specific expression patterns were shown in the tissue of heart and muscle (Necsulea and Kaessmann, 2014; Tang et al., 2017; Hao et al., 2019) . In addition, we identified numerous DEGs in gayal and yak compared with low-altitude cattle, which might have an important role in high altitude adaptation. In addition, WGCNA was further used to explore the core genes of each tissue that regulate tissue function in the high-altitude adaptation.

Hypoxia is one of the most serious challenges faced by highaltitude animals, and oxygen supply is directly related to the development of red blood cells (Zhang D. et al., 2017) . Studies in mice (Li et al., 2011) and yak (Xin et al., 2019) have shown an increase in the number of red blood cells and platelets under hypoxia, but the number of granulocyte-macrophage progenitor cells periodically declined. In this study, similar results were found in the gayal. Hematopoietic cell lineage pathway was uniquely enriched in multiple tissues (heart, kidney, lung, spleen) of gayal vs. cattle group. Among these, IL11RA is a member of the hematopoietic cytokine receptor family (Ng et al., 2021) , which was uniquely up-regulated by 1.33 times in the heart of gayal. GP1BA is the alpha subunit of platelet surface membrane glycoprotein Ib (Luo et al., 2007) , were 9.58, 3.90, 7.77, 11.41, 12.38, and 12.39-fold higher in the tissues (corresponding to heart, spleen, lung, kidney) of gayal than in the cattle, suggesting that the number of red blood cells and platelets might be increased. As for the regulatory mechanisms underlying granulocyte-macrophage progenitor cells, we observed three genes (IL1B, CSF1, IL6R) were down-regulated in multi tissues of gayal, which might be a reason for the decreased proportion of granulocytes. Among them, IL1B is produced by activated macrophages as a proprotein (Lappas, 2013) , which was down-regulated 1.94-fold in the lung of gayal. CSF1 controls the production and differentiation, and function of macrophages (Hamilton et al., 2016) , displaying log 2 fold changes of -1.06, -1.10, and -1.76 in the heart, spleen, and muscle, respectively. In addition, Interleukin 6 is involved in the maturation of B cells. IL6R encodes a subunit of the interleukin 6 (IL6) receptor complex and was found down-regulated in all tissues of gayal, which suggests lymphocytes may decrease, which is consistent with the previous results on the hypoxia adaptation of yak (Xin et al., 2019) . In addition, in the high-altitude related gene set that we searched (Supplementary Table S5 ), several crucial positive selection genes related to erythrocyte development were found differential expression in the gayal and yak. EGLN1, a key positive selection gene in the HIF pathway (Bigham et al., 2010) , was found consistently down-regulated in all tissues of gayal, which contribute to HIF1-α activity and has a role in increased Hb levels via erythropoiesis under hypoxia conditions (Epstein et al., 2001; Xiang et al., 2013) . HFE is a positively selected gene in Tibetans (Yi et al., 2010) , which indirectly regulates iron balance by reducing the entry of transferrin into cells (Salter-Cid et al., 1999) . Compared with lowland cattle, the appropriate downregulation of HFE gene expression in gayal can accelerate iron storage to facilitate the synthesis of hemoglobin to cope with the reduced oxygen partial pressure (Muckenthaler et al., 2020) . In addition, in the high-altitude related DEG set of yak, we also found EGLN1 was consistently down-regulated in all tissues of yak. In addition, HOXB6 is considered to be a marker for erythropoiesis (Magli et al., 1997) , and its expression was uniquely up-regulated 1.85-and 1.11-fold in the spleen and kidney of the yak, respectively, which may contribute to yak adaptation to the hypoxia environment (Xin et al., 2019) .

The regulation of angiogenesis by hypoxia links blood oxygen supply to metabolism (Pugh and Ratcliffe, 2003) . In this study, the shared and unique DEGs involved in blood vessel and vascular development were identified. For example, glutamyl aminopeptidase (ENPEP), a member of the M1 family of endopeptidases involved in blood pressure regulation and blood vessel formation (Holmes et al., 2017) , was up-regulated in the lungs of gayal and yak with the log 2 fold changes of 2.78 and 2.98, respectively. Additionally, we found that SRPX2 and ESM1 were up-regulated in the spleen of gayal and yak relative to lowland cattle, and SAT1 and SERPINE1 genes were up-regulated in the kidney and muscle of gayal and yak relative to lowland cattle, respectively, which are angiogenesis promoting genes (Aitkenhead et al., 2002; Mohr et al., 2017; Xiao et al., 2020) . HEY1 is the downstream effector of notch signaling required for cardiovascular development (Aitkenhead et al., 2002) with log 2 fold changes of 1.08 and 1.71 in the heart and lung of gayal, respectively, while displaying log 2 fold changes of 1.25 and 1.03 in the yak, respectively. In addition, MEOX2 has been reported to have a role in inhibited angiogenesis (Dhahri et al., 2020) , was down-regulated in all 6 tissues of gayal and the kidney of yak. Study has found VHL-deficient mice show increased vasculogenesis in embryonic cells , and VHL gene was down-regulated expression in all six tissues of yak and down-regulated in the lung of gayal. In addition, we have identified two unique DEGs (RAMP2 and TGFBI) in gayal that promote angiogenesis. RAMP2 mediates the pro-angiogenic effect (Albertin et al., 2010) , which was up-regulated in all 6 tissues of gayal. TGFBI is a matricellular protein coding gene that plays an important role in tumor angiogenesis (Fico and Santamaria-Martínez, 2020) , which was up-regulated in the liver and kidney of gayal with log 2 fold changes of 1.34 and 1.45, respectively. Additionally, DEGs (THY1 and PLCD3) that promote angiogenesis were also uniquely found in yak. THY1 was expressed on mouse tumor-associated lymphatic vessels and blood vessels (Jurisic et al., 2010) , displaying log 2 fold changes of 1.63 and 1.03 in spleen and lung of yak. PLCD3 induce angiogenesis in human endothelium (Kim et al., 2002) , which was up-regulated in all 6 tissues, displaying log 2 fold changes of 1. 77, 1.83, 2.07, 3.20, 2.12, and 1.15 in the heart, liver, spleen, lung, kidney, and muscle of yak, respectively. This suggests that angiogenesis may be conducive to the adaptation of gayal and yak to high-altitude environment.

Increasing evidence shows that the harsh environment of hypoxia may regulate the immune system (Mishra and Ganju, 2010; Gaur et al., 2020) . In the present study, we found that the immune system of the spleen in gayal and yak was highly regulated. DEGs of the spleen in gayal were significantly enriched in lymphocytic chemotaxis, inflammatory response, and cytokine receptor interaction pathways (Supplementary Table S4 ), and the DEGs of the spleen in yak were also significantly enriched in TNF signaling pathway, cytokine receptor interaction, and chemokine signaling pathway (Supplementary Table S4 ). In the spleen of gayal vs. cattle group, chemokines and their receptors like CCL4, CCL5, CCL14, CCL16, and CXCR4 displayed log 2 fold changes of 2.33, 2.78, 1.24, 1.83, and 1.18, respectively, which mediates chemokinetic inflammatory response (Loetscher et al., 1994; Ariel et al., 2000) . Interleukins like IL12RB2 promotes the proliferation of T-cells as well as NK cells (Presky et al., 1996) , IL23A contribute to the production of proinflammatory cytokines (Hoeve et al., 2006) , and IL12RB2 and IL23A with log 2 fold changes of 1.77 and 1.09, respectively; tumor necrosis factors like TNFRSF8, TNFRSF4, and LTBR displayed log 2 fold changes of 4.06, 1.41, and 1.32, respectively, which activated T and B cells expression (Nishikori et al., 2005) . CD40LG regulates B cell function by engaging CD40 on the B cell surface (Henn et al., 1998) , which was up-regulated 2.63-fold in the spleen of gayal. In the spleen of the yak vs. cattle group, CXCR4 and CD40LG were also up-regulated, with log 2 fold changes of 1.52 and 2.83, respectively. CXCL14 which fulfills a unique role in antimicrobial immunity was up-regulated 2.15-fold in the spleen of yak (Maerki et al., 2009 ). These changes might help gayal and yak resist inflammation and disease in highaltitude harsh environments.

High altitude pulmonary edema (HAPE) is a serious lifethreatening disease in humans and animals characterized by uneven vasoconstriction of pulmonary blood vessels due to hypoxia (Stream and Grissom, 2008) . It has been confirmed that HAPE is related to coagulation activation and fibrin formation enhancement (Mannucci et al., 2002) . In addition, hypoxia can lead to increased platelet activity, significantly increasing the risk of thrombosis (Shilei et al., 2016) . In this study, we found that a large number of DEGs are related to platelet activation and blood coagulation in the two comparison groups. Coagulation factors play an important role in the endogenous and exogenous blood coagulation process in animals (Thiel et al., 2017) . In this study, compared with lowaltitude cattle, the expression levels of coagulation factor II (F2) and coagulation factor V (F5) in the lung of gayal were separately reduced 1.26-and 1.64-fold, and displayed log 2 fold changes of -2.33 and -1.23 in the lung of yak, respectively. In addition, coagulation factor X (F10) was found uniquely down-regulated 5.62 times in the lungs of yak. Meanwhile, serine protease inhibitors are inhibitors of blood coagulation factors (Seixas et al., 2007) , and we have found the expression level of SERPINA5, SERPINF2, and SERPING1 were uniquely increased in the lung of gayal, and the corresponding log 2 fold changes of 1.80, 3.06, and 1.48, which may play an important role in the mechanism of pulmonary edema resistance, thus reducing the risk of pulmonary edema. In addition, impaired decomposition of bradykinin and its metabolites may be an important cause of pulmonary edema (De Maat et al., 2020) . We found that the expression of the kininogen 1 (KNG1) gene was down-regulated by 2.46 times in the lung of yak, while downregulated by 12.39-and 7.38-fold in the heart of gayal and yak than that of lowland cattle, respectively. These results suggest yak and gayal have developed a mechanism to prevent the occurrence of severe pulmonary edema.

High energy metabolism or reduced energy turnover is essential for survival under hypoxia (Qiu et al., 2012) . In this study, mitochondrion-related categories were widely enriched in the tissues of gayal and yak (Supplementary Table S4 ). Mitochondria is the main place for aerobic respiration of cells, which provides ATP for organisms through the oxidation of sugar, fats, and amino acids (Solaini et al., 2010) . Similarly, several pathways related to carbohydrate and fat energy metabolism were enriched, such as glycolysis/gluconeogenesis in yak lung, oxidative phosphorylation, and galactose metabolism in multiple tissues of gayal, and PPAR signaling pathway in the heart of gayal and yak. Among the genes related to these functions, TPI1 is essential for efficient energy production and plays an important role in glycolysis (Shimoda et al., 2012) . G6PC (glucose-6-phosphatase) is a key enzyme of gluconeogenesis (Jia et al., 2012) . We found that the gene expression of TPI1 was upregulated in all tissues of both gayal and yak, and G6PC was down-regulated in three tissues (heart, spleen, and lung) of both gayal and yak, which is conducive to the production of ATP through strengthened glycolysis and reduced gluconeogenesis. Oxidative phosphorylation is an important process for ATP synthesize in aerobic organisms. We found that most genes related to oxidative phosphorylation were up-regulated in gayal, such as COX6A2, NDUFA9, COX8B, NDUFB7, NDUFA11, NDUFB5, UQCR10, NDUFS8, UQCRQ, PPA1, and NDUFS3. In addition, the significant down-regulation of COX5B found in all tissues of gayal and yak may be beneficial to hypoxia adaptation (Trueblood et al., 1988) . Additionally, galactose metabolism (bta00052) was uniquely enriched in the tissues (lung, liver, spleen, and muscle) of gayal. Galactose is a kind of monosaccharide, most of which are converted into glucose in the liver, and then incorporated into glycogen or used for energy metabolism by glycolysis. Lactase (LCT), a key enzyme involved in the conversion of galactose to glucose was up-regulated in all 6 tissues of gayal. The PPAR signaling pathway, a classic pathway in lipid catabolism was significantly enriched in the heart of gayal and yak. EHHADH, a dehydrogenase in the β-oxidase system (Yeh et al., 2006) , was down-regulated 4.02 and 2.69 times in the heart of gayal and yaks. Additionally, the utilization of fatty acids first needs to be activated to produce fatty acyl CoA under the catalysis of the ACSL and SLC27A family (Bowman et al., 2016) . In this study, we found that SLC27A2 and SLC27A6 of the SLC27A gene family were down-regulated 6.24 and 3.24 times in the heart of gayal, and SLC27A2 and SLC27A4 were downregulated 6.81 and 2.08 times in yak heart, respectively. In addition, the ACSL3 and ACSL4 genes expression of ACSL gene family members were found uniquely decreased in the heart of gayal. Therefore, our results suggest gayal and yak may adapt to high altitude environments by strengthening carbohydrate metabolism and reducing fatty acid degradation although the related regulation genes may be different between gayal and yak.

Despite we found numerous genes that may be related to highaltitude adaptation of gayal and yak, we must admit the limitation that age differences between gayal and other species may affect the experimental results obtained to some extent. Therefore, it is necessary to conduct comparative transcriptome analysis on more bovine species of the same age that are combined with physiologic experiments to verify our findings, and better clarify the high-altitude adaptation mechanism of bovine subfamily species.

This study shows a comprehensive multi tissues transcriptome expression landscape of high-and lowaltitude bovine populations, and reveals the gene expression regulation associated with high-altitude acclimatization. The gene expression profiles of gayal, yak, and cattle showed tissue-specific expression patterns. By the comparative analysis of 6 hypoxia-related tissues of high-and low-altitude bovine subfamily species, which highlights numerous DEGs and terms underlying associated with hypoxia. Notably, categories related to angiogenesis, blood coagulation, energy metabolism, and immunity were commonly enriched in the tissue DEGs of gayal vs. cattle group and yak vs. cattle group, and we found that many expression regulatory genes related to these functions in gayal and yak are different, which may serve as an important regulatory mechanism for gayal and yak to adapt to the corresponding local environment. In addition, we also found numerous hub genes related to myocardial contraction, energy metabolism, and pulmonary vasoconstriction by WGCNA analysis. Overall, our study lays a foundation for further study on the environmental adaptability of mammals in subtropical plateau and provides a valuable basis for resource utilization in high-altitude areas.

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: NCBI PRJNA783860.

The animal study was reviewed and approved by the Scientific Research Department of the Institute of Animal Science, Chinese Academy of Agricultural Sciences, and strictly followed the guidelines formulated by the Chinese Animal Protection Commission.

JM, TZ, and XG conceived and designed the experiments. JM analyzed the data and wrote the paper. YC, WW, and BZ collected the data. WC contributed analysis tools, LZ, LX, HG, JL, and XG discussed and improved the article. All authors read and approved the final manuscript.

Identification of Endothelial Cell Genes Expressed in an In Vitro Model of Angiogenesis: Induction of ESM-1, βig-h3, and NrCAM

Transcription of the Non-coding RNA Upperhand Controls Hand2 Expression and Heart Development

Induction of Interactions between CD44 and Hyaluronic Acid by a Short Exposure of Human T Cells to Diverse Pro-inflammatory Mediators

Identifying Signatures of Natural Selection in Tibetan and Andean Populations Using Dense Genome Scan Data

A Comparison of Gene Expression and DNA Methylation Patterns across Tissues and Species

Delineation ofCCDC39/CCDC40mutation Spectrum and Associated Phenotypes in Primary Ciliary Dyskinesia

Acyl CoA Synthetase 5 (ACSL5) Ablation in Mice Increases Energy Expenditure and Insulin Sensitivity and Delays Fat Absorption

BLAST+: Architecture and Applications

Identification of Crucial Genes in Abdominal Aortic Aneurysm by WGCNA

UpSetR: an R Package for the Visualization of Intersecting Sets and Their Properties

Impaired Breakdown of Bradykinin and its Metabolites as a Possible Cause for Pulmonary Edema in COVID-19 Infection

Reduced Expression of microRNA-130a Promotes Endothelial Cell Senescence and Age-dependent Impairment of Neovascularization

Loss-of-function Mutations in the Human Ortholog of Chlamydomonas Reinhardtii ODA7 Disrupt Dynein Arm Assembly and Cause Primary Ciliary Dyskinesia

Flow Regulates Intercellular Communication in HAEC by Assembling Functional Cx40 and Cx37 gap Junctional Channels

OrthoFinder: Solving Fundamental Biases in Whole Genome Comparisons Dramatically Improves Orthogroup Inference Accuracy

C. elegans EGL-9 and Mammalian Homologs Define a Family of Dioxygenases that Regulate HIF by Prolyl Hydroxylation

Genomic Organization and Expression of KCNJ8/Kir6.1, a Gene Encoding a Subunit of an ATP-Sensitive Potassium Channel

Comparative microRNA Transcriptomes in Domestic Goats Reveal Acclimatization to High Altitude

Nuclear Gene Mutations as the Cause of Mitochondrial Complex III Deficiency

TGFBI Modulates Tumour Hypoxia and Promotes Breast Cancer Metastasis

Temporal Transcriptome Analysis Suggest Modulation of Multiple Pathways and Gene Network Involved in Cell-Cell Interaction during Early Phase of High Altitude Exposure

Metabolic Insight into Mechanisms of High-Altitude Adaptation in Tibetans

Expression Pattern of Heme Oxygenase 1 Gene and Hypoxic Adaptation in Chicken Embryos

Involvement of Vascular Endothelial Growth Factor Signaling in CLR/RAMP1 and CLR/RAMP2-mediated Pro-angiogenic Effect of Intermedin on Human Vascular Endothelial Cells

Fifteen Years of Genomewide Scans for Selection: Trends, Lessons and Unaddressed Genetic Sources of Complication

Anti-colony-stimulating Factor Therapies for Inflammatory and Autoimmune Diseases

Comparative Transcriptomics of 3 High-Altitude Passerine Birds and Their Low-Altitude Relatives

Heart Rate Recovery after Maximal Exercise Is Associated with Acetylcholine Receptor M2 (CHRM2) Gene Polymorphism

Genetic Variants at the EGLN1 Locus Associated with Highaltitude Adaptation in Tibetans Are Absent or Found at Low Frequency in highland Andeans

CD40 Ligand on Activated Platelets Triggers an Inflammatory Reaction of Endothelial Cells

Divergent Effects of IL-12 and IL-23 on the Production of IL-17 by Human T Cells

Mammalian Glutamyl Aminopeptidase Genes (ENPEP) and Proteins: Comparative Studies of a Major Contributor to Arterial Hypertension

Systematic and Integrative Analysis of Large Gene Lists Using DAVID Bioinformatics Resources

The Dual Role of AQP4 in Cytotoxic and Vasogenic Edema Following Spinal Cord Contusion and its Possible Association with Energy Metabolism via COX5A

Gene Coexpression Network Analysis Unraveling Transcriptional Regulation of High-Altitude Adaptation of Tibetan Pig

Maternal Low-Protein Diet Induces Gender-dependent Changes in Epigenetic Regulation of the Glucose-6-Phosphatase Gene in Newborn Piglet Liver

Thymus Cell Antigen 1 (Thy1, CD90) Is Expressed by Lymphatic Vessels and Mediates Cell Adhesion to Lymphatic Endothelium

Effect of Hypoxia on Integrin-Mediated Adhesion of Endothelial Progenitor Cells

x Frontiers in Genetics | www.frontiersin.org

Transcriptome Differences in Frontal Cortex between Wild Boar and Domesticated Pig

TNF-related Activation-Induced Cytokine (TRANCE) Induces Angiogenesis through the Activation of Src and Phospholipase C (PLC) in Human Endothelial Cells

Transcriptome Profile and Unique Genetic Evolution of Positively Selected Genes in Yak Lungs

WGCNA: an R Package for Weighted Correlation Network Analysis

NOD1 and NOD2 Regulate Proinflammatory and Prolabor Mediators in Human Fetal Membranes and Myometrium via Nuclear Factor-Kappa B1

Mitochondria in the Human Heart

Regulation of Bone Marrow Hematopoietic Stem Cell Is Involved in High-Altitude Erythrocytosis

featureCounts: an Efficient General Purpose Program for Assigning Sequence Reads to Genomic Features

Genomewide Association Studies for Female Fertility Traits in Chinese and Nordic Holsteins

Cloning of a Human Seven-Transmembrane Domain Receptor, LESTR, that Is Highly Expressed in Leukocytes

Glycoprotein Ibα Forms Disulfide Bonds with 2 Glycoprotein Ibβ Subunits in the Resting Platelet

Bioinformatics Analyses and Biological Function of lncRNA ZFPM2-AS1 and ZFPM2 G-ene in H-epatocellular C-arcinoma

Potent and Broad-Spectrum Antimicrobial Activity of CXCL14 Suggests an Immediate Role in Skin Infections

Effects ofHOX Homeobox Genes in Blood Cell Differentiation

Short-term Exposure to High Altitude Causes Coagulation Activation and Inhibits Fibrinolysis

The Endothelin Family: an Overview

Evolutionary Dynamics of Gene and Isoform Regulation in Mammalian Tissues

Effects of Grazing and Precipitation on Herbage Biomass, Herbage Nutritive Value, and Yak Performance in an Alpine Meadow on the Qinghai-Tibetan Plateau

Influence of High Altitude Exposure on the Immune System: a

Proteome Profiling in IL-1β and VEGF-Activated Human Umbilical Vein Endothelial Cells Delineates the Interlink between Inflammation and Angiogenesis

Mapping and Quantifying Mammalian Transcriptomes by RNA-Seq

Iron Metabolism in High-Altitude Residents

Muscle Transcriptome Signature and Gene Regulatory Network Analysis in Two Divergent Lines of a Hilly Bovine Species Mithun (Bos frontalis)

Evolutionary Dynamics of Coding and Non-coding Transcriptomes

Similarities and Differences between IL11 and IL11RA1 Knockout Mice for Lung Fibro-Inflammation

Stimulation of CD30 in Anaplastic Large Cell Lymphoma Leads to Production of Nuclear Factor-kappaB P52, Which Is Associated with Hyperphosphorylated Bcl-3

RNA Sequencing: Advances, Challenges and Opportunities

Cross-Species Insights into Genomic Adaptations to Hypoxia

CAPS Mutations Are Potentially Associated with Unexplained Recurrent Pregnancy Loss

Down-Regulation ofEPAS1Transcription and Genetic Adaptation of Tibetans to High-Altitude Hypoxia

Transcriptlevel Expression Analysis of RNA-Seq Experiments with HISAT, StringTie and Ballgown

A Functional Interleukin 12 Receptor Complex Is Composed of Two -type Cytokine Receptor Subunits

Regulation of Angiogenesis by Hypoxia: Role of the HIF System

The Transcriptomic Landscape of Yaks Reveals Molecular Pathways for High Altitude Adaptation

The Yak Genome and Adaptation to Life at High Altitude

edgeR: a Bioconductor Package for Differential Expression Analysis of Digital Gene Expression Data

RNA-seq Transcriptomics and Pathway Analyses Reveal Potential Regulatory Genes and Molecular Mechanisms in High-and Low-Residual Feed Intake in Nordic Dairy Cattle

Transferrin Receptor Is Negatively Modulated by the Hemochromatosis Protein HFE: Implications for Cellular Iron Homeostasis

Sequence Diversity at the Proximal 14q32.1 SERPIN Subcluster: Evidence for Natural Selection Favoring the Pseudogenization of SERPINA2

Inverse Relationship between Sec14l3 mRNA/protein Expression and Allergic Airway Inflammation

Frontiers in Genetics | www.frontiersin.org

Genome Wide Expression Analysis Suggests Perturbation of Vascular Homeostasis during High Altitude Pulmonary Edema

Sympathetic Stimulation Facilitates Thrombopoiesis by Promoting Megakaryocyte Adhesion, Migration, and Proplatelet Formation

Genetic Evidence for High-Altitude Adaptation in Tibet

Hypoxia and Mitochondrial Oxidative Metabolism

High-Altitude Adaptation: Mechanistic Insights from Integrated Genomics and Physiology

Update on High-Altitude Pulmonary Edema: Pathogenesis, Prevention, and Treatment

A Global View of Gene Activity and Alternative Splicing by Deep Sequencing of the Human Transcriptome

Comparative Transcriptomics of 5 High-Altitude Vertebrates and Their Low-Altitude Relatives

High Altitude Genetic Adaptation in Tibetans: No Role of Increased Hemoglobin-Oxygen Affinity

Effect of Clostridium perfringens β-Toxin on Platelets

Decoding the Long Noncoding RNA during Cardiac Maturation

Differential Regulation of the Two Genes Encoding Saccharomyces cerevisiae Cytochrome C Oxidase Subunit V by Heme and the HAP2 and REO1

Semidomesticated and Irreplaceable Genetic Resource Gayal (Bos frontalis) Needs Effective Genetic Conservation in Bangladesh: A Review

Developing Human Lung and Kidney Are Major Sites for Synthesis of Bone Morphogenetic Protein-3 (Osteogenin)

Detection and Integrated Analysis of lncRNA and mRNA Relevant to Plateau Adaptation of Yak

Different Gene Expressions between Cattle and Yak Provide Insights into High-Altitude Adaptation

KLF4, a Key Regulator of a Transitive Triplet, Acts on the TGF-β Signaling Pathway and Contributes to High-Altitude Adaptation of Tibetan Pigs

Endothelial Cells by Inactivation of VHL Gene Direct Angiogenesis, Not Vasculogenesis via Twist1 Accumulation Associated with Hemangioblastoma Neovascularization

O2 Sensing in Hypoxic Pulmonary Vasoconstriction: the Mitochondrial Door Re-opens

Comparison of Histological Structure of Pulmonary Alveoli between 180 Days Old Yak and plain Cattle

Genetic Characterisation of the Mithun (Bos frontalis) and Studies of Spermatogenesis, Blood Groups and Haemoglobins of its Hybrids with Bos indicus

Identification of a Tibetan-specific Mutation in the Hypoxic Gene EGLN1 and its Contribution to High-Altitude Adaptation

LINC00265 Targets miR-382-5p to Regulate SAT1, VAV3 and

Transcriptome Profiles Revealed the Mechanisms Underlying the Adaptation of Yak to High-Altitude Environments

Gene Co-expression Network Analysis Reveals Common System-Level Properties of Prognostic Genes across Cancer Types

Fatty Acid Metabolism Pathway Play an Important Role in Carcinogenesis of Human Colorectal Cancers by Microarray-Bioinformatics Analysis

Ubiquinol-cytochrome C Reductase Core Protein 1 Overexpression Protects H9c2 Cardiac Cells against Mimic Ischemia/ reperfusion Injury through PI3K/Akt/GSK-3β Pathway

Sequencing of 50 Human Exomes Reveals Adaptation to High Altitude

Comparative Transcriptomic and Proteomic Analyses Provide Insights into the Key Genes Involved in High-Altitude Adaptation in the Tibetan Pig

Characterization and Comparative Analyses of Muscle Transcriptomes in Dorper and Small-Tailed Han Sheep Using RNA-Seq Technique

Genetic Adaptation of Schizothoracine Fish to the Phased Uplifting of the Qinghai-Tibetan Plateau

Hypoxia Adaptations in the Grey Wolf (Canis lupus Chanco) from Qinghai-Tibet Plateau

Mammalian Homeobox B6 Expression Can Be Correlated with Erythropoietin Production Sites and Erythropoiesis during Development, but Not with Hematopoietic or Nonhematopoietic Stem Cell Populations

Conflict of Interest: The authors declare that the research was conducted in the

This study was supported by the National Natural Science Foundation of China (31572376), the Agricultural Science and Technology Innovation Program in Chinese Academy of Agricultural Sciences (ASTIP-IAS03), and the Science and Technology Project of Inner Mongolia Autonomous Region (2020GG0210), and the China Agriculture Research System of MOF and MARA (CARS-37). The funders had no role in study design, data collection, and analysis, decision to publish, or preparation of the manuscript.

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fgene.2021.778788/ full#supplementary-material