From milk to cheese: genomic background, biological pathways and latent phenotypes of bovine cheese-related traits

The aim of this PhD thesis was the study of the genomic, biological and phenotypic background of bovine cheese-related traits. The primary goal of this PhD thesis was to unravel the genomic background of bovine milk technological and cheese-related traits to specific chromosomic regions (CHAPTERS 1 to 3). To achieve this, the cow’s ability to produce cheese was decomposed into 11 milk coagulation (MCP) and curd-firming properties (CFt), and 7 cheese yield and milk component recoveries into the curd (REC) traits. Besides, to tackle the problem of the large number of variables required to describe the cow’s ability to produce cheese, posing restrictions in the construction of selection indices, and thereby selection decisions, factor analysis (FA) was used (CHAPTERS 4 and 5). The MCP traits were: 3 traditional single point lacto-dynamographic properties (RCT: rennet coagulation time, min; k20: time to a curd firmness (CF) of 20 mm, min; a30: CF 30 min after rennet addition), 6 parameters modeling 360 CF data for each milk sample (CFP: potential asymptotic CF at infinite time, mm; kCF: curd firming instant rate constant, %×min-1; kSR: syneresis instant rate constant, %×min-1; RCTeq: RCT from modeling; CFmax: maximum CF, mm; tmax: time at CFmax, min), milk- protein (%) and pH. The 3 CY traits were the weight (wt) of fresh curd (%CYCURD), curd solids (%CYSOLIDS), and curd moisture (%CYWATER) as % of wt of milk processed. The 4 REC (RECFAT, RECPROTEIN, RECSOLIDS, and RECENERGY) were calculated as the % ratio between the nutrient in curd and the corresponding nutrient in processed milk. For FA 26 traits related to milk yield and quality (including milk protein fractions), MCP-CFt and CY-REC traits were analyzed. Single marker genome-wide association analyses (GWAS) complemented by gene-set enrichment and pathway-based analyses were conducted. In total, 1,152 Italian Brown Swiss cows reared in 85 herds were genotyped with the Illumina SNP50 Beadchip v.2. Single marker regression GWAS were fitted using the GenABEL R package (GRAMMAR-GC). The Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway databases were queried for the enrichment analyses. In GWAS (CHAPTERS 1 and 2), sharp peaks were detected on Bos taurus autosome (BTA) 6, at 84 to 88 Mbp, with the highest peak detected at 87.4 Mbp in the region harboring the casein genes and more precisely of κ-CN (CSN3). Marker Hapmap52348-rs29024684 (~87.4 Mbp), closely located to the casein genes on BTA6, was strongly associated with RECFAT (P = 1.91×10-15) and CFP (P = 1.62×10-17). Evidence of quantitative trait loci at 82.6 and 88.4 Mbp on the same chromosome was found. On BTA11, marker ARS-BFGL-NGS-104610 (~104.3 Mbp) was highly associated with RECPROTEIN (P = 6.07×10-36). Apart from BTA6 and 11, SNP located in 15 more chromosomes (1, 2, 9, 12, 13, 14, 15, 16, 18, 19, 20, 23, 26, 27 and 28) were significantly associated to the MCP-CFt and CY-REC traits. The gene-set enrichment and pathway-based analysis (CHAPTER 3) revealed 21 GO and 17 KEGG categories significantly associated (false discovery rate controlled at 5%) with 7 of the traits (RCT, RCTeq, kCF, %CYSOLIDS, RECFAT, RECSOLIDS and RECENERGY), with some being in common between traits. The significantly enriched categories included calcium signaling pathway, salivary secretion, metabolic pathways, carbohydrate digestion and absorption, the tight junction and the phosphatidylinositol pathways, as well as pathways related to the bovine mammary gland health status, and contained a total of 150 genes located in all chromosomes but 9, 20, and 27. In FA (CHAPTERS 4 and 5), ten mutual orthogonal Fs were obtained using a varimax rotation. The 10 Fs explaining 74% of the original variability. Those Fs captured basic concepts of the “milk to cheese” process. More precisely, the first four Fs, sorted by variance explained, were able to capture the underlying structure of the CY percentage (F1%CY), the CF process with time (F2CFt), the milk and solids yield (F3Yield) and the presence of nitrogen (N) into the cheese (F4Cheese N). Moreover, 4 Fs (F5 αs1-β-CN, F7β-κ-CN, F8αs2-CN, F9αs1-CN-P) were related to the basic milk caseins (as1-CN, as2-CN, β-CN, κ-CN, and the phosphorylated form of as1-CN) and 1 factor was associated with the α-LA whey protein (F10α-LA). A factor describing the udder health status of a cow (F6Udder health), mainly loaded on lactose, other nitrogen compounds and SCS, was also obtained. In general, FA results were coherent to the given name of the factor. Stage of lactation had a significant effect for the majority of the Fs, followed by parity. Moreover, considerable genetic correlations existed among the Fs (CHAPTER 4). All Fs showed significant associations (P < 5 ×10-5) in GWAS, but F5Yield. High peaks on BTA6 (~87Mbp) and at the tail of BTA11 (~104Mbp) were mainly associated to F6β-κ-CN and F1Cheese N, respectively. In addition, 33 GO terms and 6 KEGG categories were mainly enriched for F8αs2-CN, but also for F1%CY, F4Cheese N, and F10α-LA. Biological pathways were mainly related to the broader categories of ion activity, neurons and the tight junction. Moreover, the considerably large number of enriched GO and KEGG terms for F8αs2-CN suggests that, perhaps, more focus should be given on αs2-CN (CHAPTER 5)

Keep Garfagnina alive. An integrated study on patterns of homozygosity, genomic inbreeding, admixture and breed traceability of the Italian Garfagnina goat breed

The objective of this study was to investigate the genetic diversity of the Garfagnina (GRF) goat, a breed that currently risks extinction. For this purpose, 48 goats were genotyped with the Illumina CaprineSNP50 BeadChip and analyzed together with 214 goats belonging to 9 other Italian breeds (~25 goats/breed), whose genotypes were available from the AdaptMap project [Argentata (ARG), Bionda dell'Adamello (BIO), Ciociara Grigia (CCG), Di Teramo (DIT), Garganica (GAR), Girgentana (GGT), Orobica (ORO), Valdostana (VAL) and Valpassiria (VSS)]. Comparative analyses were conducted on i) runs of homozygosity (ROH), ii) admixture ancestries and iii) the accuracy of breed traceability via discriminant analysis on principal components (DAPC) based on cross-validation. ROH analyses was used to assess the genetic diversity of GRF, while admixture and DAPC to evaluate its relationship to the other breeds. For GRF, common ROH (more than 45% in GRF samples) was detected on CHR 12 at, roughly 50.25-50.94Mbp (ARS1 assembly), which spans the CENPJ (centromere protein) and IL17D (interleukin 17D) genes. The same area of common ROH was also present in DIT, while a broader region (~49.25-51.94Mbp) was shared among the ARG, CCG, and GGT. Admixture analysis revealed a small region of common ancestry from GRF shared by BIO, VSS, ARG and CCG breeds. The DAPC model yielded 100% assignment success for GRF. Overall, our results support the identification of GRF as a distinct native Italian goat breed. This work can contribute to planning conservation programmes to save GRF from extinction and will improve the understanding of the socio-agro-economic factors related with the farming of GRF

Predictive formulas for different measures of cheese yield using milk composition from individual goat samples

The objective of this study was to develop formulas based on milk composition of individual goat samples for predicting cheese yield (%CY) traits (fresh curd, milk solids, and water retained in the curd). The specific aims were to assess and quantify (1) the contribution of major milk components (fat, protein, and casein) and udder health indicators (lactose, somatic cell count, pH, and bacterial count) on %CY traits (fresh curd, milk solids, and water retained in the curd); (2) the cheese-making method; and (3) goat breed effects on prediction accuracy of the %CY formulas. The %CY traits were analyzed in duplicate from 600 goats, using an individual laboratory cheese-making procedure (9-MilCA method; 9 mL of milk per observation) for a total of 1,200 observations. Goats were reared in 36 herds and belonged to 6 breeds (Saanen, Murciano-Granadina, Camosciata delle Alpi, Maltese, Sarda, and Sarda Primitiva). Fresh %CY (%CYCURD), total solids (%CYSOLIDS), and water retained (%CYWATER) in the curd were used as response variables. Single and multiple linear regression models were tested via different combinations of standard milk components (fat, protein, casein) and indirect udder health indicators (UHI; lactose, somatic cell count, pH, and bacterial count). The 2 %CY observations within animal were averaged, and a cross-validation (CrV) scheme was adopted, in which 80% of observations were randomly assigned to the calibration (CAL) set and 20% to the validation (VAL) set. The procedure was repeated 10 times to account for sampling variability. Further, the model presenting the best prediction accuracy in CrV (i.e., comprehensive formula) was used in a secondary analysis to assess the accuracy of the %CY predictive formulas as part of the laboratory cheese-making procedure (within-animal validation, WAV), in which the first %CY observation within animal was assigned to CAL, and the second to the VAL set. Finally, a stratified CrV (SCrV) was adopted to assess the %CY traits prediction accuracy across goat breeds, again using the best model, in which 5 breeds were included in CAL and the remaining one in the VAL set. Fitting statistics of the formulas were assessed by coefficient of determination of validation (R2VAL) and the root mean square error of validation (RMSEVAL). In CrV, the formula with the best prediction accuracy for all %CY traits included fat, casein, and UHI (R2VAL = 0.65, 0.96, and 0.23 for %CYCURD, %CYSOLIDS, and %CYWATER, respectively). The WAV procedure showed R2VAL higher than those obtained in CrV, evidencing a low effect of the 9-MilCA method and, indirectly, its high repeatability. In the SCrV, large differences for %CYCURD and %CYWATER among breeds evidenced that the breed is a fundamental factor to consider in %CY predictive formulas. These results may be useful to monitor milk composition and quantify the influence of milk traits in the composite selection indices of specific breeds, and for the direct genetic improvement of cheese production

Population Structure and Genetic Diversity of Italian Beef Breeds as a Tool for Planning Conservation and Selection Strategies

The aim was to investigate the population structure of eight beef breeds: three local Tuscan breeds under extinction, Calvana (CAL), Mucca Pisana (MUP), and Pontremolese (PON); three local unselected breeds reared in Sardinia, Sarda (SAR), Sardo Bruna (SAB), and Sardo Modicana (SAM); and two cosmopolitan breeds, Charolais (CHA) and Limousine (LIM), reared in the same regions. An effective population size ranges between 14.62 (PON) to 39.79 (SAM) in local breeds, 90.29 for CHA, and 135.65 for LIM. The average inbreeding coefficients were higher in Tuscan breeds (7.25%, 5.10%, and 3.64% for MUP, CAL, and PON, respectively) compared to the Sardinian breeds (1.23%, 1.66%, and 1.90% in SAB, SAM, and SAR, respectively), while for CHA and LIM they were &lt;1%. The highest rates of mating between half-siblings were observed for CAL and MUP (~9% and 6.5%, respectively), while the highest rate of parent&ndash;offspring mating was ~8% for MUP. Our findings describe the urgent situation of the three Tuscan breeds and support the application of conservation measures and/or the development of breeding programs. Development of breeding strategies is suggested for the Sardinian breeds

Prediction accuracies of cheese-making traits using Fourier-transform infrared spectra in goat milk

The objectives of this study were to explore the use of Fourier-transform infrared (FITR) spectroscopy on 458 goat milk samples for predicting cheese-making traits, and to test the effect of the farm variability on their prediction accuracy. Calibration equations were developed using a Bayesian approach with three different scenarios: i) a random cross-validation (CV) [80% calibration (CAL); 20% validation (VAL) set], ii) a stratified CV [(SCV), 13 farms used as CAL, and the remaining one as VAL set], and iii) a SCV where 20% of the goats randomly selected from the VAL farm were included in the CAL set (SCV80). The best prediction performance was obtained for cheese yield solids, justifying for its practical application at population level. Overall results were similar to or outperformed those reported for bovine milk. Our results suggest considering specific procedures for calibration development to propose reliable tools applicable along the dairy goat chain

Genetic Diversity in the Italian Holstein Dairy Cattle Based on Pedigree and SNP Data Prior and After Genomic Selection

Genetic diversity has become an urgent matter not only in small local breeds but also in more specialized ones. While the use of genomic data in livestock breeding programs increased genetic gain, there is increasing evidence that this benefit may be counterbalanced by the potential loss of genetic variability. Thus, in this study, we aimed to investigate the genetic diversity in the Italian Holstein dairy cattle using pedigree and genomic data from cows born between 2002 and 2020. We estimated variation in inbreeding, effective population size, and generation interval and compared those aspects prior to and after the introduction of genomic selection in the breed. The dataset contained 84,443 single-nucleotide polymorphisms (SNPs), and 74,485 cows were analyzed. Pedigree depth based on complete generation equivalent was equal to 10.67. A run of homozygosity (ROH) analysis was adopted to estimate SNP-based inbreeding (FROH). The average pedigree inbreeding was 0.07, while the average FROH was more than double, being equal to 0.17. The pattern of the effective population size based on pedigree and SNP data was similar although different in scale, with a constant decrease within the last five generations. The overall inbreeding rate (ΔF) per year was equal to +0.27% and +0.44% for Fped and FROH throughout the studied period, which corresponded to about +1.35% and +2.2% per generation, respectively. A significant increase in the ΔF was found since the introduction of genomic selection in the breed. This study in the Italian Holstein dairy cattle showed the importance of controlling the loss of genetic diversity to ensure the long-term sustainability of this breed, as well as to guarantee future market demands

A genome-wide association analysis for body weight at 35 days measured on 137,343 broiler chickens

Background: Body weight (BW) is an economically important trait in the broiler (meat-type chickens) industry. Under the assumption of polygenicity, a "large" number of genes with "small" effects is expected to control BW. To detect such effects, a large sample size is required in genome-wide association studies (GWAS). Our objective was to conduct a GWAS for BW measured at 35 days of age with a large sample size.Methods: The GWAS included 137,343 broilers spanning 15 pedigree generations and 392,295 imputed single nucleotide polymorphisms (SNPs). A false discovery rate of 1% was adopted to account for multiple testing when declaring significant SNPs. A Bayesian ridge regression model was implemented, using AlphaBayes, to estimate the contribution to the total genetic variance of each region harbouring significant SNPs (1 Mb up/downstream) and the combined regions harbouring non-significant SNPs.Results: GWAS revealed 25 genomic regions harbouring 96 significant SNPs on 13 Gallus gallus autosomes (GGA1 to 4, 8, 10 to 15, 19 and 27), with the strongest associations on GGA4 at 65.67-66.31 Mb (Galgal4 assembly). The association of these regions points to several strong candidate genes including: (i) growth factors (GGA1, 4, 8, 13 and 14); (ii) leptin receptor overlapping transcript (LEPROT)/leptin receptor (LEPR) locus (GGA8), and the STAT3/STAT5B locus (GGA27), in connection with the JAK/STAT signalling pathway; (iii) T-box gene (TBX3/TBX5) on GGA15 and CHST11 (GGA1), which are both related to heart/skeleton development); and (iv) PLAG1 (GGA2). Combined together, these 25 genomic regions explained similar to 30% of the total genetic variance. The region harbouring significant SNPs that explained the largest portion of the total genetic variance (4.37%) was on GGA4 (similar to 65.67-66.31 Mb).Conclusions: To the best of our knowledge, this is the largest GWAS that has been conducted for BW in chicken to date. In spite of the identified regions, which showed a strong association with BW, the high proportion of genetic variance attributed to regions harbouring non-significant SNPs supports the hypothesis that the genetic architecture of BW35 is polygenic and complex. Our results also suggest that a large sample size will be required for future GWAS of BW35

Goat farm variability affects milk Fourier-transform infrared spectra used for predicting coagulation properties

Peer-ReviewedDriven by the large amount of goat milk destined for cheese production, and to pioneer the goat cheese industry, the objective of this study was to assess the effect of farm in predicting goat milk-coagulation and curd-firmness traits via Fourier-transform infrared spectroscopy. Spectra from 452 Sarda goats belonging to 14 farms in central and southeast Sardinia (Italy) were collected. A Bayesian linear regression model was used, estimating all spectral wavelengths' effects simultaneously. Three traditional milk-coagulation properties [rennet coagulation time (min), time to curd firmness of 20 mm (min), and curd firmness 30 min after rennet addition (mm)] and 3 curd-firmness measures modeled over time [rennet coagulation time estimated according to curd firmness change over time (RCTeq), instant curd-firming rate constant, and asymptotical curd firmness] were considered. A stratified cross validation (SCV) was assigned, evaluating each farm separately (validation set; VAL) and keeping the remaining farms to train (calibration set) the statistical model. Moreover, a SCV, where 20% of the goats randomly taken (10 replicates per farm) from the VAL farm entered the calibration set, was also considered (SCV80). To assess model performance, coefficient of determination (R2VAL) and the root mean squared error of validation were recorded. The R2VAL varied between 0.14 and 0.45 (instant curd-firming rate constant and RCTeq, respectively), albeit the standard deviation was approximating half of the mean for all the traits. Although average results of the 2 SCV procedures were similar, in SCV80, the maximum R2VAL increased at about 15% across traits, with the highest observed for time to curd firmness of 20 mm (20%) and the lowest for RCTeq (6%). Further investigation evidenced important variability among farms, with R2VAL for some of them being close to 0. Our work outlined the importance of considering the effect of farm when developing Fourier-transform infrared spectroscopy prediction equations for coagulation and curd-firmness traits in goats.Università degli Studi di Sassar

Tracing Cattle Breeds with Principal Components Analysis Ancestry Informative SNPs

The recent release of the Bovine HapMap dataset represents the most detailed survey of bovine genetic diversity to date, providing an important resource for the design and development of livestock production. We studied this dataset, comprising more than 30,000 Single Nucleotide Polymorphisms (SNPs) for 19 breeds (13 taurine, three zebu, and three hybrid breeds), seeking to identify small panels of genetic markers that can be used to trace the breed of unknown cattle samples. Taking advantage of the power of Principal Components Analysis and algorithms that we have recently described for the selection of Ancestry Informative Markers from genomewide datasets, we present a decision-tree which can be used to accurately infer the origin of individual cattle. In doing so, we present a thorough examination of population genetic structure in modern bovine breeds. Performing extensive cross-validation experiments, we demonstrate that 250-500 carefully selected SNPs suffice in order to achieve close to 100% prediction accuracy of individual ancestry, when this particular set of 19 breeds is considered. Our methods, coupled with the dense genotypic data that is becoming increasingly available, have the potential to become a valuable tool and have considerable impact in worldwide livestock production. They can be used to inform the design of studies of the genetic basis of economically important traits in cattle, as well as breeding programs and efforts to conserve biodiversity. Furthermore, the SNPs that we have identified can provide a reliable solution for the traceability of breed-specific branded products