DNA methylation may affect beef tenderness through signal transduction in Bos indicus

Abstract

Background: Beef tenderness is a complex trait of economic importance for the beef industry. Understanding the epigenetic mechanisms underlying this trait may help improve the accuracy of breeding programs. However, little is known about epigenetic effects on Bos taurus muscle and their implications in tenderness, and no studies have been conducted in Bos indicus.

Results: Comparing methylation profile of Bos indicus skeletal muscle with contrasting beef tenderness at 14 days after slaughter, we identified differentially methylated cytosines and regions associated with this trait. Interestingly, muscle that became tender beef had higher levels of hypermethylation compared to the tough group. Enrichment analysis of predicted target genes suggested that differences in methylation between tender and tough beef may affect signal transduction pathways, among which G protein signaling was a key pathway. In addition, different methylation levels were found associated with expression levels of GNAS, PDE4B, EPCAM and EBF3 genes. The differentially methylated elements correlated with EBF3 and GNAS genes overlapped CpG islands and regulatory elements. GNAS, a complex imprinted gene, has a key role on G protein signaling pathways. Moreover, both G protein signaling pathway and the EBF3 gene regulate muscle homeostasis, relaxation, and muscle cell-specificity.

Conclusions: We present differentially methylated loci that may be of interest to decipher the epigenetic mechanisms affecting tenderness. Supported by the previous knowledge about regulatory elements and gene function, the methylation data suggests EBF3 and GNAS as potential candidate genes and G protein signaling as potential candidate pathway associated with beef tenderness via methylation.

Keywords: Cattle; EBF3; Epigenome; GNAS; Methylation; Muscle; Nelore; RRBS; Shear force.

Fig. 1

Chromosome distribution of A 123 differentially methylated cytosines (DMCs) and B 42 differentially methylated regions (DMRs) between two Nelore groups of divergent phenotypes for shear force at 14 days of aging. Only DMCs and DMRs showing q value < 0.05 and differences in methylation > 25% between groups are shown. The asterisks at the top of the represented DMCs and DMRs show which methylation feature overlapped with tenderness-related QTLs found in our population or QTL database. Blue bars represent DMCs or DMRs that were hypermethylated in animals with tender beef, while red bars represent those hypomethylated in the same group. The color gradient represents the difference of methylation (%) between the divergent groups

Fig. 2

Overview of the over-represented biological processes (BP) annotations for the predicted target genes of differentially methylated CpGs (DMCs) and regions (DMRs). The bar plot displays result from a functional enrichment analysis to identify BP gene ontology (GO) terms. Significant ontology terms were identified at an estimated false discovery rate (FDR) < 0.05. The colored bars represent the percentage of genes per GO BP term and the number of genes associated with the term is shown as a label. Ontology terms that were closely related based on similar genes/biological functions were grouped by kappa score (kappa score > 0.4) and represented within the figure as the same color

Fig. 3

Significant Pearson coefficient correlation between the methylation percentage of DMCs and DMRs with the expression of its target genes (p < 0.05), using 11 RNA-Seq Nelore muscle samples. A Percent of methylation of DMC23 was positively correlated with the expression of gene GNAS. B Percent of methylation of DMC98 was negatively correlated with the expression of gene PDE4B. Percent of methylation of DMC89 (C), DMC90 (D) and DMR40 (E) were negatively correlated with expression of EBF3 gene. F Percent of methylation of DMC12 was negatively correlated with the expression of gene EPCAM

Fig. 4

Regulatory element features of the region overlapped by the differentially methylated cytosine DMC23, which was correlated with GNAS expression, in adipose (A) and muscle (M). DMC23 (represented by the vertical black line) overlapped a CpG island (CpG358, represented by the dark green solid horizontal bar) in the intron 1 of GNAS isoform 001,271,771 and less than 1 kb from the start site of isoform 181,021. Histone marks and ATAC peaks enrichment suggested that this region was classified as bivalent/poised TSS (represented by the orange solid horizontal bar #12) in muscle (M) and strongly active promoter/transcript in adipose (A; red solid horizontal bar #1) in two male Bos taurus [23]. The GNAS gene is represented by the blue line (introns) and blue blocks (exons). Pink tracks represent the ATAC peaks. The peaks of histone marks H3K27ac, H3K27me, H3K4me1 and H3K4me3 are represented by the red, black, yellow, and green tracks, respectively. The image was obtained from UCSC Genome Browser and edited by the authors

Fig. 5

Regulatory element feature of the region overlapped by the two differentially methylated cytosines (DMC89 and DMC90) and region (DMR40) which were correlated with EBF3 expression. DMCs and DMRs region (represented by the vertical blue line) overlapped a CpG island (CpG1144, represented by the green solid horizontal bar). The DMCs and the DMR overlap a region enriched with histone marks and ATAC peaks, which characterizes the region as strongly active promotors/transcripts in muscle and bivalent/poised TSS in adipose in two male Bos taurus [23]. The EBF3 gene is represented by the dark red line (introns) and dark red blocks (exons). Pink tracks represent the ATAC peaks. The peaks of histone marks H3K27ac, H3K27me, H3K4me1 and H3K4me3 are represented by the red, black, yellow, and green tracks, respectively. The image was obtained from UCSC Genome Browser and edited by the authors