Integrated analyses of the methylome and transcriptome to unravel sex differences in the perirenal fat from suckling lambs

Abstract

In sheep, differences were observed regarding fat accumulation and fatty acid (FA) composition between males and females, which may impact the quality and organoleptic characteristics of the meat. The integration of different omics technologies is a relevant approach for investigating biological and genetic mechanisms associated with complex traits. Here, the perirenal tissue of six male and six female Assaf suckling lambs was evaluated using RNA sequencing and whole-genome bisulfite sequencing (WGBS). A multiomic discriminant analysis using multiblock (s)PLS-DA allowed the identification of 314 genes and 627 differentially methylated regions (within these genes), which perfectly discriminate between males and females. These candidate genes overlapped with previously reported QTLs for carcass fat volume and percentage of different FAs in milk and meat from sheep. Additionally, differentially coexpressed (DcoExp) modules of genes between males (nine) and females (three) were identified that harbour 22 of these selected genes. Interestingly, these DcoExp were significantly correlated with fat percentage in different deposits (renal, pelvic, subcutaneous and intramuscular) and were associated with relevant biological processes for adipogenesis, adipocyte differentiation, fat volume and FA composition. Consequently, these genes may potentially impact adiposity and meat quality traits in a sex-specific manner, such as juiciness, tenderness and flavour.

Keywords: RNA sequencing; fat deposits; omics integration; perirenal fat; suckling lambs; whole-genome bisulfite sequencing.

FIGURE 1

Flowchart summarizing the methods applied and the respective objectives reached in the current study. The grey boxes highlight the main results obtained in the methodological pipeline designed and applied in the present study.

FIGURE 2

Distribution and genomic density of differentially methylated loci (DMLs) across the genome of Assaf lambs in the three evaluated methylation contexts CG (A), CHG (B), and CHH (C). On the right-hand side, the circular Manhattan plots show the distribution of–DMLs in the comparison between males and females (log₁₀ (p-values) after false discovery rate adjustment for the detected). On the left-hand side, the bar plots show the density of DMLs within 1 Mb windows for each sheep chromosome. A darker red shade indicates a higher number of DMLs within the 1 Mb windows.

FIGURE 3

Distribution of (A) mean methylation level (B) length, and (C) genomic context (exon, intergenic, intron or promoter region) for the differential methylation regions (DMRs) detected between males and females in Assaf suckling lambs. (A) Violin plots showing the distribution of mean methylation levels within DRMs in the three methylation contexts (CG, CHG and CHH) and between males (in blue) and females (in red). (B) Violin plots showing the distribution of the DMR lengths in the three methylation contexts (CG, CHG and CHH) and between males (in blue) and females (in red). The length of DMRs was truncated at 500 bp for all three contexts to obtain a clear visualization of the distributions. (C) Pie plots showing the percentage of DMRs detected within each genomic context (intergenic, exon, intron and promoter).

FIGURE 4

Multi-omics discriminant analysis between males (blue) and females (orange) samples using multiblock (s)PLS-DA. (A) Scatterplot from plotDiablo displaying the first principal component in RNA-Seq and WGBS datasets (upper diagonal plot), thus, constructed with expression values of 314 discriminant genes between males and females (RNA-Seq data) and 627 mean methylation levels within DMRs mapped within the same genes (WGBS data), and Pearson correlation (0.98) between each component (lower diagonal plot) (B) Loading values for the top 20 DMRs (left-hand side) and gene expression values (right-hand side) obtained in the multiblock (s)PLS-DA used to discriminate male and female samples.

FIGURE 5

Functional analysis for the 314 discriminant genes obtained in the (s)PLS-DA discriminant analysis between male and female samples. (A) Functional grouping tree diagram for the Gene Ontology (GO) terms annotated for the 314 genes selected after the multiblock (s)PLS-DA discriminant analysis. Each color in the dendrogram represents a functional group obtained after estimating the Jaccard correlation coefficient. The area of the circles represents the number of genes assigned to each GO term, and the color of the circle indicates the p-value estimated for each GO term. (B) Pie plot showing the percentage of each QTL type annotated within the genomic coordinate of the 314 discriminant genes and bubble plot with QTL enrichment analysis results, where the darker the shade of red, the smallest is the enrichment p-value and the area of the circle represents the number of annotated QTLs.

FIGURE 6

Pearson correlation between module eigengenes and fat percentage in different adipose deposits (Renal, Pelvic, Renal and Pelvic, Leg subcutaneous, and Leg intramuscular). The modules shown in the figure are the differentially coxpressed modules harboring at least one of the 314 discriminant genes between males and females. For each module, after the color identification of the module, it is indicated if the coexpressed gene module was identified in males or females. The Pearson correlation between the fat percentage and the module is shown in the figure for each specific fat deposit, and the p-value assigned for each correlation coefficient is shown between parenthesis (^** p-value<0.05). The color scale represents the signal of the correlation coefficient, where green shades represent negative correlations and red shades positive correlations.