Multi-omic data integration for the study of production, carcass, and meat quality traits in Nellore cattle

Abstract

Data integration using hierarchical analysis based on the central dogma or common pathway enrichment analysis may not reveal non-obvious relationships among omic data. Here, we applied factor analysis (FA) and Bayesian network (BN) modeling to integrate different omic data and complex traits by latent variables (production, carcass, and meat quality traits). A total of 14 latent variables were identified: five for phenotype, three for miRNA, four for protein, and two for mRNA data. Pearson correlation coefficients showed negative correlations between latent variables miRNA 1 (mirna1) and miRNA 2 (mirna2) (-0.47), ribeye area (REA) and protein 4 (prot4) (-0.33), REA and protein 2 (prot2) (-0.3), carcass and prot4 (-0.31), carcass and prot2 (-0.28), and backfat thickness (BFT) and miRNA 3 (mirna3) (-0.25). Positive correlations were observed among the four protein factors (0.45-0.83): between meat quality and fat content (0.71), fat content and carcass (0.74), fat content and REA (0.76), and REA and carcass (0.99). BN presented arcs from the carcass, meat quality, prot2, and prot4 latent variables to REA; from meat quality, REA, mirna2, and gene expression mRNA1 to fat content; from protein 1 (prot1) and mirna2 to protein 5 (prot5); and from prot5 and carcass to prot2. The relations of protein latent variables suggest new hypotheses about the impact of these proteins on REA. The network also showed relationships among miRNAs and nebulin proteins. REA seems to be the central node in the network, influencing carcass, prot2, prot4, mRNA1, and meat quality, suggesting that REA is a good indicator of meat quality. The connection among miRNA latent variables, BFT, and fat content relates to the influence of miRNAs on lipid metabolism. The relationship between mirna1 and prot5 composed of isoforms of nebulin needs further investigation. The FA identified latent variables, decreasing the dimensionality and complexity of the data. The BN was capable of generating interrelationships among latent variables from different types of data, allowing the integration of omics and complex traits and identifying conditional independencies. Our framework based on FA and BN is capable of generating new hypotheses for molecular research, by integrating different types of data and exploring non-obvious relationships.

Keywords: Bayesian network; factor analysis; latent variables; meat quality; omics data.

FIGURE 1

Correlation plot of 22 phenotypes. The degree of shading and the value reported correspond to the correlations among the traits. BWi: pre-feedlot body weight; REAi: initial ribeye area by ultrasonography; REAf: final ribeye area on steak; BFTi: initial backfat thickness by ultrasonography; BFTf: final backfat thickness by carcass; fat_pelvis: pelvis fat content at carcass; fat_kidney: kidney fat content; carcass_hot: hot carcass weight; carcass_cold: cold carcass weight; carcass_depth: carcass depth; pH: pH at 24 h; water_free: free water; w_ret_cap: water holding capacity; moisture: meat moisture; SF: shear-force; MFI: miofibrilar fragmentation index; L*, a*, b*: color parameters; and IMF: intramuscular fat.

FIGURE 2

Principal component analysis of total counts as a batch effect in miRNAs. Principal component analysis of miRNAs before (A) and after (B) the limma batch effect normalization. The total counts refer to the total number of reads per sample. Three colors were used to represent 1) samples with a higher total number of reads: higher than mean +standard deviation (345,281 reads) (blue); 2) samples with a lower total number of reads: lower than mean—standard deviation (125,193 reads) (red); 3) samples with the average total number of reads: between mean > + standard deviation and mean < + standard deviation (green).

FIGURE 3

Principal component analysis of proteomic data. (A) Principal component analysis with all animals normalized together. (B) Principal component analysis when animals were normalized separately by equipment acquisition. The colors denote the equipment batch effect. The samples in red, green, dark blue, and blue colors are from equipment 1. The samples in pink are from equipment 2.

FIGURE 4

Final underlying latent structures of phenotypes generated by exploratory factor analysis. BWi: initial body weight on feedlot trial; REAi: initial ribeye area by ultrasonography; REAf: final ribeye area on steak; BFTi: initial backfat thickness by ultrasonography; BFTf: final backfat thickness by carcass; SF: shear-force; MFI: myofibrillar fragmentation index; and L*, a*, b*: color parameters.

FIGURE 5

Final underlying latent structures of miRNA (yellow), proteins (green), and mRNA (blue) generated by exploratory factor analysis. N denotes the total number of features in each latent variable.

FIGURE 6

Correlation plot of 14 factor scores. The degree of shading and the value reported correspond to the correlation between each pair of latent variables.

FIGURE 7

Bayesian network between latent variables based on the score-based (hill climbing and tabu) algorithms. The quality of the structure was evaluated by bootstrap resampling and model averaging across 1,000 replications. Orange nodes: phenotype latent variables; yellow nodes: miRNA latent variables; green nodes: protein latent variables; blue node: gene expression of mRNA1 module (WGCNA); and white node: gene expression of mRNA2 module (WGCNA). The labels of the arcs correspond to the strength and direction (in parenthesis).