Associations between Milk Coagulation Properties and Microbiological Quality in Sheep Bulk Tank Milk

This study conducted a seasonal analysis of bulk tank milk from 77 sheep farms to establish relationships between the concentration of major microbial groups and milk coagulation properties. The investigated milk traits included composition (pH, fat, casein, lactose), coagulation properties (curd firmness: A60-, rennet clotting time: RCT-, curd firming time: k20-, curd yield: CY-), and somatic cell score (SCS). The main microbial groups analyzed were total mesophilic bacteria (SPC), thermodurics (THERMO), psychrotrophs (PSYCHRO), Pseudomonas spp. (PSEUDO), lactic acid bacteria (LAB), catalase-negative gram-positive cocci (GPCNC), Escherichia coli (ECOLI), coliforms other than Escherichia coli (COLI), coagulase-positive staphylococci (CPS), coagulase-negative staphylococci (CNS), and spores of lactate-fermenting Clostridium (BAB). Mixed linear models were used to explore associations between coagulation properties and the aforementioned variables. Results demonstrated that incorporating microbial loads into the models improves their fit and the relative quality of the outcomes. An important seasonality is demonstrated by an increase in CY and A60, along with a decrease in RCT and k20 during autumn and winter, contrasting with spring and summer. BAB concentration resulted in a reduction of A60 and an increase in RCT, whereas SPC concentration led to an enhancement of A60 and a reduction in RCT. An increase in GPCNC concentration was associated with an increase in k20 and a decrease in CY.


Introduction
Sheep milk production holds significant importance in countries across the Mediterranean basin.Traditionally, native breeds have been raised, with many now facing the threat of extinction.These breeds have played a crucial role in sustaining livestock activity and associated industries in vast rural areas.This not only helps to maintain local populations but also contributes to environmental protection [1].Sheep milk production is primarily directed towards cheese and dairy product processing, with many of these products holding protected designation of origin (PDO) status [2].Many of these cheeses are known for their unique sensory characteristics, which result from traditional manufacturing techniques.The quality of sheep milk is a key factor in creating exceptional cheeses and is subject to rigorous quality control systems in the industry [3].Traditionally, milk quality assessment has focused on physicochemical aspects, mainly fat and protein content.Recently, hygienic parameters have been integrated into these quality systems, including analyses like total mesophilic microorganism and somatic cell counts (SCC).However, practices considering aspects such as differential milk microbiology are still limited [3].
Research in dairy cattle has provided valuable insights into the microbial composition of milk and its relationship with several sources of contamination, emphasizing the importance of this aspect in milk quality [4].These contaminants introduce a variety of microorganisms associated with farm environments, such as lactic acid bacteria, hygiene practices on farms leading to the presence of coliforms, and intramammary infections caused by staphylococci [5].Other researchers [6,7] have also highlighted the complex microbial composition found in milk from different species, emphasizing the importance of its contamination sources and its impact on cheese-making processes, especially concerning lactic acid bacteria.Some studies have even suggested that stringent hygiene measures have significantly reduced microbial populations, potentially adversely affecting the production of traditional and artisanal cheeses [8].However, there is limited research specifically addressing the differential microbial load and microbiota of sheep milk [9,10].
Several studies have underscored the significance of milk coagulation properties (MCP) as complementary traits alongside the inherent quality of raw milk [11,12].For that reason, there is well-documented literature on the relevance and significance of MCP, mainly in cattle [13,14] and, to a lesser extent, in sheep and goats [15,16].Also, these studies underline the importance of improving milk MCP and explore the relationship with several factors that contribute its modification, such as breeding, feeding, farm management system, season, or milking system.Factors such as bacterial counts and the intricate interplay among microbiological characteristics and technological parameters can influence MCP and cheese production.However, research in this area, particularly regarding its impact on milk from small ruminants, remains limited [5].
Therefore, this study aims to assess the impact of composition and microorganisms present in sheep milk on its coagulation properties.Additionally, it seeks to compare the predictive performance of models incorporating microbiological variables with reference models that do not include these variables.The ultimate goal is to reach a better understanding of the microbiological factors influencing milk coagulation, providing valuable information for quality control and dairy production management.

Study Design and Sample Determination
The present study included 77 Manchega dairy sheep flocks from the region of Castilla-La Mancha (Spain), registered under Protected Designation of Origin "Queso Manchego", representing 15% of the farms listed in PDO [17].The average flock size was 1144 ewes, ranging from 150 to 5500 ewes, with an average production of 143 kg/ewe/year, varying from 91 to 200 kg, over a mean lactation period of 121 days, ranging from 100 to 150 days.All farms employ mechanical milking systems, and ewes are milked twice daily.The milk is then filtered and stored in refrigerated tanks, maintaining a maximum temperature of 6 • C. The production system is semi-intensive, traditionally associated with grazing on natural pastures and utilizing the residues and remains of cereal crops.On average, each sheep has access to 1.44 ha, ranging from 0.5 to 1.94 ha.
The flocks were visited four times a year (once during each season) to collect bulk tank samples, which included milk from the morning milking of the visit day and the previous evening milking.Samples were stored in hermetically sealed sterilized containers and kept at 4 • C until the arrival at the Experimental Dairy Science and Cheese Laboratory at the Regional Centre for Animal Selection and Reproduction (CERSYRA-IRIAF, Valdepeñas, Ciudad Real, Spain), where they were aliquoted and sent to other facilities for further analyses.

Laboratory Analysis
Microbiological analyses were conducted at CERSYRA-IRIAF within the first 24 h following sample collection.The techniques and methods used to determine the micro-Foods 2024, 13, 886 3 of 13 biological counts are described in detail in a previous paper from our research group [5].Serial dilutions were performed from each sample to inoculate 0.1 mL onto different culture media.Bacterial counts for the following groups of microorganisms were determined on PCA media (Panreac, Barcelona, Barcelon, Spain): standard plate count (SPC) was incubated in aerobic conditions at 30 • C for 72 h [18]; thermoduric bacteria (THERMO, Panreac, Barcelona, Spain) were incubated under the same conditions as SPC, after pasteurizing milk at 62.8 • C for 30 min [19]; and psychrotrophic bacteria (PSYCHRO, Panreac, Barcelona, Spain) were incubated at 6.5 • C for 10 days [20].Pseudomonas spp.(PSEUDO, Panreac, Barcelona, Spain) was cultured on Cetrimide agar (Panreac, Spain) and incubated at 35 • C for 48 h [21].Determination of Escherichia coli (ECOLI, Basingstoke Hampshire, UK) and other coliforms (COLI, Basingstoke Hampshire, UK) was conducted using CromoIDTM Coli medium (bioMérieux, Madrid, Spain), incubated at 37 • C for 24 h [22].Gram-positive catalase-negative cocci count (GPCNC) was determined in modified Edwards medium with colistin and oxolinic acid supplement (Oxoid, Basingstoke Hampshire, UK), incubated at 35 • C for 48 h [9].Lactic acid bacteria (LAB) were seeded on MRS medium (Panreac, Barcelona, Spain) acidified to pH 5.7, and incubated at 30 • C for 72 h [23].Lactatefermenting Clostridium spore (BAB, Basingstoke Hampshire, UK) count was performed using the most probable number (MPN) technique in Bryant and Burkey Broth (BBB, Merck, Darmstadt, Germany) [24].For the enumeration of coagulase-positive staphylococci (CPS) and coagulase-negative staphylococci (CNS), Baird-Parker RPF agar medium (bioMérieux, Madrid, Spain) was used, and incubated at 37 • C for 24 h [4].The microbial counts were expressed in colony-forming units per mL of milk (CFU/mL), and Clostridium spores were expressed as spores/mL.A logarithmic transformation (log 10 ) was applied to both microbial counts to normalize their distribution.
Milk composition, pH, and somatic cell count (SCC) analyses were conducted at the Interprofessional Dairy Laboratory of Castilla-La Mancha (LILCAM, Talavera de la Reina, Toledo, Spain) within 48 h from sample collection.The pH was measured using a Crison Basic 20 pH meter (Crison Instruments, Barcelona, Spain).Milk composition, including fat content (FAT), total casein (CAS), and lactose (LAC) was determined using a Milkoscan 6000 FT device (Foss Electric, Hillerød, Denmark), which operates on mid-infrared spectroscopy.The total SCC was quantified by flow cytometry on a Fossomatic FC (Foss Electric, Hillerød, Denmark).The SCC was subsequently expressed as somatic cell scores (SCS) by applying the following formula [25]: SCS = log 2 [SCC cells/mL × 10 −5 ] + 3 Milk coagulation properties were monitored at the Small Ruminant Dairy Laboratory (University of Córdoba, Spain) within the first 24 h following sample collection.Using a Formagraph viscosimeter (Foss Electric A/S, Hillerød, Denmark), key coagulation parameters were determined based on the methodology adapted from McMahon and Brown [26].These parameters included rennet coagulation time (RCT), curd-firming time (k 20 ), curd firmness at 60 min after rennet addition (A 60 ), and curd yield (CY). Curd yield was subsequently calculated following draining of the fresh curds [5].

Statistical Analysis
Data analysis was performed on SAS version 9.3 (SAS Institute Inc., Cary, NC, USA).Preliminary testing of the data was carried out to determine outliers to be discarded before further analysis.After the editing, a total of 306 data points were used (Table 1).The explanatory variables (i.e., SCS, BAB, SPC, THERMO, PSYCHRO, PSEUDO, LAB, GPCNC, ECOLI, COLI, CPS, CNS) were discretized into 5 classes based on the 20th, 40th, 60th, and 80th percentiles.This approach was taken to ensure a better assessment of the pattern of microbiological count effects and to guarantee the classes were balanced [27].A linear mixed model (MIXED) procedure was used to investigate the associations between microbiological quality and milk coagulation properties.The model is represented as follows: where Y denotes milk coagulation traits (CY, A 60 , RCT, k 20 ), β i are the unknown parameters to be estimated, f i denotes the explanatory variables, and δ is the error term.
The WSM model incorporating microbiological counts was adjusted using the same variables as BM, completing it with microbiological counts discretized into five classes as fixed factors (Table 2).
To address challenges arising from the limited sample size (308 complete cases), a step-by-step sequence was undertaken, considering the primary effects of all fixed, random, and covariate factors.A significance level of p < 0.05 was used for the retention criterion.In the initial phase, all conceivable single-variable models were assessed using the Akaike Information Criterion (AIC) value.Subsequently, each remaining factor was incrementally added, one by one, to the model with the lowest AIC value, and each predictor was compared based on the AIC value.This iterative process continued until the model with the lowest AIC was identified [28].This model was deemed the most plausible and was selected as the final model.After defining the model, an analysis of collinearity was conducted using the variance inflation factor (VIF) of the regression coefficients.A severe multicollinearity issue was present if any of the VIFs exceeded 10 [28].The normal distribution of residuals was verified using Kolmogorov-Smirnov test, while the absence of autocorrelation of residuals was assessed through the Durbin-Watson test.Heteroscedasticity was evaluated using White's test.The model fit was assessed using the adjusted coefficient of determination (R 2 ) and the root mean squared error (RMSE) [29].

Results and Discussion
The most plausible models identified for each of the analyzed coagulation characteristics are shown in Table 3. Broadly, the models exhibited poor fit (R 2 < 50%), except for CY.In the case of curd yield, the coefficient of determination ranged between 54.8% for the baseline model and 57.9% for the model incorporating microbiological counts.This emphasizes the fact that the inclusion of microbiological counts not only significantly improved the model fit but also enhanced relative quality across all evaluated outcomes.A significant association was observed between season and all coagulation traits.This association persisted even in models that incorporated microbiological counts, albeit with a slightly reduced contribution.As expected, season was a significant source of variation for all investigated traits, underscoring the importance of considering this effect in the statistical model for a better evaluation of the effect of microbiological counts.It is important to note that Manchega does not exhibit as pronounced a seasonality in milk production as other Mediterranean breeds, such as Sarda [30] or Latxa [31].Therefore, the seasonal differences found are arguably more interesting in terms of their relevance to the cheese industry.
The impact of flock on the variability of the examined coagulation traits was modest, as indicated by previous research [11,32].Although the AIC favored models excluding the flock effect, this does not completely rule out its influence.These results are explained by the fact that the evaluated flocks follow a similar farming system, with other variables within the scope of this study exerting a more pronounced influence on coagulation characteristics than flock variability.The level of somatic cells was only significantly associated with RCT.This result is in line with a previous study by Pazzola et al. [33] focused on the Sarda breed, in which higher somatic cell counts were linked to slower coagulations.Moreover, Caballero-Villalobos et al. [34] observed in individual samples from Manchega that an increase in SCC led to a rise in RCT and k 20 and a decrease in curd firmness measured at 30 min.This could be attributed to an increase in casein breakdown as a consequence of increased plasmin activity, the main enzyme associated with somatic cells in milk [35].
Regarding CY, the most plausible baseline model was adjusted considering season, CAS, FAT, and pH, in that order of importance.The potential inclusion of microbiological counts led to a better model, which specifically incorporated the level of gram-positive catalase-negative cocci (GPCNC), thus excluding pH.As expected, CY was positively associated with FAT and CAS, which is consistent with the correlations found for these quality parameters and CY in Manchega sheep milk [5].In other breeds, such as Sarda, Assaf, or Churra, this relationship has also been documented [30,36,37].CY was also negatively associated with the level of GPCNC, with the linear pattern of the relationship being confirmed across the classes of GPCNC (Figure 1).
It has been described that high microbial counts could be associated with an enzymatic load acting proteolytically and lipolytically on the milk, leading to a decrease in its fat and protein content, indirectly affecting coagulation properties [38].In terms of season, CY was higher in autumn and winter milk than in spring and summer.The higher yields obtained in these seasons could be related to the increased levels of the major components of the milk, as has also been shown in various studies [39,40].
sociated with FAT and CAS, which is consistent with the correlations found for these quality parameters and CY in Manchega sheep milk [5].In other breeds, such as Sarda, Assaf, or Churra, this relationship has also been documented [30,36,37].CY was also negatively associated with the level of GPCNC, with the linear pattern of the relationship being confirmed across the classes of GPCNC (Figure 1).It has been described that high microbial counts could be associated with an enzymatic load acting proteolytically and lipolytically on the milk, leading to a decrease in its fat and protein content, indirectly affecting coagulation properties [38].In terms of season, CY was higher in autumn and winter milk than in spring and summer.The higher yields obtained in these seasons could be related to the increased levels of the major components of the milk, as has also been shown in various studies [39,40].
Regarding A60, the best baseline model was adjusted with season and CAS, the latter being the most important factor.However, the model fit was the poorest of all, with a coefficient of determination of 8.5%.The possibility of including microbiological counts led to improving the model, incorporating BAB and SPC, and also including LAC.The level of BAB acquired a similar importance to CAS for A60.The association of A60 was positive with CAS and SPC and negative with LAC and BAB (Figure 2).SPC and BAB levels showed a good fit with linear and quadratic patterns.The associations found are consistent with the correlations for these quality parameters and A60 in Manchega [5] and Sarda [11,12] sheep milk, which indicate that the decrease in curd firmness is inversely y = 0.36x + 0.28 R² = 0.29 Regarding A 60 , the best baseline model was adjusted with season and CAS, the latter being the most important factor.However, the model fit was the poorest of all, with a coefficient of determination of 8.5%.The possibility of including microbiological counts led to improving the model, incorporating BAB and SPC, and also including LAC.The level of BAB acquired a similar importance to CAS for A 60 .The association of A 60 was positive with CAS and SPC and negative with LAC and BAB (Figure 2).SPC and BAB levels showed a good fit with linear and quadratic patterns.The associations found are consistent with the correlations for these quality parameters and A 60 in Manchega [5] and Sarda [11,12] sheep milk, which indicate that the decrease in curd firmness is inversely proportional to the increase in lactose content, with the relationship between firmer curds and the level of CAS being noteworthy.An opposite effect of the vegetative forms of total mesophilic microorganisms (SPC) and spore-forming forms in milk (BAB) on curd hardness was observed, for which no previous references have been found, and which will be appropriate to address in future studies.A 60 follows a seasonal pattern similar to CY, with increased curd firmness in autumn and winter.Todaro et al. [41] found greater hardness in spring and autumn, possibly related to the production system itself, which conditions milk quality.
The most plausible baseline model for RCT was adjusted with LAC, season, and somatic cell level, achieving a coefficient of determination of 39.1%.LAC and season were the most relevant factors for the onset of coagulation.In contrast, the model that incorporated microbiological levels was more plausible, obtaining a coefficient of determination of 44.4%.This model included the levels of SPC and BAB and also incorporated pH, this time excluding somatic cell levels.The association of RCT was positive with pH, LAC, and BAB and negative with the SPC.These results are in line with those reported by other authors [11,12], who reported an increase in coagulation times when individual sheep milk samples have a more alkaline pH, although this same study shows a negative relationship between lactose content and RCT.The mechanisms that exactly explain the association between lactose and technological traits are still unknown, although it could be associated with a modification of the percentage and composition of minerals and proteins, in accordance with the role of lactose as an osmotic regulator of milk [11,12].BAB followed a linear pattern, while linear and quadratic patterns were confirmed for SPC (Figure 3).It is confirmed that different types of microorganisms do not have the same effect on RCT, as was also observed for curd hardness; on one hand, active vegetative cells, such as SPC, consume lactose to assure growth and, consequently, decrease RCT, while spore-forming microorganisms like BAB are associated with a longer duration of coagulation, which would have certain negative connotations in optimizing the cheese-making process.Regarding season, RCT is significantly higher in spring and lower in autumn, similar to results obtained in Latxa sheep, with increased RCT from April to July [31], and in dairy cattle, with the lower value in autumn [42].Conversely, other studies reported the lowest coagulation times in spring and the highest in winter for bulk tank sheep milk [41].In any case, in our study, differences in coagulation time depending on season are much higher than those of the cited studies, which, as previously mentioned, may be related to the production system of this particular breed.proportional to the increase in lactose content, with the relationship between firmer curds and the level of CAS being noteworthy.An opposite effect of the vegetative forms of total mesophilic microorganisms (SPC) and spore-forming forms in milk (BAB) on curd hardness was observed, for which no previous references have been found, and which will be appropriate to address in future studies.A60 follows a seasonal pattern similar to CY, with increased curd firmness in autumn and winter.Todaro et al. [41] found greater hardness in spring and autumn, possibly related to the production system itself, which conditions milk quality.The most plausible baseline model for RCT was adjusted with LAC, season, and somatic cell level, achieving a coefficient of determination of 39.1%.LAC and season were the most relevant factors for the onset of coagulation.In contrast, the model that incorporated microbiological levels was more plausible, obtaining a coefficient of determination of 44.4%.This model included the levels of SPC and BAB and also incorporated pH, this time excluding somatic cell levels.The association of RCT was positive with pH, LAC, Foods 2024, 13, 886 9 of 13 [31], and in dairy cattle, with the lower value in autumn [42].Conversely, other studies reported the lowest coagulation times in spring and the highest in winter for bulk tank sheep milk [41].In any case, in our study, differences in coagulation time depending on season are much higher than those of the cited studies, which, as previously mentioned, may be related to the production system of this particular breed.The best baseline model for k 20 was adjusted with season, CAS, and pH, achieving a poor coefficient of determination of 10.4%.The model with microbiological counts added the GPCNC level and increased the coefficient of determination to 13.1%.The association of k 20 was positively associated with pH, similar to what was observed for RCT, and negatively associated with CAS, similar to GPCNC, which showed a negative linear pattern (Figure 4).The relationship between k 20 and CAS follows the trend already highlighted in previous studies, particularly when coagulation times are high (RCT > 30 min) [5].According to some authors, high microorganism counts, as in this case, GPCNC, affect coagulation by causing acidification of the milk and a faster firming of the curd [38].Similar to RCT, k 20 was significantly higher in spring and lower in autumn.Others, on the contrary, reported lower levels for k 20 in spring and winter [41], differences that could be due, as discussed for other parameters, to either the production system of Valle del Belice sheep, which is more extensive and dependent on the availability of pastures, or to differences in the stage of lactation [43,44].
previous studies, particularly when coagulation times are high (RCT > 30 min) [5].According to some authors, high microorganism counts, as in this case, GPCNC, affect coagulation by causing acidification of the milk and a faster firming of the curd [38].Similar to RCT, k20 was significantly higher in spring and lower in autumn.Others, on the contrary, reported lower levels for k20 in spring and winter [41], differences that could be due, as discussed for other parameters, to either the production system of Valle del Belice sheep, which is more extensive and dependent on the availability of pastures, or to differences in the stage of lactation [43,44].

Conclusions
The applied models demonstrated limited fitting for most milk coagulation properties, except for curd yield (CY). Incorporating microbiological analysis into the models significantly enhanced both the fit and the relative quality of the results, highlighting the importance of microbiological factors in dairy industry quality control systems.Significant associations were found between seasonal variations and all coagulation properties, even when microbiological counts were considered.During autumn and winter, the coagulation process initiates earlier (RCT) and proceeds at a slower pace (k20), resulting in firmer curds (A60) and yielding better curd production (CY).Both lactate-fermenting Clostridium spores (BAB) and total mesophilic bacteria (SPC) showed significant associations with A60 and RCT, albeit with opposing effects.While BAB concentration reduced curd firmness and increased coagulation onset, SPC concentration enhanced curd firmness and promoted coagulation onset.Coagulation speed and curd yield were significantly associated with gram-positive catalase-negative cocci (GPCNC).As GPCNC concentration increased, coagulation speed accelerated, and curd yield decreased.These results emphasize the complexity of interactions among microbiological factors and milk composition

Conclusions
The applied models demonstrated limited fitting for most milk coagulation properties, except for curd yield (CY). Incorporating microbiological analysis into the models significantly enhanced both the fit and the relative quality of the results, highlighting the importance of microbiological factors in dairy industry quality control systems.Significant associations were found between seasonal variations and all coagulation properties, even when microbiological counts were considered.During autumn and winter, the coagulation process initiates earlier (RCT) and proceeds at a slower pace (k 20 ), resulting in firmer curds (A 60 ) and yielding better curd production (CY).Both lactate-fermenting Clostridium spores (BAB) and total mesophilic bacteria (SPC) showed significant associations with A 60 and RCT, albeit with opposing effects.While BAB concentration reduced curd firmness and increased coagulation onset, SPC concentration enhanced curd firmness and promoted coagulation onset.Coagulation speed and curd yield were significantly associated with gram-positive catalase-negative cocci (GPCNC).As GPCNC concentration increased, coagulation speed accelerated, and curd yield decreased.These results emphasize the complexity of interactions among microbiological factors and milk composition from farms and their collective impact on coagulation characteristics and seasonal variations, a crucial aspect for optimizing cheese production.and J.M.P.; project administration, R.A.; funding acquisition, R.A. and J.M.P.All authors have read and agreed to the published version of the manuscript.
Funding: The present study was funded by the National Institute of Research and Agrarian and Food Technology-INIA (Government of Spain) as part of the project RTA2011-00057-C02-01.
Institutional Review Board Statement: Ethical approval was not required for this project, since all samples were bulk tank milk.
Informed Consent Statement: Not applicable.

Table 1 .
Descriptive statistics of the variables for composition, microbiological quality, and coagulation properties of sheep milk from bulk tank (n = 308).

Table 2 .
Classes of microbial counts used in linear mixed models as fixed factors (n = 308).

Table 3 .
More plausible linear mixed models (F-value and significance) based on the composition of the milk (BM), and the incorporation of microorganism counts (WSM) (n = 308).