Descriptive statistics and genomic heritability estimates
Summary statistics and genomic heritability estimates for the 81 fatty acid traits
of SQ and 83 fatty acid traits of LL are presented in Table 1. In general, the estimates of heritability for the same fatty acids are comparable
in both the adipose and muscle tissues, with a correlation coefficient of 0.61. Relatively
higher (0.40) heritability estimates were found for 10:0, 12:0, 18:0, ai15:0, 9c-14:1,
9c-16:1, 13c-18:1, 18:3n-3, 18:2n-6, n-3, n-6, sumtrans 18:1, total PUFA, P/S, and
P/(S?+?B) in the SQ tissue, and for 12:0, 14:0, 16:0, 9c-14:1, 9c-16:1, 12c-16:1,
9c-18:1, 13c-18:1, SFA, SFA?+?BFA, MUFA, n-6/n-3, and health index (HI) in the LL
muscle, which suggests greater direct host genetic effects on these traits in the
corresponding tissues. Very low (0.05) or zero heritability were observed for 22:0,
7c-17:1, 12Â t-18:1, 15Â t-18:1, 6Â t,8Â t:18:2, 7Â t,9Â t-18:2, 9Â t,11Â t-18:2, 10Â t,12Â t-18:2,
12Â t,14Â t-18:2, and n-6/n-3 in the SQ tissue, and for 7c-17:1, 15Â t-18:1, 6Â t,8Â t-18:2,
7Â t,9Â t-18:2, 7Â t,9c-18:2, 8Â t,10Â t-18:2, 12Â t,14Â t-18:2 in the LL muscle, which indicates
weak host direct genetic control on these traits. In general, the heritability estimates
for the fatty acids in this study are in line with those reported in other studies
10], 12], 15]. Therefore, the genomic estimated additive genetic variance and heritability for
the fatty acid traits were further used in the calculation of total genetic variance
explained by significant markers identified in GWAS and in the derivation of realised
accuracy of genomic prediction in this study.
Table 1. Summary statistics of mean, standard deviation (SD), additive genetic variance (? a2
), and heritability estimates (h
2
?±?SE)
Genome-wide association study
In total, 302 and 360 significant SNPs spanning all autosomal chromosomes were identified
to be associated with one or more fatty acid traits in the SQ and LL tissues, respectively,
at the genome-wise empirical significance threshold at ??=?0.05. Significant SNPs
and their distributions over the genome varied for different fatty acid traits. Manhattan
plots of posterior probability of inclusion (PPI) were provided in Additional file
1 for all fatty acid traits in the two tissues. Proportions of genotypic variance explained
by individual significant SNPs ranged from 0.03 to 11.06Â % in SQ, and from 0.005 to
24.28Â % in LL. Among these, 28 and 41 SNPs individually explained greater than 1Â %
of total genetic variance for at least one fatty acid trait in the SQ and LL tissues,
respectively. Figures 1 and 2 showed these SNPs and their associated traits in SQ and LL, respectively. Of these
SNPs, SNP rs41921177 at the location of BTA19:51326750 had the largest effects on multiple fatty acid
traits in both tissues, followed by SNP rs42714483 at BTA29:18090509 and SNP rs42090719 at BTA26:20903573. Details of all significant SNPs including SNP name, chromosome
position, allele substitution effect, percentage of total genetic variance explained,
and PPI were also provided in additional files (Additional file 2 for SQ and Additional file 3 for LL). Candidate genes within 1 mega base pair (Mb) region centering the significant
SNPs were provided separately (Additional file 4).
Fig. 1. Summary of fatty acid trait associations across genomic regions (SNPs) and percentage
of genetic variance explained by significant SNPs in the subcutaneous adipose tissue
(SQ). Each row represents a trait and each column represents a SNP. Only traits with
at least one significant SNP explaining greater than 1Â % of genetic variance were
listed, and only SNPs that explain greater than 1Â % of genetic variance for at least
one trait were shown. The top of the figure shows the chromosome and position and
the bottom shows the name of the SNP. Va %: percentage of genetic variance
Fig. 2. Summary of fatty acid trait associations across genomic regions (SNPs) and percentage
of genetic variance explained by significant SNPs in longissimus lumborum muscle (LL). Each row represents a trait and each column represents a SNP. Only traits
with at least one significant SNP explaining greater than 1Â % of genetic variance
were listed, and only SNPs that explain greater than 1Â % of genetic variance for at
least one trait were shown. The top of the figure shows the chromosome and position
and the bottom shows the name of the SNP. Va %: percentage of genetic variance
SNP rs41921177 was significantly associated with 19 individual and grouped fatty acids in the LL
muscle including SFAs 10:0, 12:0, 13:0, 14:0, 15:0, 16:0, 18:0, branched fatty acids
(BFAs) ai 15:0 and iso 18:0, MUFAs 9c-14:1, 9c-15:1, 9c-16:1, 12c-16:1, 9c-18:1, 11c-20:1,
grouped fatty acids total SFA, SFA?+?BFA, total MUFA and HI, with genetic variance
explained from 1.37Â % (18:0) to 24.28Â % (14:0). The same SNP also showed significant
associations with 11 of the above fatty acids in SQ including saturated fatty acids
12:0, 14:0, 15:0, 16:0, SFA, monounsaturated fatty acids 9c-14:1, 9c-15:1, 9c-16:1,
9c-18:1, 12c-16:1 and HI, explaining 0.33Â % (SFA) to 11:06Â % (14:0) of the genetic
variance (also see Additional file 2). This chromosomal region was previously identified to be associated with 14:0, 16:0,
16:1 and 18:1 in adipose and muscle tissues of a Jersey and Limousin crossbred beef
cattle 21], with 9c-18:1, and 14:0, 14:1, 16:0, 16:1 in intramuscular fat of Japanese Black
cattle 7], 37], with 14:0, 16:0, 9c-14:1 and 9c-18:1 in adipose tissue of an Australian multi-breed
beef population 38], with 14:0, 14:1, 16:0, 16:1, 9c-18:1, MUFA, SFA, and Atherogenic index (AI, the
inverse of HI) in muscle of American Angus beef cattle 16]. The association of this chromosomal region with the fatty acid traits was therefore
confirmed in both the SQ and LL tissues of a Canadian beef population of diverse breed
compositions, indicating a strong host genetic effect on the fatty acid composition
in beef tissues. Multiple genes are within 1Â Mb region centering the SNP (see Additional
file 4), with FASN being a strong candidate gene due to its function in fatty acid synthesis 40], 41]. Different SNPs of the FASN gene have also been reported to be associated with concentrations of saturated and
monounsaturated fatty acids in various beef and dairy cattle populations 7], 16], 18], 22], 37], 38], 40]–45].
SNP rs42714483 showed significant associations with concentrations of 15 fatty acids in the LL tissue
and 10 fatty acids in SQ including 10:0, 12:0, 13:0, 14:0, 15:0, 9c-14:1, 12c-16:1,
13c-18:1, 9c,15c-18:2, and HI in both the tissues, and 16:0, 18:0, 9c-15:1, 9c-16:1,
and 9c-18:1 in the LL tissue. Saatchi et al.16] also identified the same chromosomal region associated with fatty acids 14:0, 9c-14:1,
16:0, 16:1, 18:0, 9c-18:1, and AI, and Kelly et al.38] found SNPs in the same chromosomal region that were associated with fatty acids 14:0,
9c-14:1 in subcutaneous adipose tissue of an Australian multi-breed beef population
38]. These results strongly support that the chromosome region on BTA 29 harbors host
genes that influence fatty acid composition of beef tissues. In this study, the SNP
at BTA29:18090509 is a missense mutation (T/C) of the thyroid hormone responsive gene (THRSP), causing amino acid change from isoleucine to valine (I16V). Recently THRSP has been considered as a candidate gene for fatty acid composition in beef 16], 46]. Substitution of allele T with C of this missense mutation was associated with decrease of 10:0, 12:0, 13:0, 14:0,
15:0, 16:0, 9c-14:1, 9c-15:1, 9c-16:1, 12c-16:1, 13c-18:1, 9c,15c-18:2, and increase
of 18:0, 9c-18:1, and HI (see Additional files 2 and 3). The direction of the allele substitution effect on different fatty acid traits
also coincided with that of SNP rs41921177, which is close to FASN gene, suggesting possible co-ordinations between THRSP and FASN genes in fatty acid synthesis.
SNP rs42090719 at BTA26: 20903573 was found to be significantly associated with 9c-14:1, 12c-16:1,
13c-18:1 in both the LL and SQ tissues, 9c,15c-18:2 and CLA isomers 11Â t,13c?+?11c,13Â t-18:2
(also see Additional file 3) in the LL tissue. In addition, SNP rs41646463 at BTA26:21258113 also showed significant associations with 13c-18:1 in both the
LL and SQ tissues. In the nearby chromosomal region of BTA26:18994785, SNP rs109465094 was significantly associated with 9c-14:1, 12c-16:1, 13c-18:1 and 9c,15c-18:2 in
LL. The chromosomal regions on BTA 26 were previously found associated with a variety
of fatty acids in muscle of American Angus and in adipose tissue in an Australian
multi-breed beef population, and stearoyl-CoA desaturase gene (SCD) was suggested as a candidate gene 16], 38]. The two SNPs, rs42090719 and rs41646463, are within 250 kilo base pairs (Kb) of SCD. The other SNP rs109465094 is more than 2Â Mb distant from SCD, indicating a possible alternative candidate gene or its association could merely
be due to LD with SCD. Linkage disequilibrium between SNPs around the SCD gene were analysed and visualised
using the Haploview software 47] and results are shown in Additional file 5. Indeed, the three significant SNPs are in moderate to high LD with SNPs in a LD
block containing the SCD gene. The SCD gene is involved in the synthesis of particular MUFA and CLA isomers, in creating
a double bond at the ?
9
position of fatty-acyl CoA 48], 49]. SCD has been reported to be associated with both meat and milk fatty acid composition
in cattle 7], 16], 18], 32], 34], 37], 38], 42], 50]–56]. The present study showed that SNPs close to the SCD gene were associated with many MUFAs and several CLA isomers but none of the SFAs,
which supports the proposed role of SCD in fatty acid composition in beef. However, in this study none of the SNPs around
SCD were associated with oleic acid, 9c-18:1, the most abundant MUFA in beef. This could
be partly due to lack of SNPs in the current panel that are in a high LD with SCD to capture all its effects. Interestingly, several other studies also showed no associations
between SCD and oleic acid in various beef and dairy cattle populations, using different SNP
panels or SCD gene SNP 7], 16], 38], 50], 51], 57], although two other studies have reported significant associations between SCD SNP variants and oleic acid concentrations in Japanese Black cattle 18], 32]. The role of SCD on the concentration of oleic acid in beef is worthy of further investigation.
Other SNPs on BTA 1, 3, 4, 5, 6, 7, 10, 16, 20, 23, 24, 25, 28 and on 1, 2, 4, 5,
6, 8, 10, 13, 14, 16, 23, 24, 25, 27, 28 were found significantly associated with
one or more fatty acid concentrations in the SQ and LL, respectively, but with relatively
smaller effects (Figs. 1 and 2). The SNP rs41642879 at BTA23:6760915 was associated with 12:0, 14:0, and HI of both the LL and SQ, and
with 16:0, 9c-14:1 in LL tissue. There are several genes within the 1Â Mb window centering
the SNP. One possible candidate gene is the glutamate-cysteine ligase catalytic subunit
gene (GCLC), which is involved in the synthesis of glutathione (GSH) 58]. Glutathione has an antioxidant function by oxidizing itself into Glutathione disulfide
(GSSG), which in turn is reduced to GSH at the expense of nicotinamide adenine dinucleotide
phosphate (NADPH) oxidase 58]. The latter is essential for fatty acid synthesis 40]. The SNP rs41574597 at BTA23:25956421 was found to be associated with 11c-20:1 in both tissues. Several
genes belonging to the butyrophilin family are located nearby. Butyrophilin is the
major protein associated with milk fat droplets and has been reported to be related
to milk quality in cattle 59]. However, it was suggested that butyrophilin is specific to mammary tissue 60] hence its role in meat fatty acid production remains unclear.
Fewer SNPs were identified for PUFAs (90 in SQ and 87 in LL) in comparison to the
number of significant SNPs for SFAs (121 in SQ and 117 in LL) and MUFAs (174 in SQ
and 120 in LL). One SNP rs41582945 at BTA2:50069820 explained 2.45Â % of genetic variance for dihomo-gamma-linolenic
acid (Dihomo-GLA, 20:3n-6) in LL. However, no known genes exist in the 1Â Mb region
of this SNP. Several SNPs were also found to be associated with the intermediate product
of Dihomo-GLA, arachidomic acid in LL (20:4n-6). The most significant SNP rs110776216 on BTA2:131968682 explained 1.91Â % of total genetic variance of 20:4n-6 and was located
within the endothelin converting enzyme 1 gene (ECE1) which encodes the enzyme that converts big endothelin-1 to endothelin-1. Endothelin-1
was previously found to stimulate arachidonic acid release in human pericardial smooth
muscle cells 61], 62]. In this study, no significant SNPs were found to be associated with iso14:0, 7c-17:1,
15Â t-18:1, 9c-20:1, 6Â t,8Â t-18:2, 7Â t,9c-18:2, 9Â t,11Â t-18:2, 8Â t,10Â t-18:2, 7Â t,9Â t-18:2,
20:5n3 in LL and 22:0, 7c-17:1, 11c-18:1, 12Â t-18:1, 15Â t-18:1, 6Â t,8Â t-18:2, 10Â t,12Â t-18:2,
9Â t,11Â t-18:2, 8Â t,10Â t-18:2, 7Â t,9Â t-18:2, 18:3n6, 20:3n9, and n-6/n-3 in SQ. These
fatty acid traits had very low or near zero heritability estimates (Table 1), therefore their concentrations were less likely influenced by host direct genetic
effects.
Genomic prediction
Realized accuracies of genomic prediction measured as the Pearson’s correlation coefficients
between genomic estimated breeding values (GEBV) and adjusted phenotypic values of
fatty acid traits divided by square root of heritability are presented in Table 2. Accuracies of breeding values estimated from the pedigree-based BLUP method (PBLUP)
are also presented in Table 2 as comparisons. The realized accuracy of genomic prediction ranged from ?0.05 for
15Â t-18:1 to 0.73 for 14:0 in the LL muscle, and varied from ?0.05 for 20:3n-9 to
0.65 for 16Â t-18:1 in the SQ tissue. Averaged across all traits, accuracies from PBLUP,
GBLUP, and the Bayesian method were 0.23, 0.32, and 0.35, respectively, in SQ, and
0.17, 0.39, and 0.46, respectively, in LL. These results suggested the effectiveness
of genomic prediction using either GBLUP or the Bayesian method. However, the incompleteness
of the pedigree (only one generation) may largely contribute to the low accuracy for
the PBLUP method. It should be noted that the realized accuracy could be overestimated
when heritability is underestimated as pointed out by Lourenco et al.63]. Accuracies that were substantially overestimated tended to have relatively large
SE (0.10) as shown in Table 2. Additionally, Pearson’s correlation coefficient between estimated breeding values
and adjusted phenotypes, and regression coefficient by regressing adjusted phenotypes
on estimated breeding values were also calculated and provided in Additional file
6. The correlation coefficients averaged 0.11, 0.15, and 0.15 for PBLUP, GBLUP, and
the Bayesian method, respectively, in SQ, and averaged 0.08, 0.14, and 0.16, respectively
in LL. The average regression coefficients in SQ were 1.02, 0.77, and 0.92, and were
1.04, 0.90, and 0.83 in LL for PBLUP, GBLUP, and the Bayesian method, respectively.
The regression coefficient is expected to be 1 if the estimated breeding values were
unbiased predictions of the true breeding values. Nevertheless, for most of the fatty
acid traits, the accuracy of genomic prediction were relatively low (0.40), which
was expected given the low heritability estimates and the small sample size used in
this study 64]. Relatively higher accuracy (r(GEBV,y)
/h???0.50 with SE??0.10) were achieved for 10:0 (0.53), 12:0 (0.53), 14:0 (0.73), 15:0
(0.69), 16:0 (0.50), 9c-14:1 (0.55), 12c-16:1 (0.55), 13c-18:1 (0.51), and HI (0.59)
in LL, and for 12:0 (0.58), 14:0 (0.61), 15:0 (0.62), 10Â t,12c-18:2 (0.52), and 11Â t,13c?+?11c,13Â t-18:2
(0.56) in SQ. The relatively higher accuracy for certain saturated and monounsaturated
fatty acids, and HI, and relatively lower accuracy for CLAs and other PUFAs in muscle
were compatible with the magnitude of their estimated heritability (Table 1). The correlations between heritability estimates and realised accuracy of genomic
prediction in LL were 0.61 and 0.39 for Bayesian and GBLUP methods, respectively.
However, in SQ such correlations were only 0.10 for GBLUP and 0.23 for the Bayesian
method, which is likely due to many overestimations of realised accuracy for traits
with low and inaccurate heritability estimates. Genomic prediction from the Bayesian
method performed similarly as GBLUP for most of the traits, but substantially better
for several traits in LL muscle such as 10:0 (0.37 for GBLUP vs 0.53Â for BayesC?),
12:0 (0.31 vs 0.53), 14:0 (0.45 vs 0.73), 15:0 (0.57 vs 0.69), 16:0 (0.36 vs 0.50),
9c-14:1 (0.34 vs 0.55), 9c-16:1 (0.37 vs 0.49), 12c-16:1 (0.32 vs 0.55), 9c-18:1 (0.27
vs 0.37), 13c-18:1 (0.36 vs 0.51), and HI (0.41 vs 0.59), and for traits in SQ including
12:0 (0.42 vs 0.58), 14:0 (0.39 vs 0.61), and 9c-14:1 (0.31 vs 0.43). These traits
have been shown to have SNPs with larger effects from GWAS results (Figs. 1 and 2). The Bayesian method adopted in this study allows a fraction of SNPs to take relatively
large effects, which may better characterize the genetic architecture of traits that
have QTL of larger effects than the GBLUP method 65], which assumes all SNPs have the same genetic variance.
Table 2. Realised accuracy (±SE) of breeding value prediction for fatty acid traits in the
subcutaneous adipose and longissimus lumborum musle
Fatty acid composition is a complex trait and it is difficult and expensive to measure,
making it a good candidate trait for genomic selection. To date, genomic prediction
for fatty acid composition in beef cattle has only been reported by Saatchi et al.16] for 24 individual and grouped/ratio of fatty acids in steaks of American Angus beef
cattle, and by Onogi et al.39] for 8 fatty acid traits in musculus trapezius of Japanese Black cattle. Relatively higher prediction accuracies were found for
14:0 (0.57), 16:0 (0.53), total long chain saturated fatty acids (0.57), total medium
chain saturated fatty acids (0.57), 9c-18:1 (0.35), 12c-18:1 (0.35), total MUFA (0.38),
(14:0?+ 16:0)/all (0.55), and AI (0.56) in Saatchi’s study, compared to other fatty
acid traits. In this study, relatively higher accuracies were also obtained for SFAs
12:0, 14:0, and 15:0 in both the LL and SQ tissues, and for 10:0, 16:0, 9c-14:1, 12c-16:1,
13c-18:1, and HI in LL (Table 2), suggesting strong host genetic controls on synthesis of these SFAs and MUFAs. Saatchi
et al.16] reported genomic prediction accuracies for 12 PUFAs and all were very low (0.30).
In this study, we analyzed 32 and 34 PUFAs and PUFA-BHI (including CLAs and 11Â t-18:1)
in the adipose and muscle tissues, respectively, and found moderate accuracies (between
0.30 and 0.45) for 11Â t-18:1, 9c,13Â t?+?8Â t,12c-18:2, 9c,15c-18:2, 8Â t,13c-18:2, 11Â t,15c-18:2,
18:2n-6, 18:3n-3, n-3, n-6, and total PUFA in both the adipose and muscle tissues,
moderate accuracies for 20:3n-6, 20:3n-9, 22:4n-6, 22:6n-3 in the muscle, and relatively
high accuracies (0.50) for 12Â t,14c?+?12c,14Â t-18:2, and 11Â t,13c?+?11c,13Â t-18:2
in the adipose tissue, suggesting considerable host genetic influence on these fatty
acids. Different beef cattle populations, environments where the animals were raised,
sample sizes and statistical models may also contribute to the differences of genomic
prediction accuracy observed between different studies. Although most dietary PUFAs
are biohydrogenated by rumen bacteria 66], a portion of PUFAs and PUFA-BHI may escape and deposit into body fat of beef. In
addition, some PUFAs can be endogenously synthesized, for example CLAs can be synthesized
from one of the PUFA-BHI, vaccenic acid (11Â t-18:1) by the host 67]. Therefore, contents of both PUFAs and PUFA-BHI are potentially influenced by host
genetics and thus predictable by genomic prediction. Onogi et al. 39] also reported a relatively high accuracy (0.56) for PUFA C18:2 in Japanese Black cattle.
Although it would be worthwhile to further verify the genomic prediction accuracy
in other beef cattle populations, the moderate to relatively high genomic prediction
accuracies achieved in this study for the HI, several individual SFAs, MUFAs, PUFAs
and PUFA-BHI suggest that genomic selection is a promising tool for genetic improvement
of fatty acid profiles in beef cattle to produce healthier meat. Therefore, as consumers’
demand for healthier meat continues to grow, beef producers may get more premiums
by producing meat with enhanced fatty acid profiles, which can be achieved by incorporating
fatty acid composition traits into a multi-trait selection index for selection and/or
by genetic based diet management.