DNA methylation, microRNA expression profiles and their relationships with transcriptome in grass-fed and grain-fed Angus cattle rumen tissue
Authors:
Yaokun Li aff001; José A. Carrillo aff002; Yi Ding aff002; Yanghua He aff002; Chunping Zhao aff003; Jianan Liu aff002; Linsen Zan aff003; Jiuzhou Song aff002
Authors place of work:
College of Animal Science, South China Agricultural University, Guangzhou, Guangdong, P.R. China
aff001; Department of Animal & Avian Sciences, University of Maryland, College Park, Maryland, United States of America
aff002; College of Animal Science and Technology, Northwest A&F University, Yangling, Shaanxi, P.R. China
aff003
Published in the journal:
PLoS ONE 14(10)
Category:
Research Article
doi:
https://doi.org/10.1371/journal.pone.0214559
Summary
Rumen is an organ for supplying nutrients for the growth and production of bovine, which might function differently under grass-fed and grain-fed regimens considering the association of gene expression, DNA methylation, and microRNA expression. The objective of this study was to explore the potential mechanism influencing rumen function of grass-fed and grain-fed animals. Methylated DNA binding domain sequencing (MBD-Seq) and microRNA-Seq were respectively utilized to detect the DNA methylation and microRNA expression in rumen tissue of grass-fed and grain-fed Angus cattle. Combined analysis revealed that the expression of the differentially expressed genes ADAMTS3 and ENPP3 was correlated with the methylation abundance of the corresponding differentially methylated regions (DMRs) inside these two genes, and these two genes were reported to be respectively involved in biosynthesis and regulation of glycosyltransferase activity; the differentially expressed microRNA bta-mir-122 was predicted to possibly target the differentially expressed genes OCLN and RBM47, potentially affecting the rumen function; the microRNA bta-mir-655 was exclusively detected in grain-fed group; its targets were significantly enriched in insulin and TGF-beta signaling pathways, which might worked together to regulate the function of rumen, resulting in different characteristics between grass-fed and grain-fed cattle. Collectively, our results provided insights into understanding the mechanisms determining rumen function and unraveled the biological basis underlying the economic traits to improve the productivity of animals.
Keywords:
DNA methylation – Gene expression – Gene regulation – Polymerase chain reaction – MicroRNAs – Cattle – Gene ontologies – Beef
Introduction
Historically, most of beef products were from grass-finished cattle. Since the 1950’s, numerous studies were performed to improve the efficiency of beef production; meanwhile, the beef cattle feedlot industry began to emerge, where high-energy grains were utilized to improve the productivity of cattle. And, consumers have been accustomed to grain-fed beef, considering the flavor and palatability resulted from the diets with large proportion of high-energy grain [1]. However, due to new research on the effects of diverse feeding regimens, consumers’ preferences for beef quality have changed, making certain producers to regress to the pastoral beef production in spite of the feeding inefficiency. Recently, growing number of consumers are interested in products obtained from grass-finished animals, raising concerns about the quality difference between grass-fed and grain-fed beef [2]. During the past few decades, it was elucidated that altered feeding regimes could bring about changes in the nutritional quality of beef. The essential fatty acids omega-3 and omega-6 are critical for animal health; however, they cannot be produced by the host and must be obtained from food [2]. Studies illustrated that grass-fed beef had a significantly higher content of omega-3 and displayed a more favorable omega-3/omega-6 ratio compared with grain-fed beef [3–4]. Beta-carotene could modulate immune reaction and prevent animals from bacterial and viral infection. Compared with grain-fed animals, beta-carotene was more abundant in the muscle of grass-fed animals [5–8]. Diterpenoids and derivatives of chlorophyll were also higher in grass-fed beef than in grain-fed beef, which played important roles in changing the flavor and aroma of the cooked beef [9]. Additionally, higher level of vitamin E was found in the grass-fed beef than the beef products from concentrate diets [10]. Besides nutritional components, the transcriptomic difference of cattle could also be detected due to varied feeding regimens [11–12].
In the field of genetics, epigenetics refers to cellular and physiological phenotypic trait alterations because of external factors that can influence genes expression without changing the DNA sequence [13]. As one important type of the epigenetic mechanism, DNA methylation was generally used to regulate gene expression; disturbing it could cause abnormal expression of certain genes [14]. DNA methylation has been largely studied due to its involvement in regulating most of the biological processes, including embryonic development, chromosome stability, chromatin structure, and transcription [15–19]. DNA methylation can interfere with transcription factor binding, indirectly repressing the gene activity through the recruitment of methyl-CpG binding domain (MBD) protein that can alter chromatin structure [20]. The repression strength could be determined by the concentration of DNA methylation. In porcine, differentially methylated regions in the promoter have been found to be associated with the repression of genes related to obesity and some novel genes [21]. Up to now, genome-wide DNA methylation profiles in many animal breeds, including chicken, pig, arabidopsis, human and bovine, have been reported [21–26]. However, the DNA methylation pattern of bovine rumen tissue was rarely explored.
Additionally, the gene expression could also be regulated by miRNAs, which could specifically bind to corresponding mRNAs, leading to the mRNA silencing and translation repression. It is known that each miRNA might have many mRNA targets; and each mRNA could also be regulated by more than one miRNA [27]. In mammals, miRNAs are predicted to regulate the expression of ~60% of all protein-coding genes, and found to be associated with many biological processes, including cell apoptosis, cell growth and differentiation, and embryo development [28–29].
As the largest part of the stomach, rumen constitutes the important workplace for feed digestion and microorganism fermentation, which is critical for supplying nutrients for animals’ growth, development and production. Given the emerging roles of miRNAs and DNA methylation in gene regulation, identifying the DNA methylation profile and expression pattern of miRNAs is critical to understand their functions in the process of bovine development. Therefore, we utilized Illumina sequencing technology to explore the miRNA expression pattern and DNA methylation profile in the rumen tissue of Angus cattle raised under grass-fed and grain-fed regimens. Further, we conducted the combined analysis of miRNAs, DNA methylation, and mRNA expression, which was important for further functional study and valuable biomarkers detection. Accordingly, our results would be significant for exploring potential mechanisms that led to differences observed between animals under altered feeding regimens.
Materials and methods
Ethics statement
We conducted the animal experiments according to the NIH guidelines for housing and care of laboratory animals. The Institutional Animal Care and Use Committee (IACUC) of University of Maryland at College Park reviewed and authorized the protocols (permit number R-08-62).
Sample collection
The studied Angus cattle were born and raised at the Wye Angus farm. This farm has produced many progenies with similar genetic background. All experimental animals in this study received the same diet until weaning. After weaning, the calves were randomly assigned to one diet (grass-fed or grain-fed) and followed that regimen until slaughter. The grain-fed steers received conventional diet comprised of trace minerals, shelled corn, corn silage and soybean. The grass-fed steers normally ate grazed alfalfa; in winter, we provided bailage for the animals. The grass-fed animals consumed no industrial or agricultural byproducts and did not receive any type of grain. The grain-fed individuals met the market weight around the age of 14 months; grass-fed steers reached the similar weight value with the age of 20 months. After slaughter at the Old Line Custom Meat Company (Baltimore, MD), the rumen samples from all experimental animals were immediately excised at the same location around the cardiac ostium, and then the samples were rinsed and kept at -80°C for subsequent analysis. For MBD-Seq, two animals were contained per experimental group. Based on the MBD-Seq experiment, we respectively added one animal per group in miRNA-Seq; thus, three animals were respectively utilized in grass-fed and grain-fed group for miRNA-Seq analysis.
DNA extraction and MBD-Seq library construction
Firstly, we extracted the genomic DNA from the four experimental animals through the Wizard Genomic DNA purification kit (Promega, A1120), and measured the concentration via the Qubit dsDNA Broad-Range Assay (Invitrogen, Q32850). The genomic DNA from each individual was adjusted to 0.1μg/μl for a final volume of 55 μl and sheared into 300–500 bp fragments. Secondly, MethylCap kit (Diagenode, C02020010) was employed to obtain the methylated DNA (meDNA). The 141.8 μl of capture reaction mix containing 12 μl of sheared DNA was prepared according to the instruction of the kit. Then, the capture reaction mix (119 μl) was incubated with 1 μl of diluted MethylCap protein at 40 rpm for 2 hours at 4°C. The rest 22.8 μl was utilized as input sample. Next, meDNA Capture beads (coated with GSH) provided by the kit were utilized to capture methylated DNA. The bead pellet combined by the methylated DNA was collected. After washing the bead pellet, we performed the elution of the captured DNA; 150 μl of low, medium, and high concentration elution buffer were sequentially utilized per reaction. All fractions were treated by the MiniElute PCR Purification Kit (QIAGEN, 28006).
We repaired the end of the fragmented methylated DNA by using the NEBNext End Repair Module (NEB, E6050S, USA). Then the 3’ poly “A” was added through DNA Polymerase I, Large (Klenow) Fragment (NEB, M0210L, USA). Next, we ligated a pair of Solexa adaptors (Illumina) to the repaired ends by T4 ligase (Promega, M1801, USA). The fragments of the ligated products (DNA plus adaptors) from 200 to 500bp were selected from 2% agarose gels and purified by QIA- quick Gel Extraction Kit (QIAGEN, USA). Then, the Phusion Hot Start High-Fidelity DNA Polymerase (NEB, M0530S, USA) was used to enrich the purified DNA templates through PCR. After purification of the PCR products (MinElute PCR Purification Kit, QIAGEN, USA), the concentration of the DNA library was measured through the Qubit assay (Life Technology, Q32850). Finally, the DNA sequencing was performed in the Solexa 1G Genome Analyzer (Illumina) according to the specification provided by the manufacturer.
MBD-Seq data analysis
We evaluated the quality of the raw reads through FastQC, which is a web-based software. It could thoroughly examine the reads and create a detailed quality assurance report containing “per base sequence quality” and “per sequence GC content”, etc. Then, Bowtie was used to align the sequencing reads to the reference genome (Bos_taurus_UMD3.1/bosTau6) [30], which was obtained from the iGenomes web site (http://support.illumina.com/sequencing/sequencing_software/igenome.html). During this step, according to the indicators generated by FastQC, we trimmed the first 10 bases and the last 5 bases of each read (50 bp) to keep high sequence quality score; and 35 bp tags were then obtained. For data format conversion, SAMtools and BEDtools were applied in our analysis.
For peaks identification, the Model Based Analysis of ChIP-Seq (MACS) was employed [31]. Identification of the DMRs was performed by the R package Diffbind; it calculates the differentially bound sites using affinity data. The input data for DiffBind was the bam file containing the aligned reads and the peaks set identified by MACS. The default method TMM (Trimmed Mean of M-values) considering the effective library size was used for normalization. After creating a contrast between conditions, DiffBind carried out a DESeq2 analysis with a false discovery rate (FDR) < 0.1 to call the DMRs. DESeq2 approach employed shrinkage estimation for fold changes and dispersions, making the estimated values more stable and interpretable [32].
The ChIPpeakAnno package was utilized for the genomic annotation of the previously identified DMRs [33]. This software provides DMRs information about the distance, relative position and overlaps for the inquired features, including the information relative to the exon and the transcriptional start site (TSS) of the corresponding genes. The annotation information was obtained from BioMart, using Ensembl 80 in the archive site; dataset “Bos taurus_genes_ensembl (UMD3.1)” corresponding to the bosTau6 reference genome was used for alignment. The CpG islands annotation was retrieved from the UCSC web browser. Finally, the DMRs were annotated based on specific genomic features. In addition, we combined the annotation results of DMRs and differentially expressed genes (DEGs) which has been published in 2015 [34], and performed the association analysis between the methylation abundance of the DMRs and the expression level of their corresponding genes including the DEGs between grass-fed and grain-fed Angus cattle through software SPSS Statistics 17.0, from which we hypothesized the relationship between the methylation abundance of the DMRs and the expression level of their corresponding genes.
Bisulfite sequencing for MBD-Seq Validation
After quality evaluation and quantification, equal amounts of genomic DNA from two experimental animals of each group were pooled together, which were then served as the template for the bisulfite conversion. Firstly, 500 ng of each DNA pool was treated by the sodium bisulfite conversion reagents (Methyl EdgeTM Bisulfite Conversion System, Promega, USA). We randomly selected the DMRs from the bioinformatics analysis results for validation. The PCR primers were designed via MethPrimer (http://www.urogene.org/methprimer/) and shown in S1 Table. Then, PCR products were purified using QIA- quick Gel Extraction Kit (QIAGEN, USA), which were subsequently ligated to pGEM-T Vector (pGEM-T Vector System I, Promega, USA), and transformed to DH5α competent cells (Z-Competent E. Coli Cells—Strain Zymo 5α, ZYMO Research, USA) for screening successful insertions (blue-white selection). Next, ten white colonies from each culture plate were cultured overnight at 37°C in a shaker. The plasmid DNA was isolated utilizing Zyppy Plasmid Miniprep Kit (ZYMO Research, USA). BigDye Terminator v3.1 Cycle Sequencing Kit (Applied Biosystems, 4337456) was used for the sequencing in the ABI 3730 machine. Bisulfite sequencing results were analyzed through QUMA (http://quma.cdb.riken.jp).
RNA extraction and microRNA-Seq library construction
We extracted the total RNA from the six experimental individuals (three animals per group) by Trizol (Invitrogen, Carlsbad, CA, USA), the RNA was then treated by DNase digestion and purified by Qiagen RNeasy column (Invitrogen, Carlsbad, CA, USA), as previously described [35]. We checked the integrity and quality by a NanoDrop 1000 spectrophotometer and by visualization on a 1.5% agarose gel. The microRNA-Seq library was built utilizing NEBNext Multiplex Small RNA Library Prep Set for Illumina (Set 1) following the instructions provided by the manufacturer (NEB, E7300S/L, USA). The library construction was started with 800ng total RNA. Firstly, the 3’ adaptor and 5’ adaptor were sequentially ligated to the RNA. Then, we performed the reverse transcription to synthesize the first strand. Next, we executed the PCR amplification and added the 6-bp index to the DNA products; different libraries were assigned different indexes. Subsequently, we selected and recycled the DNA fragments from 140 to 150 bp using 6% polyacrylamide gel (6% Novex® TBE PAGE gel, Life Technology, USA), which were corresponding to RNAs from 21 to 30 bp. After purification, the concentration was measured by the Qubit assay (Life Technology, Q32850). Finally, we sequenced the libraries identified by the 6-bp index through an Illumina HiSeq 2000 sequencer, as described previously [36].
microRNA-Seq data analysis
After quality evaluation of the tags performed by FastQC, all samples were analyzed employing miRDeep* software, which could quantify known and novel miRNAs [37]. This program had a user-friendly graphic interface implemented totally in java and accepted raw data in FastQ and SAM/BAM format as input. The default length of miRNAs was set as 18–23 nucleotides. Low-quality reads were filtered out at the alignment stage; the read with a Phred quality score of 20 or higher was considered as good read. Meanwhile, multi-mapping reads with alignments to more than 100 genomic loci were also filtered out; the score based on the probabilistic score of the potential miRNA precursor was set as -10 as default [37]. In the further miRNA expression analysis, we only considered known miRNAs. Expression values for known miRNAs were estimated individually by miRDeep* in each sample. Subsequently, these expression levels were extracted and recorded in a matrix, which was used later as input for the edgeR package to call the differentially expressed miRNAs.
After identification of the differentially expressed miRNAs, the target genes of those miRNAs were predicted by TargetScan [38], which were then utilized for the combined analysis with transcriptome analysis results (including the DEGs) published in 2015 [34]. We screened the overlapped genes of the differentially expressed miRNAs targets and the DEGs; based on the discovered overlapped DEGs, we then searched the corresponding miRNAs and performed the association analysis between the expression of the differentially expressed miRNAs and the corresponding DEGs through software SPSS Statistics 17.0, from which the relationship between the differentially expressed miRNAs and the DEGs could be analyzed.
Quantitative real-time PCR (qRT-PCR) analysis
qRT-PCR was performed to verify the differentially expressed miRNAs found in the microRNA-Seq analysis on the iCycler iQ PCR system (Bio-Rad, Hercules, CA, USA). The qRT-PCR reactions were performed with a QuantiTect SYBR Green PCR Kit (Qiagen, Valencia, CA, USA) according to the manufacturer’s instructions. We conducted three technical replicates and two independent biological replicates for each product. RPS18 was chosen as the control gene [39]. The primer sequences were listed in S2 Table.
Results
Landscape of the DNA methylomes
For the methylome study, we had four experimental samples in total; the alignment levels were 94.92%, 82.98%, 95.54% and 96.10%, respectively (S1 Fig). We finally found 217 DMRs between the two groups (S3 Table), which were presented as red dots in Fig 1. Based on the identified DMRs, the cluster analysis of our experimental samples was performed (Fig 2), the results showed consistency with the group assignment of the samples, suggesting that these DNA fragments carrying DMRs could be used to predict the potential mechanism causing the difference between our two groups. From those DMRs, we found that only two DMRs were highly methylated and all the others low methylated in the grain-fed bovines compared with the grass-fed group. For the distribution of the DMR length, the average was 3760 bp with extreme value of 485 bp and 14,387 bp. Approximately one percent of the DMRs were less than 1,000 bp and 37.8% of the DMRs accounted for fragments longer than 4,000 bp (Fig 3). Fig 4 showed the number of DMRs per chromosome, chromosome 8 accounted for 18 DMRs indicating the largest number, following by chromosome 1 with 15 DMRs. Chromosomes 19, 21, 23, and 15 had the lowest number of DMRs with values of 2, 2, 1 and 1, respectively.
The binding affinity between grass-fed and grain-fed cattle was shown in Fig 5. We found that global DNA methylation level decreased in grain-fed individuals. Then, we detected the DNA methylation level in the 2 kb region upstream of the transcription start site (TSS), the 2 kb region downstream of the transcription end site (TES) and the gene body region from TSS to TES. In our experimental groups, the DNA methylation level declined significantly before the TSS and increased notably towards the gene body region with slight changes before the TES, followed by a sharp decrease in the downstream of TES (Fig 6). Compared with grass-fed individuals, the grain-fed steers showed a higher level of DNA methylation around the gene body region.
The DMR annotation was subsequently performed (S4 Table). The distribution of DMRs’ distance to the nearest TSS could be found in Fig 7, suggesting that approximately 14.3% of the DMRs are located within a range of 10,000 bp from the TSS. Moreover, DMRs’ location related to genes and CpG islands were summarized in Fig 8. We found that almost a quarter of the DMRs were contained within genes and half in regions upstream of the TSS. The DMRs downstream of the genes accounted for 25%. Only four genes were found in the DMRs. Regarding CpG islands, approximately 8.8% of the annotated CpG islands were included in the DMRs. Most of the DMRs were outside of the CpG island boundaries with 42.4% upstream and 50% downstream of the CpG islands.
Combined analysis of MBD-Seq with transcriptome (RNA-Seq)
In this study, we screened the genes with promoters containing DMRs, using 10,000 bp as the parameter for maximum distance. From the 217 identified DMRs, 21 were located within the promoters, corresponding to 21 different genes. Our previous study explored the transcriptome profiling in grass-fed and grain-fed animals, and 342 DEGs (FDR < 0.1) were detected between the two groups [34]. We found that, of the above 21 genes, gene SEC22C (Ensemble gene ID: ENSBTAG00000006570) was positively correlated with the corresponding DMR (r > 0.9, p = 0.055), gene HGSNAT (Ensemble gene ID: ENSBTAG00000033182) was negatively correlated with the corresponding DMR (r < -0.9, p = -0.029), however, none of the them could be found in the DEGs list; the p-value of the other genes were all larger than 0.1 (S5 Table). Additionally, 57 DMRs were found inside 52 genes, and all of the DMRs were highly methylated in grass-fed steers. And two of the corresponding genes were detected in the DEGs list which were ADAMTS3 and ENPP3, respectively; ADAMTS3 (Ensemble gene ID: ENSBTAG00000006507) was negatively correlated with the corresponding DMRs (r < -0.6), ENPP3 (Ensemble gene ID: ENSBTAG00000020196) was also negatively correlated with the DMR (r = -0.143), however, the p-value of these two genes was all larger that 0.05; generally, we found that 55.8% of the corresponding 52 genes were negatively correlated with the DMR, only one showed significant correlation (r = -0.967, p < 0.05) (S6 Table). Meanwhile, we found that 4 genes were inside the DMRs; however, none of them were differentially expressed. Based on the annotation regarding the TSS, another two genes could also be discovered in the DEGs list, which were CRISPLD1 and PRR5; CRISPLD1 was negatively correlated with the corresponding DMR and PRR5 was positively correlated with the DMR; however, the difference was not significant (p > 0.05); the nearest distance of the corresponding DMR to the TSS was respectively 14,973 and 45,887 bp.
MBD-Seq and microRNA-Seq data validation
In order to assess the reliability and accuracy of DMRs detection from MBD-Seq, 5 regions were randomly selected for validation. The calculated methylation levels are shown in the Fig 9. We found that the methylation level of all the 5 regions agreed with MBD-Seq analysis. The bisulfite sequencing employed only a segment of the methylation region to perform the validation, thus the magnitude of methylation difference of the validation results could not be exactly the same as MBD-Seq analysis.
For the differentially expressed miRNA bta-mir-122, the expression level in grain-fed ruminal wall was significantly higher than in grass-fed ruminal wall, and the results of qPCR and RNA-Seq suggested the same direction (S2 Fig).
Expression profiles of microRNAs in grass-fed and grain-fed groups
To identify the microRNAs involved in bovine rumen function, total RNAs from three biological replicates for each condition were used to construct small RNA libraries. Via high-throughput sequencing, we obtained an average of 15 million reads per sample. The percentage of alignment for all samples exceeded 90% as shown in Fig 10. Based on the aligned sequencing reads, we totally identified 321 known miRNAs (S7 Table), of which 72.9% were high expressed in grain-fed individuals compared with grass-feed ones. Between the two groups, we identified only one differentially expressed miRNA (FDR < 0.1), which was bta-mir-122 with improved expression in grain-fed group as compared with grass-fed animas. Additionally, we found that microRNA bta-mir-655 was exclusively expressed in the grain-fed group.
Target gene prediction, Gene Ontology and pathways enriched by the differentially expressed microRNAs
Prediction of target genes of the miRNA was performed using TargetScan (http://www.targetscan.org/). For bta-mir-122 and bta-mir-655, we respectively found 145 and 749 target genes. To explore the specific functional features shared by the targets of the above miRNA, online software David Bioinformatics Resources 6.7 was used to perform the GO enrichment analysis, mainly regarding the biological processes, cellular components and molecular functions. For the targets of bta-mir-122, the most significant GO terms were: regulation of GTPase activity, regulation of protein signal transduction, intracellular, GTPase activator activity and GTPase regulator activity (Table 1). To some extent, the targets were also involved in the pathway glycolysis/gluconeogenesis (P = 7.9×10−2), which could produce energy and is an critical metabolic pathway employed by a host of organisms. Therefore, this might be an interesting pathway to study in rumen. For the targets of bta-mir-655, the most significant GO terms of the three categories were: regulation of cellular process, regulation of biological process, biological regulation (Table 2). The complete list of GO terms enriched with the target genes of bta-mir-655 was provided in S8 Table. Then, we detected the pathways that involved the targets of bta-mir-655, the most significantly enriched pathways were shown in Table 3. Several interesting pathways were found including adherens junction, insulin signaling pathway, TGF-beta signaling pathway and neurotrophin signaling. These results would provide prior knowledge to explain the difference between grass-fed and grain-feed Angus cattle.
Combined analysis of microRNA-Seq with transcriptome (RNA-Seq)
For bta-mir-122, two of the 145 target genes were found in the DEGs list from RNA-Seq, which were OCLN and RBM47, respectively. These two genes were all highly expressed in grass-fed group, showing a negative correlation with bta-mir-122, while the difference was not significant for the correlation (p > 0.05). For bta-mir-655 exclusively expressed in grain-fed animals, we found that 14 of the 749 targets were overlapped with the DEGs. Among the 14 overlapped genes, four genes displayed increased expression in grain-fed group, among which GAS1, ANO1 and NFASC displayed significantly positive correlation with bta-mir-655 (p < 0.05); the other nine genes (MYO10, SASH1, EMP1, SLC14A1, PCDH19, IRX5, MAL2, FAM84A, AHCYL2 and DSG3) exhibited enhanced expression level in grass-fed group, which demonstrated negative correlation with bta-mir-655, but none of them showed significant difference (p > 0.05) (S9 Table).
Discussion
In addition to genetic makeup, age, species, gender and environment, grass and grain ratio of the diet could also contribute to significant difference in the general profile of fatty acid, animals’ growth rate, the content of vitamin A and E in the meat, the antioxidant enzyme concentration in body tissues and lipid depots [10, 40–41]. In the past several years, miRNA and DNA methylation have been studied widely. However, limited studies have focused on the rumen tissue of bovine, one of the most significant livestock in the world. Rumen was the most important workshop for the digestion of nutritional substance of bovine. The potential changes of rumen metabolism may have effects on the quality and quantity of protein, affecting other digestive organs, such as reticulum, small intestine and large intestine. Therefore, rumen function was critical for animals’ growth, health and productivity. In the present study, we mainly focused on the rumen under different feeding diets, detecting the genome-wide DNA methylation profiles and miRNA expression profiles, and we also performed combined analysis of these with transcriptomic results. The objective was to identify the DMRs and miRNAs that might potentially affect the rumen function. This study performed a comprehensive analysis of DNA methylation profiles in the rumen tissues of grass-fed and grain-fed animals, and totally revealed 217 DMRs across the whole genome. In our analysis, we found that 21 DMRs were located in the promoter regions; however, none of their corresponding genes were overlapped with the DEGs. For these DMRs, although the methylation level was different, the concentration might not be high enough to cause the expression difference of their targets. Or, these DMRs could actually exert no effects on the genes. Additionally, we detected 57 DMRs inside 52 genes, of which two genes (ADAMTS3 and ENPP3) were overlapped with the DEGs. The methylation abundance of the DMRs exhibited negative and positive correlation with the expression of the corresponding gene ADAMTS3 and ENPP3, respectively. Studies suggested that DNA methylation in the gene body region might change the chromatin structure and alter the transcription elongation efficiency [18]. Gene body methylation was more prevalent than promoter; its relationship with gene expression levels was very complex and not monotonic. Normally, gene body DNA methylation showed positive correlation with gene expression [42]. However, the information about the role of DNA methylation in the gene body region was still insufficient. For the two DEGs predicted to be correlated with the DMR in the gene body region, ADAMTS3 belonged to the ADAMTS (a disintegrin and metalloproteinase domain with thrombospondin motifs) metalloproteinase family mediating cartilage aggrecan degradation as well as collagen biosynthesis [43]. It was involved in the biosynthesis of type II procollagen, the main collagen of articular cartilage [44]. ENPP3 was known to be a typical ectoenzyme localized to the cell surface, playing a role in metabolizing extracellular nucleotides and their derivatives [45]. And, it was suggested that ENPP3 could modulate nucleotide-mediated signal transduction, which was known as purinergic signaling [46–47]. Recent study demonstrated that ENPP3 could catalyze the hydrolysis of the nucleotide sugar, and the levels of several intracellular sugars, including UDP-Fuc, UDP-GalNAc and UDP-GlcA, were significantly influenced by knocking down endogenous ENPP3 [48]. There was evidence showing that intracellular sugars played a role in regulating glycosyltransferase activity and controlling the total cellular glycosylation profile [49–51]. Therefore, through DNA methylation regulation, the expression of ENPP3 would be changed, which might indirectly regulate the activity of a broad range of glycosyltransferase and consequently influence the total cellular glycosylation pattern.
The miRNAs could negatively regulate gene expression through degrading the target mRNAs or inhibiting the translation [52]. However, the understanding of the contribution of miRNA to rumen function was still limited. In this study, we aimed to get a deeper insight into underlying mechanism of miRNA in rumen function. Between grass-fed and grain-fed group, we found one differentially expressed miRNA (bta-mir-122) with higher expression level in grain-fed steers. Bta-mir-122 was a conserved miRNA between vertebrate species. It was also found differentially expressed between normal and aberrant placental samples, suggesting that bta-mir-122 might be involved in the development of placentae [53]. In liver, bta-mir-122 was relevant in maintenance of homeostasis and had critical metabolic and anti-inflammatory functions [54]. In our study, we detected 145 target genes for bta-mir-122 in total. To explore the specific functional features shared by the 145 targets, we performed GO enrichment analysis. Results showed that the targets were mainly involved in the GTPase activity and regulation of protein signal transduction. GTP hydrolysis played an essential role in controlling numerous biological processes, including protein biosynthesis, growth control and differentiation, and various transport processes. Additionally, two of the 145 target genes were found in the DEGs list, which were OCLN and RBM47, respectively. OCLN was a key tight junction protein that could interact with intracellular signaling pathways, which played roles in regulating intestinal function [55]. RBM47 was a novel RNA-binding protein; it could contribute to the basic machinery for C to U RNA editing in intestine, influencing the expression of some genes [56]; however, the function of this gene has not been studied in rumen. We made the hypothesis that RBM47 would alter the expression of some genes in rumen, consequently regulating rumen function. Meanwhile, we detected one exclusively expressed miRNA (bta-mir-655) in grain-fed group; it could target 749 different genes. After the GO analysis based on the targets, the most significant GO terms included regulation of cellular process, biological regulation, developmental process, regulation of metabolic process, and transcription factor activity, which played a critical role in rumen function. The most interesting pathways included adherens junction, insulin signaling pathway and TGF-beta signaling pathway, which were related to animals’ growth, survival and development. In addition, 14 of the 749 targets were overlapped with the DEGs; however, none of them was discovered in the pathways. Thus, we hypothesized that bta-mir-655 just slightly changed the expression of the targets without significant expression change; but the target genes working together in the corresponding pathway could exert significant effects on the rumen function. Among the 14 overlapped genes, FAM84A, localized in the subcellular membrane region, was involved in invasion and/or metastasis of colon cancer cells influencing colorectal cancer [57]; SASH1 was a member of the SH3-domain containing expressed in lymphocytes (SLY1) gene family, it encodes signal adapter proteins which were composed of certain protein–protein interaction domains, showing prognostic significance in human cancer [58]. A previous study suggested that epithelial membrane protein 1 (EMP1) gene could prevent tumor proliferation and was associated with gastric carcinoma [59]. DSG3 was predicted to be involved in the GO term cell adhesion; it was reported to be associated with oncogenesis [60]. However, information about the function of these genes in rumen was still limited; they might play a potential role in animal development through affecting the gastrointestinal function. Accordingly, it might be of great interest to perform functional experiment of these genes to better understand the mechanisms causing the varied performance.
Our results provided evidence for explaining the molecular mechanisms leading to the differences between grass-fed and grain-fed cattle. For the function of the DEGs containing the DMR and the differentially expressed miRNAs, extensive experimental validation work was still needed. Thus, overexpression and inhibition of our identified DEGs and miRNAs could be considered for the functional validation, which would provide more supportive information for our findings.
Conclusions
In this study, we combined DNA methylation and miRNA expression with transcriptome analysis to explore the potential mechanism influencing rumen function of grass-fed and grain-fed animals. We found that the expression of ADAMTS3 and ENPP3 might be altered by the corresponding DMR inside these two genes, through which rumen function of grass-fed and grain-fed cattle could be regulated. For the differentially expressed miRNA bta-mir-122, it might modulate the rumen function by targeting the two DEGs OCLN and RBM47, which were possibly associated with gastrointestinal function. While expanding the scope of future studies with putative genes relevant to bovine growth and meat quality traits, our analysis provided biological insights into the mechanisms regulating rumen function and uncovered the molecular basis underlying the economic traits to enhance the productivity of animals.
Supporting information
S1 Fig [tif]
Alignment level of MBD-Seq reads to the bovine genome.
S2 Fig [tif]
Validation of differentially expressed miRNA.
S1 Table [xls]
Primers used for validating the randomly selected DMRs.
S2 Table [xls]
Primers used for validating the differentially expressed miRNA.
S3 Table [xls]
The DMRs between grass-fed and grain-fed Angus cattle.
S4 Table [xls]
Annotation of the identified DMRs.
S5 Table [xlsx]
The association analysis of the identified DMRs located within the genes promoters and their corresponding genes.
S6 Table [xlsx]
The association analysis of the identified DMRs inside the genes and their corresponding genes.
S7 Table [xlsx]
The known miRNAs in grass-fed and grain-fed animals.
S8 Table [xls]
The GO terms enriched with the target genes of bta-mir-655.
S9 Table [xlsx]
The association analysis between bta-mir-122 and bta-mir-655 and their targeted DEGs.
Zdroje
1. Wood JD, Richardson RI, Nute GR, Fisher AV, Campo MM, Kasapidou E, et al. Effects of fatty acids on meat quality: a review. Meat Sci 2004,66:21–32. doi: 10.1016/S0309-1740(03)00022-6 22063928
2. Daley CA, Abbott A, Doyle PS, Nader GA, Larson S. A review of fatty acid profiles and antioxidant content in grass-fed and grain-fed beef. Nutr J 2010,9:10. doi: 10.1186/1475-2891-9-10 20219103
3. French P, Stanton C, Lawless F, O’Riordan EG, Monahan FJ, Caffrey PJ, et al. Fatty acid composition, including conjugated linoleic acid, of intramuscular fat from steers offered grazed grass, grass silage, or concentrate-based diets. J Anim Sci 2000,78:2849–2855. doi: 10.2527/2000.78112849x 11063308
4. Wood JD, Enser M. Factors influencing fatty acids in meat and the role of antioxidants in improving meat quality. Br J Nutr 1997,78 Suppl 1:S49–60.
5. Descalzo AM, Insani EM, Biolatto A, Sancho AM, Garcia PT, Pensel NA, et al. Influence of pasture or grain-based diets supplemented with vitamin E on antioxidant/oxidative balance of Argentine beef. Meat Sci 2005,70:35–44. doi: 10.1016/j.meatsci.2004.11.018 22063278
6. Beauchesne-Rondeau E, Gascon A, Bergeron J, Jacques H. Plasma lipids and lipoproteins in hypercholesterolemic men fed a lipid-lowering diet containing lean beef, lean fish, or poultry. Am J Clin Nutr 2003,77:587–593. doi: 10.1093/ajcn/77.3.587 12600847
7. Smith DR, Wood R, Tseng S, Smith SB. Increased beef consumption increases apolipoprotein A-I but not serum cholesterol of mildly hypercholesterolemic men with different levels of habitual beef intake. Exp Biol Med (Maywood) 2002,227:266–275.
8. Hunninghake DB, Maki KC, Kwiterovich PO, Davidson MH, Dicklin MR, Kafonek SD. Incorporation of lean red meat into a National Cholesterol Education Program Step I diet: A long-term, randomized clinical trial in free-living persons with hypercholesterolemia. Journal of the American College of Nutrition 2000,19:351–360. doi: 10.1080/07315724.2000.10718931 10872897
9. Elmore JS, Warren HE, Mottram DS, Scollan ND, Enser M, Richardson RI, et al. A comparison of the aroma volatiles and fatty acid compositions of grilled beef muscle from Aberdeen Angus and Holstein-Friesian steers fed diets based on silage or concentrates. Meat Sci 2004,68:27–33. doi: 10.1016/j.meatsci.2004.01.010 22062004
10. De la Fuente J, Diaz MT, Alvarez I, Oliver MA, Font IFM, Sanudo C, et al. Fatty acid and vitamin E composition of intramuscular fat in cattle reared in different production systems. Meat Sci 2009,82:331–337. doi: 10.1016/j.meatsci.2009.02.002 20416720
11. Bondzio A, Gabler C, Badewien-Rentzsch B, Schulze P, Martens H, Einspanier R. Identification of differentially expressed proteins in ruminal epithelium in response to a concentrate-supplemented diet. Am J Physiol Gastrointest Liver Physiol 2011,301:G260–268. doi: 10.1152/ajpgi.00304.2010 21566014
12. Naeem A, Drackley JK, Lanier JS, Everts RE, Rodriguez-Zas SL, Loor JJ. Ruminal epithelium transcriptome dynamics in response to plane of nutrition and age in young Holstein calves. Funct Integr Genomics 2014,14:261–273. doi: 10.1007/s10142-013-0351-2 24318765
13. Egger G, Liang G, Aparicio A, Jones PA. Epigenetics in human disease and prospects for epigenetic therapy. Nature 2004,429:457–463. doi: 10.1038/nature02625 15164071
14. Su ZX, Xia JF, Zhao ZM. Functional complementation between transcriptional methylation regulation and post-transcriptional microRNA regulation in the human genome. Bmc Genomics 2011,12.
15. Li E, Bestor TH, Jaenisch R. Targeted mutation of the DNA methyltransferase gene results in embryonic lethality. Cell 1992,69:915–926. doi: 10.1016/0092-8674(92)90611-f 1606615
16. Li E, Beard C, Jaenisch R. Role for DNA methylation in genomic imprinting. Nature 1993,366:362–365. doi: 10.1038/366362a0 8247133
17. Boyes J, Bird A. DNA methylation inhibits transcription indirectly via a methyl-CpG binding protein. Cell 1991,64:1123–1134. doi: 10.1016/0092-8674(91)90267-3 2004419
18. Lorincz MC, Dickerson DR, Schmitt M, Groudine M. Intragenic DNA methylation alters chromatin structure and elongation efficiency in mammalian cells. Nat Struct Mol Biol 2004,11:1068–1075. doi: 10.1038/nsmb840 15467727
19. Eden A, Gaudet F, Waghmare A, Jaenisch R. Chromosomal instability and tumors promoted by DNA hypomethylation. Science 2003,300:455. doi: 10.1126/science.1083557 12702868
20. Klose RJ, Bird AP. Genomic DNA methylation: the mark and its mediators. Trends Biochem Sci 2006,31:89–97. doi: 10.1016/j.tibs.2005.12.008 16403636
21. Li M, Wu H, Luo Z, Xia Y, Guan J, Wang T, et al. An atlas of DNA methylomes in porcine adipose and muscle tissues. Nat Commun 2012,3:850. doi: 10.1038/ncomms1854 22617290
22. Li Q, Li N, Hu X, Li J, Du Z, Chen L, et al. Genome-wide mapping of DNA methylation in chicken. PLoS One 2011,6:e19428. doi: 10.1371/journal.pone.0019428 21573164
23. Zhang X, Yazaki J, Sundaresan A, Cokus S, Chan SW, Chen H, et al. Genome-wide high-resolution mapping and functional analysis of DNA methylation in arabidopsis. Cell 2006,126:1189–1201. doi: 10.1016/j.cell.2006.08.003 16949657
24. Eckhardt F, Lewin J, Cortese R, Rakyan VK, Attwood J, Burger M, et al. DNA methylation profiling of human chromosomes 6, 20 and 22. Nat Genet 2006,38:1378–1385. doi: 10.1038/ng1909 17072317
25. Diederich M, Hansmann T, Heinzmann J, Barg-Kues B, Herrmann D, Aldag P, et al. DNA methylation and mRNA expression profiles in bovine oocytes derived from prepubertal and adult donors. Reproduction 2012,144:319–330. doi: 10.1530/REP-12-0134 22733804
26. Ponsuksili S, Murani E, Schwerin M, Schellander K, Tesfaye D, Wimmers K. Gene expression and DNA-methylation of bovine pretransfer endometrium depending on its receptivity after in vitro-produced embryo transfer. PLoS One 2012,7:e42402. doi: 10.1371/journal.pone.0042402 22952593
27. Yang BF, Lu YJ, Wang ZG. MicroRNAs and apoptosis: implications in the molecular therapy of human disease. Clinical and Experimental Pharmacology & Physiology 2009,36:951–960.
28. Filipowicz W, Bhattacharyya SN, Sonenberg N. Mechanisms of post-transcriptional regulation by microRNAs: are the answers in sight? Nat Rev Genet 2008,9:102–114. doi: 10.1038/nrg2290 18197166
29. van Kouwenhove M, Kedde M, Agami R. MicroRNA regulation by RNA-binding proteins and its implications for cancer. Nat Rev Cancer 2011,11:644–656. doi: 10.1038/nrc3107 21822212
30. Langmead B. Aligning short sequencing reads with Bowtie. Current Protocols in Bioinformatics 2010,Chapter 11:Unit 11 17.
31. Liu T. Use model-based Analysis of ChIP-Seq (MACS) to analyze short reads generated by sequencing protein-DNA interactions in embryonic stem cells. Methods in Molecular Biology 2014,1150:81–95 doi: 10.1007/978-1-4939-0512-6_4 24743991
32. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 2014,15:550. doi: 10.1186/s13059-014-0550-8 25516281
33. Zhu LJ, Gazin C, Lawson ND, Pages H, Lin SM, Lapointe DS, et al. ChIPpeakAnno: a Bioconductor package to annotate ChIP-seq and ChIP-chip data. BMC Bioinformatics 2010,11:237. doi: 10.1186/1471-2105-11-237 20459804
34. Li Y, Carrillo JA, Ding Y, He Y, Zhao C, Zan L, et al. Ruminal Transcriptomic Analysis of Grass-Fed and Grain-Fed Angus Beef Cattle. PLoS One 2015,10:e0116437. doi: 10.1371/journal.pone.0116437 26090810
35. Li RW, Li C. Butyrate induces profound changes in gene expression related to multiple signal pathways in bovine kidney epithelial cells. BMC Genomics 2006,7:234. doi: 10.1186/1471-2164-7-234 16972989
36. Li RW, Schroeder SG. Cytoskeleton remodeling and alterations in smooth muscle contractility in the bovine jejunum during nematode infection. Funct Integr Genomics 2012,12:35–44. doi: 10.1007/s10142-011-0259-7 22203460
37. An J, Lai J, Lehman ML, Nelson CC. miRDeep*: an integrated application tool for miRNA identification from RNA sequencing data. Nucleic Acids Res 2013,41:727–737. doi: 10.1093/nar/gks1187 23221645
38. Wang HJ, Wang BZ, Zhang PJ, Deng J, Zhao ZR, Zhang X, et al. Identification of four novel serum protein biomarkers in sepsis patients encoded by target genes of sepsis-related miRNAs. Clinical Science 2014,126:857–867. doi: 10.1042/CS20130301 24303815
39. Sun J, Zhang B, Lan X, Zhang C, Lei C, Chen H. Comparative transcriptome analysis reveals significant differences in MicroRNA expression and their target genes between adipose and muscular tissues in cattle. PLoS One 2014,9:e102142. doi: 10.1371/journal.pone.0102142 25006962
40. De Smet S, Raes K, Demeyer D. Meat fatty acid composition as affected by fatness and genetic factors: a review. Animal Research 2004,53:81–98.
41. Garcia PT, Pensel NA, Sancho AM, Latimori NJ, Kloster AM, Amigone MA, et al. Beef lipids in relation to animal breed and nutrition in Argentina. Meat Science 2008,79:500–508. doi: 10.1016/j.meatsci.2007.10.019 22062910
42. Ball MP, Li JB, Gao Y, Lee JH, LeProust EM, Park IH, et al. Targeted and genome-scale strategies reveal gene-body methylation signatures in human cells. Nature Biotechnology 2009,27:361–368. doi: 10.1038/nbt.1533 19329998
43. Clark IM, Parker AE. Metalloproteinases: their role in arthritis and potential as therapeutic targets. Expert Opinion on Therapeutic Targets 2003,7:19–34. doi: 10.1517/14728222.7.1.19 12556200
44. Fernandes RJ, Hirohata S, Engle JM, Colige A, Cohn DH, Eyre DR, et al. Procollagen II amino propeptide processing by ADAMTS-3. Insights on dermatosparaxis. J Biol Chem 2001,276:31502–31509. doi: 10.1074/jbc.M103466200 11408482
45. Stefan C, Jansen S, Bollen M. NPP-type ectophosphodiesterases: unity in diversity. Trends Biochem Sci 2005,30:542–550. doi: 10.1016/j.tibs.2005.08.005 16125936
46. Stefan C, Jansen S, Bollen M. Modulation of purinergic signaling by NPP-type ectophosphodiesterases. Purinergic Signal 2006,2:361–370. doi: 10.1007/s11302-005-5303-4 18404476
47. Zimmermann H. Nucleotide signaling in nervous system development. Pflugers Arch 2006,452:573–588. doi: 10.1007/s00424-006-0067-4 16639549
48. Korekane H, Park JY, Matsumoto A, Nakajima K, Takamatsu S, Ohtsubo K, et al. Identification of Ectonucleotide Pyrophosphatase/Phosphodiesterase 3 (ENPP3) as a Regulator of N-Acetylglucosaminyltransferase GnT-IX (GnT-Vb). Journal of Biological Chemistry 2013,288:27912–27926. doi: 10.1074/jbc.M113.474304 23960081
49. Lubke T, Marquardt T, Etzioni A, Hartmann E, von Figura K, Korner C. Complementation cloning identifies CDG-IIc, a new type of congenital disorders of glycosylation, as a GDP-fucose transporter deficiency. Nat Genet 2001,28:73–76. doi: 10.1038/88299 11326280
50. Noda K, Miyoshi E, Gu JG, Gao CX, Nakahara S, Kitada T, et al. Relationship between elevated FX expression and increased production of GDP-L-fucose, a common donor substrate for fucosylation in human hepatocellular carcinoma and hepatoma cell lines. Cancer Research 2003,63:6282–6289. 14559815
51. Taniguchi N. A sugar-coated switch for cellular growth and arrest. Nature Chemical Biology 2007,3:307–309. doi: 10.1038/nchembio0607-307 17510646
52. Lee YH, Na HS, Jeong SY, Jeong SH, Park HR, Chung J. Comparison of inflammatory microRNA expression in healthy and periodontitis tissues. Biocell 2011,35: 43–49. 22128589
53. Liu FJ, Jin LJ, Ma XG, Zhang YL, Zhai XW, Chen JJ, et al. Differentially expressed microRNAs and affected signaling pathways in placentae of transgenic cloned cattle. Theriogenology 2014,82:338–346. doi: 10.1016/j.theriogenology.2014.04.010 24853279
54. Hsu SH, Wang B, Kota J, Yu JH, Costinean S, Kutay H, et al. Essential metabolic, anti-inflammatory, and anti-tumorigenic functions of miR-122 in liver. Journal of Clinical Investigation 2012,122:2871–2883. doi: 10.1172/JCI63539 22820288
55. Hou QK, Huang YQ, Zhu SL, Li PW, Chen XL, Hou ZK, et al. MiR-144 Increases Intestinal Permeability in IBS-D Rats by Targeting OCLN and ZO1. Cellular Physiology and Biochemistry 2017,44:2256–2268. doi: 10.1159/000486059 29258088
56. Blanc V, Xie Y, Kennedy S, Riordan JD, Rubin DC, Madison BB, et al. Apobec1 complementation factor (A1CF) and RBM47 interact in tissue-specific regulation of C to U RNA editing in mouse intestine and liver. Rna 2019,25:70–81. doi: 10.1261/rna.068395.118 30309881
57. Kobayashi T, Masaki T, Sugiyama M, Atomi Y, Furukawa Y, Nakamura Y. A gene encoding a family with sequence similarity 84, member A (FAM84A) enhanced migration of human colon cancer cells. International Journal of Oncology 2006,29:341–347. 16820875
58. Rimkus C, Martini M, Friederichs J, Rosenberg R, Doll D, Siewert J, et al. Prognostic significance of downregulated expression of the candidate tumour suppressor gene SASH1 in colon cancer. British Journal of Cancer 2006,95:1419–1423. doi: 10.1038/sj.bjc.6603452 17088907
59. Sun G, Zhao G, Lu Y, Wang Y, Yang C. Association of EMP1 with gastric carcinoma invasion, survival and prognosis. Int J Oncol 2014,45:1091–1098. doi: 10.3892/ijo.2014.2488 24920167
60. Chen YJ, Chang JT, Lee L, Wang HM, Liao CT, Chiu CC, et al. DSG3 is overexpressed in head neck cancer and is a potential molecular target for inhibition of oncogenesis. Oncogene 2007,26:467–476. doi: 10.1038/sj.onc.1209802 16878157
Článek vyšel v časopise
PLOS One
2019 Číslo 10
- Řekněte mi, co jíte, a já vám řeknu, jaká je úroveň vaší kognice
- Původní český test ALBA pro odhalení mírných kognitivních deficitů získal certifikaci MZ ČR
- Duševní onemocnění se u teenagerů přenáší v rámci školní třídy
- Od kostiček k lidské kůži. Jak si postavit 3D biotiskárnu?
- Prof. Libor Grubhoffer: Nová varianta mpox je příležitostí pokročit v poznání původce varioly
Nejčtenější v tomto čísle
- Correction: Low dose naltrexone: Effects on medication in rheumatoid and seropositive arthritis. A nationwide register-based controlled quasi-experimental before-after study
- Combining CDK4/6 inhibitors ribociclib and palbociclib with cytotoxic agents does not enhance cytotoxicity
- Risk factors associated with IgA vasculitis with nephritis (Henoch–Schönlein purpura nephritis) progressing to unfavorable outcomes: A meta-analysis
- Prevalence of pectus excavatum (PE), pectus carinatum (PC), tracheal hypoplasia, thoracic spine deformities and lateral heart displacement in thoracic radiographs of screw-tailed brachycephalic dogs