Introduction
We have recently witnessed a growing interest in genomic selection programmes aiming to alter milk-fat composition, due to the growing awareness that fatty acid intake has a negative effect on serum lipids, and that milk and milk-derived foods are the main source of dietary saturated fatty acids across Europe (Hulshof et al. Reference Hulshof, van Erp-Baart, Anttolainen, Becker, Church, Couet, Hermann-Kunz, Kesteloot, Leth, Martins, Moreiras, Moschandreas, Pizzoferrato, Rimestad, Thorgeirsdottir, van Amelsvoort, Aro, Kafatos, Lanzmann-Petithory and van Poppel1999). For making predictions, in the context of marker assisted selection, results of the genomic association studies can only be applied to the population under analysis. As far as we know, genomic association studies of milk fatty acids (FA) have, to date, been performed only on first lactation Holstein cows. Shennink et al. (Reference Schennink, Stoop, Visker, Heck, Bovenhuis, van der Poel, van Valenberg and van Arendonk2007, Reference Schennink, Heck, Bovenhuis, Visker, van Valenberg and van Arendonk2008) and Conte et al. (Reference Conte, Mele, Chessa, Castiglioni, Serra, Pagnacco and Secchiari2010) showed that single nucleotide polymorphisms (SNP) in the stearoyl-CoA desaturase (SCD), diacylglycerol acyltransferase 1 (DGAT1), and sterol regulatory element binding protein (SREBP-1) genes explained part of the genetic variation of milk fat unsaturation indices. Stoop et al. (Reference Stoop, Schennink, Visker, Mullaart, van Arendonk and Bovenhuis2009) and Bouwman et al. (Reference Bouwman, Bovenhuis, Visker and van Arendok2011) performed a genome-wide scan to detect Quantitative Trait Loci (QTL) for milk FA, and identified two QTL, for short and medium-chain FA, on BTA14 and BTA26, most likely caused by mutations in the SCD and DGAT1 genes. The same authors also reported suggestive QTL for short and medium-chain FA on 21 autosomes. In the mentioned examples, QTL were detected in first lactation cows of the Holstein breed, where one single milk sample per lactation had been collected. QTL studies try to identify candidate polymorphic genes that explain the variability of the trait under study. Once the candidate gene has been identified, and once the effect of the polymorphic alleles has been assessed in different populations, it can be used to make predictions in the context of marker assisted selection.
The purpose of the present research was therefore to investigate the effect of some candidate genes on the variation of milk fat composition. Because the breed and the lactation stage are important sources of variability of milk yield and fat content (Karijord et al. Reference Karijord, Standal and Syrstad1982; Lawless et al. Reference Lawless, Stanton, Escop, Devery, Dillon and Murphy1999), three different breeds of cattle were included in the present study, and FA analyses were performed at different lactations stages for each cow, so as to evaluate the effect of the candidate genes after taking into account the major sources of variation of milk FA.
Materials and methods
Animals
The study included lactating cows of three breeds: Piedmontese and Valdostana are indigenous dairy breeds, and are still the most popular cattle of the lowlands of the North–West of Italy; the Jersey was introduced to Italy only a few decades ago, finding favour with farmers wanting to increase the milk's total solid content. The cows were reared in three nearby farms in the province of Turin (North–West of Italy), and were traditionally fed with ‘Unifeed’ (corn silage and concentrates) like all the dairy cows in that area. The farms were registered in the National Herdbook of each breed and performed the official milk recording of productivity following the regulation of the International Committee for Animal recording (ICAR). The study was conducted on 90 lactating cows (30 Jersey, 30 Piedmontese, 30 Valdostana) by collecting milk samples over three milk recordings of the same lactation, so as to have the first record before the 60th day from calving, the second record between the 100th and the 140th, and the third record on around the 210th day. Animals that calved between April and May, 2005, were chosen; they were of either second or third calving, in order to avoid seasonal and age effect on the parameters to be estimated. Analysis was carried out by using individual results of the milk recording (milk yield, fat and protein content) on the test day, the milk sample being collected in order to measure the FA. Animals of three different breeds were preferentially selected, particularly breeds which differed considerably in milk yield and composition, our hypothesis being that the effect of the candidate genes could be better estimated by taking into account, through the breed effect, the common genetic effects of all other genes, whose most common alleles are shared by each breed.
Analysis of milk
Milk samples (50 ml) were digested with 10 ml NH3 (25% v/v) and mixed with 40 ml ethanol (96% v/v). The extraction was performed with 100 ml of a mixture of diethyl ether and pentane (1 : 1 v/v). The solvent phase was filtered through 25 g anhydrous Na2SO4 and evaporated under vacuum. This procedure was based on the ISO method (2001). Milk fat (100 mg) was diluted in 5 ml hexane and derivatised as methylesters by adding 0·25 ml 2M-KOH in methanol (ISO, 2002). One millilitre of the upper phase containing fatty acid methyl esters (FAME) was diluted with 7 ml diethyl ether and 2 ml hexane for on-column injection. Gas chromatography analysis of FAME was performed by HP6890 (Agilent Technologies, Palo Alto, CA) and DB23 low bleed (J&W, Agilent Technologies, Palo Alto, CA) capillary column (30 m length, 0·32 mm I.D., 0·25 μm film thickness). On-column injection was adopted and hydrogen (1 ml/min) was used as a carrier gas. The temperature programme was as follows: 40 °C for 3 min, 25 °C/min up to 120 °C for 1 min, 4 °C/min up to 162 °C for 2 min, 8 °C/min up to 220 °C for 3 min; FID detector used, held at 250 °C. The characteristics of the capillary column did not allow the trans fatty acids to be separated. As a consequence the peak indicated as C18 : 1 includes trans-11, together with cis-9 and cis-7 isomers.
Statistical analysis
To investigate the variability of fatty acid methyl esters throughout the lactation, and in the three breeds, an ANOVA was performed with the following model (1):
$$Y_{ijkl} = \mu + B_i + B_i *b\left( {{\rm milk}} \right) + C_j + D_k + e_{ijkl\;} $$where: Y ijklm=Phenotype(1); μ=overall mean; B i=fixed effect of the breed; B i * b(milk)=covariate of the milk yield within breed; C j=fixed effect of the month of calving; D k=fixed effect of the days in milk (〈60; 100–140〉 141); b(milk)=covariate of daily milk yield within the breed; e ijk=residual.
(1)see Table 3.
Genes and genotyping
This study analysed three genes, currently considered, the world over, as strong positional and functional candidate genes in cattle: DGAT1 (Grisart et al. Reference Grisart, Farnir, Karim, Cambisano, Kim, Kvasz, Mni, Simon, Frère, Coppieters and Georges2004); SCD (Taniguchi et al. Reference Taniguchi, Utsugi, Oyama, Mannen, Kobayashi, Tanabe, Ogino and Tsuji2004; Stoop et al. Reference Stoop, Schennink, Visker, Mullaart, van Arendonk and Bovenhuis2009); and growth hormone receptor (GHR) (Blott et al. Reference Blott, Kim, Moisio, Schmidt-Kuntzel, Cornet, Berzi, Cambisano, Ford, Grisart, Johnson, Karim, Simon, Snell, Spelman, Wong, Vilkki, Georges, Farnir and Coppieters2003; Viitala et al. Reference Viitala, Szyda, Blott, Schulman, Lidauer, Maki-Tanila, Georges and Vilkki2006). Six more candidate genes were also selected for their potential effect on milk fat quality, and their involvement at various levels in fatty acid biosynthesis and desaturation as regulators, receptors, and encoding enzymes: diacylglycerol acyltransferase 2 (DGAT2), a candidate gene for quantitative traits related to dietary fat uptake, lipid synthesis and storage (Winter et al. Reference Winter, van Eckeveld, Bininda-Emonds, Habermann and Fries2003); fatty acid synthase (FASN), which is directly involved in FA synthesis and has been associated with milk fat content in cattle (Roy et al. Reference Roy, Ordovas, Zaragoza, Romero, Moreno, Altarriba and Rodellar2006; Bouwman et al. Reference Bouwman, Bovenhuis, Visker and van Arendok2011; Kaneda et al. Reference Kaneda, Linz, Sasazaki, Oyama and Mannen2011); acyl-CoA dehydrogenase (ACAD), which encodes a key enzyme in the β-oxidation of unsaturated fatty acids (Lea et al. Reference Lea, Abbas, Sprecher, Vockley and Schulz2000); fatty acid binding protein (FABP4) which falls into a significant QTL interval for beef marbling on bovine chromosome 14 (Michal et al. Reference Michal, Zhang, Gaskins and Jiang2006); lipoprotein lipase (LPL), considered a major candidate gene for lipid deposition in the yak (Ding et al. Reference Ding, Liang, Guo, Xing, Bao, Chu, Pei, Zhu and Yan2012); leptin gene (LEP), associated with production traits in dairy cattle (Buchanan et al. Reference Buchanan, Van Kessel, Waldner, Christensen, Laarveld and Schumutz2003; Madeja et al. Reference Madeja, Adamowicz, Chmurzynska, Jankowski, Melonek, Switonski and Strabel2004).
All animals of this trial were genotyped for the non-synonymous mutations of the DGAT1, SCD and LEP genes reported in Table 1; for the remaining genes, all animals were genotyped for the mutation previously described by Williams et al. (Reference Williams, Dunner, Valentini, Mazza, Amarger, Checa, Crisà, Razzaq, Delourme, Grandjean, Marchitelli, Garcia, Perez Gomez, Negrini, Ajmone Marsan and Leveziel2009), regardless of their location either in the coding sequence or within the intronic and untranslated regions (UTR).
Genes and SNP analysed in this study, location and AA change if present, and frequency of the minor allele in each breed (J=Jersey; P=Piedmontese; V=Valdostana)
A total of 21 SNP were then genotyped for this study at KBioscience Ltd. (Hoddesdon, UK) using their own novel fluorescence-based competitive allele-specific PCR (KASPar) assay. The locations of the SNP in each gene, their position relative to their GenBank accession number, and the AA changes, if present, are listed in Table 1.
Allele frequency and association analysis between alleles and milk traits
Allele frequency of each SNP was calculated with PROC FREQ in SAS (2007).
The allele substitution effect of each SNP was estimated, similarly to Barendse et al. (Reference Barendse, Harrison, Bunch and Thomas2008) and Sherman et al. (Reference Sherman, Nkrumah, Murdoch, Li, Wang, Fu and Moore2008), by regressing the number of copies of each allele against each of the traits separately, using the following Linear Mixed Model in SAS (2007) (2):
$$Y_{ijklm} = \mu + B_i + B_i *b\left( {{\rm milk}} \right) + C_j + B_i *G_k + A_m + e_{ijklmn} $$where: Y ijklm=Phenotype (as in model 1); μ=overall mean; B i=fixed effect of the breed; C j=fixed effect of the month of calving; B i*b(milk)=covariate of the milk yield within breed; G k=allele substitution effect, of the kth SNP within breed; A m=random animal effect; e ijklmn=residual.
A preliminary single locus association analysis was performed, in order to exclude, for each trait, the genes that had no influence (P<0·05), assuming that, if a marker had any effect on one trait, this would also become evident if the marker was considered on its own. A final association analysis was undertaken including, for each trait, all the possible allele combinations of the genes not excluded from the preliminary analysis. As in Manuguerra et al. (Reference Manuguerra, Matullo, Veglia, Autrup, Dunning, Garte, Gormally, Malaveille, Guarrera, Polidoro, Saletta, Peluso, Airoldi, Overvad, Raaschou-Nielsen, Clavel-Chapelon, Linseisen, Boeing, Trichopoulos, Kalandidi, Palli, Krogh, Tumino, Panico, Bueno De-Mesquita, Peeters, Lund, Pera, Martinez, Amiano, Barricarte, Tormo, Quiros, Berglund, Janzon, Jarvholm, Day, Allen, Saracci, Kaaks, Ferrari, Riboli and Vineis2007), all genes were considered in the same model, for two reasons: firstly, single locus effects cannot explain complex multifactorial traits; secondly, the sparseness of data in high dimensions, which often occurs when data are available for a large number of SNP for a relatively small number of samples, as in the present trial, could lead to false-positive effects (Mechanic et al. Reference Mechanic, Luke, Goodman, Chanock and Harris2008).
Results and discussion
Breed effect
The proportion of FA of different carbon chains was the same for all the breeds with the minimum percentage for the short chain (from 7·266% in the Valdostana to 8·737 in the Piedmontese), and the maximum percentage for the long chain (from 67·482% in the Piedmontese to 71·984% in the Valdostana) (Table 2). The breed effect was highly significant for most of the traits; in particular, saturated FA registered the highest level in the Jersey and the lowest in the Valdostana (Table 2); an opposite trend was noted for the monounsaturated FA (MUFA) (21·329% in the Jersey and 29·290% in the Valdostana). The single desaturation indices, given by the ratio of C14 : 1; C16 : 1 and C18 : 1 out of the corresponding saturated FA, confirmed this trend. No significant differences were noted between the Jersey and the Piedmontese for the FA from C6 : 0 to C12 : 0, that were significantly lower in the Valdostana. Polyunsaturated FA (PUFA) were highest in the Valdostana, due to the high content of C18 : 2 cis-9,trans-11 and C18 : 3 but not of the remaining C18 : 2.
Effect of the breed (LS Means) on each analysed variable and P of differences between breeds
Effect of the lactation stage
The effects of the lactation stage, reported in Table 3, showed that the days in milk affected the traits at various levels. Although all the breeds showed a similar behaviour in the FA variation during the lactation, with a decrease in the short-chain, accompanied by a general increase in the long chain FA, the variability between breeds in this respect was also evident. The decrease in milk yield was significant after the second milk record, in all the breeds, and particularly in the Jersey.
Fat content showed a significant increase throughout the lactation only in the Jersey and the Piedmontese. As for the short chain FA, C4 : 0 and C6 : 0 significantly decreased in the three breeds throughout the lactation; while the decrease of C8 : 0 was significant only in the Piedmontese and the Valdostana. Medium chain saturated FA up to 14 carbons significantly decreased through the lactation only in the Piedmontese. PUFA decreased through the lactation only in the Jersey, and C18 : 2 cis-9,trans-11 increased at late lactation only in the Valdostana.
Variation (LS Means) between lactation stages of each analysed variable and P of difference
The most evident differences between breeds were noted for the saturated FA. These FA increased throughout the lactation in the Jersey (from 72·7 to 76·8%; P<0·0001), while they decreased in the Piedmontese (from 72·8 to 68·7%; P=0·0002; see Table 3). The trend of the single saturated FA was however different. In fact, C16 : 0 decreased throughout the lactation in the Piedmontese, but registered a strong increase in the Jersey. On the other hand, C18 : 0 decreased in the Jersey but remained constant in the other two breeds. In the Valdostana no significant difference in the saturated FA throughout the lactation was registered.
The opposite behaviour of the saturated FA between the Jersey and the Piedmontese was accompanied, in both breeds, by the opposite behaviour of the MUFA, particularly due to C18 : 1, this decreasing from 20·8 to 17·4 (P<0·0001) in the Jersey while increasing from 21·3 to 24·6% (P=0·0002) in the Piedmontese (see Table 3).
Allele frequency
The five SNP of the FABP4, and the three SNP of the LPL showed a full within-breed linkage disequilibrium (LD; Table 1); therefore only one of them was used in the association analysis. The same thing occurred for two SNP of the GHR, and for two SNP of the DGAT2; the fact that one SNP of the GHR and one SNP of the DGAT2 were not fully correlated with the previous SNP might indicate that the genetic variability could be considered sufficient to regard the selected group of cows as a random sample of each breed under consideration, despite the relatively small sample size of each breed included in this trial.
Allele frequencies differed considerably between breeds (Table 1). Each breed sample was monomorphic at two SNP, but not at the same SNP as in the other breeds: the Jersey sample was monomorphic at the SNP of the promoter of the GHR, and at one of the SNP in the 3′UTR of the DGAT2. The Piedmontese sample was monomorphic at two consecutive SNP in the 3′UTR of the DGAT2. The Valdostana sample was monomorphic at one of the SNP in the 3′UTR of the DGAT2, and at the non-synonymous SNP of the DGAT1. Moreover, the SNP where MAF was ≦0·10, were not the same for the three breeds. These results might also indicate a sufficient genetic variability in the selected sample.
The four SNP of the DGAT2 are located in the 3′UTR at a relatively short distance (60 bp) between the two farthest of the four mutations. The first two SNP, two bp distant from each other, showed a full LD in all the breeds. The third was monomorphic in the Jersey and the Piedmontese samples, while the last SNP was monomorphic in the Piedmontese and the Valdostana samples. This unexpected behaviour of the DGAT2 alleles cannot be explained here, one reason being that the physiological role of the DGAT2 gene has not, to date, been investigated in depth: only Cases et al. (Reference Cases, Stone, Zhou, Yen, Tow, Lardizabal, Voelker and Farese2001) have reported that this gene has a function in tryglyceride storage, because of the high gene expression in adipose tissue. On the other hand, the functional causality of the DGAT1 K232A mutation has been confirmed (Grisart et al. Reference Grisart, Farnir, Karim, Cambisano, Kim, Kvasz, Mni, Simon, Frère, Coppieters and Georges2004).
The allele frequencies of the SNP under consideration, together with the smallness of the analysed sample, required a careful evaluation of the statistical model used to assess the potential influence of one SNP on each trait. After performing the preliminary SNP association analysis on each single trait, only the following SNP, which affected at least one trait, were included in the final model:
$$\eqalign{&\vskip 10pt {\rm FABP4 \,g} \,. \,{\rm 2989A} \gt {\rm G};{\rm FASN \ g}\, . \,{\rm 14726C} \gt {\rm A}; \cr & {\rm GHR \, g} \, .\,{\rm 4962A} \gt {\rm T};{\rm LPL \, g} \,. \,{\rm 25732G} \gt {\rm A};{\rm} \cr & {\rm DGAT2\, g}\,.\,{\rm 318}0{\rm 1T} \gt {\rm C};{\rm DGAT1\,g}\,.\,{\rm 1}0{\rm 434A} \gt {\rm G};\cr &{\rm SCD \,g}\,.\,{\rm 1}0{\rm 329T} \gt {\rm C}.} $$Allele substitution effect
Table 4 lists only the SNP showing a significant (P<0·05) effect on at least one trait, in more than one breed, and, obviously, if the effect of the same allele proved either positive or negative in all the breeds considered. In fact, a core assumption of the approach was that, if one gene affected one particular trait, this should be evident in different breeds, since this effect would overshadow the common genetic effects of all other genes whose alleles each breed has in common. Due to the relatively small sample size included in this trial, we also decided to take into account, in the discussion, the SNP that affected only two breeds, fearing that the possibly low frequency of one allele, in one breed, might prove insufficient to make the effect evident. For the same reason, the effects that tended to be significant (P<0·12) in one breed were also considered in the discussion, if these showed significant similar effects in the other breeds.
Genes and SNP for which an allele showed a significant effect on the indicated trait
As is evident from the table, the most important gene to affect milk fat composition was FABP4, mainly due to its influencing the medium chain FA in all the breeds; in particular, the allele which had a negative effect on medium chain FA (from 8 to 14 carbons), showed a positive effect on long chain FA. This suggests that this gene affected the partitioning of FA more than the desaturation. In fact, although the gene also influenced the ratio monounsaturated/saturated FA, this was because it affected the amount of the saturated C14 : 0, but not the homologous unsaturated FA. These results seemed particularly interesting, since, in the present study, a wide portion of the FABP4 gene has been considered (Table 1), this including three different exons and the 3′UTR. The four different SNP of this gene, which were genotyped, showed a full LD within breed, but not between breeds: this suggests the presence of a selection signature. Moreover, the SNP in exon 2 encoded a non-synonymous mutation (Ile/Val), located at the splice site, thus corroborating the finding of a new candidate gene for milk FA composition.
The effect of the SCD Val/Ala SNP on milk fat desaturation index was evident only for the C14 index (Table 4), although Stoop et al. (2009), in the Dutch Holstein Friesian population, reported that the SCD Val allele is associated with lower C10, C12, and C14 indices, and with higher C16, C18, and CLA indices in comparison with the SCD Ala allele. Because long-chain FA originate mainly from diet, it is possible that, across breeds, the effect of the SCD genotype can be perceived only on the short-chain FA, these FA being synthesised de novo in the mammary gland, and less dependent on the diet. Moioli et al. (Reference Moioli, Contarini, Avalli, Catillo, Orru, De Matteis, Masoero and Napolitano2007) had also reported that this SNP affected only the short-chain FA (C10 and C14 indices); on the other hand, in beef of the Japanese Black breed, Matsuhashi et al. (Reference Matsuhashi, Maruyama, Uemoto, Kobayashi, Mannen, Abe, Sakaguchi and Kobayashi2011) reported that this SNP affected only the C14 and C18 indices.
FASN affected PUFA, but not in the Jersey breed. Although no QTL effect has to date been reported for PUFA in cattle, Garcia-Fernandez et al. (Reference Garcia-Fernandez, Gutierrez-Gil, Garcia-Gamez, Sanchez and Arranz2010) suggested that the peaks of the QTL affecting PUFA in sheep milk map close to the FASN gene. These results are corroborated by the fact that the FASN gene encodes a key enzyme in the de novo synthesis of FA in the mammary gland, being responsible for the chain elongation of FA.
As regards the DGAT1 gene, only the Jersey breed showed the polymorphism at the K232A mutation in exon 7, with a significant effect of the g.10434A allele on fat content (+0·350; P=0·02) (not reported in Table 4), which confirms what was detected in the Holstein cows (Schennink et al. Reference Schennink, Stoop, Visker, Heck, Bovenhuis, van der Poel, van Valenberg and van Arendonk2007). By contrast, no effect was detected on any FA, although Schennink et al. (Reference Schennink, Heck, Bovenhuis, Visker, van Valenberg and van Arendonk2008) reported that this polymorphism influenced milk desaturation indices, the effect having been detected in a population of Dutch Holstein cows. Finally, although not reported in Table 4, the g.25732A allele of LPL significantly affected fat percentage in the Piedmontese (+0·312; P=0·02) to an extent comparable with the effect of the DGAT1 SNP in the Jersey. No direct effect, of any of the considered genes, on milk fat content, was observed in the Valdostana.
It should be noted that only the mutations and the genes, as described by Williams et al. (Reference Williams, Dunner, Valentini, Mazza, Amarger, Checa, Crisà, Razzaq, Delourme, Grandjean, Marchitelli, Garcia, Perez Gomez, Negrini, Ajmone Marsan and Leveziel2009), were taken into account in the present study. It is therefore possible that other SNP in the same genes have an influence on FA composition; what is most interesting, however, is that a further gene, the FABP4, might be added to the list of the known candidate genes for milk quality.
Although the animal sample was small, this work has presented some new features unmentioned by previous studies aimed at evaluating the candidate gene effects on the variability of milk FA. On the one hand, a greater number of candidate genes have here been analysed; this is contradistinction to previous studies which limited themselves to the DGAT1 and the SCD genes. On the other hand, three breeds, with markedly different milk yields, have been considered: this has enabled us to take into account, through the breed effect, the alleles of all other genes common to each breed. Finally, some of our results have shed new light on the variability of FA composition throughout the lactation, demonstrating that this variability was not the same for all the breeds.
This study is part of the ‘GENZOOT’ research programme, funded by the Italian Ministry of Agriculture.
