Dairy Sheep Grazing Management and Pasture Botanical Composition Affect Milk Macro and Micro Components: A Methodological Approach to Assess the Main Managerial Factors at Farm Level

Simple Summary Studies on the management factors that affect milk components at the farm level are important for understanding how to transfer the results from experimental study. Plant phenological stages and partially fresh herbage intakes affect the lactose and milk fatty acid profile. The botanical composition of the grassland partially affects the milk’s phenol content. A few small relationships between plant phenols and milk colour could be of interest to explain the changes in milk colour parameters. Abstract The fatty acid profile, vitamins A and E, cholesterol, antioxidant power colour and the phenols profile of Sarda sheep milk from 11 commercial sheep flocks managed under permanent grassland were investigated. In each farm, the structural and managerial data and milk samples were collected during four periods (sampling dates, SD): January, March, May, and July. Data from the milk composition (fat, protein, casein, lactose, and somatic cell count), 68 fatty acids, 7 phenols, 1 total gallocatechin equivalent, ferric reducing antioxidant power, vitamins A and E, cholesterol, degree of antioxidant protection, and the colour (b *, a * and L *) were analyzed by multivariate factorial analysis using a principal component analysis approach. A proc mixed model for repeated measurement to point out the studied factors affecting significant macro and micro milk composition was also used. Only the first five components were detailed in this paper, with approximately 70% of the explained variance detected. PC1 presented the highest positive loadings for milk lactose, de novo FA synthesis and the BH intermediate, whereas OBCFA had negative loadings values. The PC2, LCFA, UFA, MUFA, vitamins E, and DAP showed positive loadings values, while SFA had a negative value. The PC3 showed a high positive loading for total phenols and non-flavonoids. PC4 presented a high positive loading for the milk macro-composition and negative values for n-3 FAs. The PC5 is characterized by high positive loadings for the a * and L * colour parameters whereas negative loadings were detected for the milk flavonoids content. These preliminary results could help to establish future threshold values for the biomarkers in milk sourced from grazing dairy sheep in natural, permanent pasture-based diets.


Introduction
Consumers are showing growing interest in the methods of production of their food as well as in environmental issues. Dairy and meat products sourced from pasture-based diets are associated with positive characteristics since it is well known that grazing confers specific organoleptic features on the dairy products with high healthiness and sensorial value [1]. Of course, some ambiguity around the terminology of "grass-fed" vs. "pasturebased" has occurred: whereby "grass-fed" is an indication of the proportion of grass in fresh or ensiled feed, whereas "pasture-based" implies that the sheep are outdoors grazing [2]. In addition, while food safety concerns increase with the food supply chain length, the sharing of knowledge and interest among stakeholders decreases, with deleterious effects on trust and the understanding between farmers, processors, retailers, and consumers. In sheep farms, an increase in grazing and self-produced hay and feedstuff could reduce the feeding cost for purchased feed. Pasture-based feeding is sustainable, safe and delivers high quality products and a range of ecosystem services [3], related to milk production in local small ruminant breeds with high levels of grazing. In this context, farmers consider permanent pastures to have a strong link with the territory origin (heritage and identity) as also reported by [4]. Finally, consumers are prepared to pay more money for milk from pasture-based production systems, but dairy companies often say that it is difficult to obtain the expected extra value from the market [5].
In this context dairy products sourced from a pasture-based diet, could represent extra added value. In Sardinia, natural pasturelands represent the main source of animal feed, covering more than one million hectares, about 87% of the total forage area of the island [6,7]. The dairy sheep farming system is mainly characterized by semi-extensive farms with a high use of natural pasture and forage crops, and contributes 66% of the total national milk sheep production [6]. Sardinian sheep farms provide cheese with a protected denomination origin (e.g., pecorino Romano, pecorino Sardo, and Fiore Sardo) as also reported by [1] in a similar context. However, in the last decades, in the face of low sheep milk prices, more and more intensive systems are being adopted in sheep farms with a high selection rate; using large amounts of concentrate and dried fodder, with less and less grazing activity, resulting in natural pasture abandonment, particularly in less favored areas [8,9].
Previous papers have demonstrated that in Sardinia, sheep diets based on forage crops and natural pasture with moderate supplementation can increase milk polyunsaturated fatty acid (PUFA) and conjugated linoleic acid (CLA) content and, consequently, the dairy products' healthiness [10][11][12][13]. This is strongly impacted by the plants' phenological stage and season, since herbage fatty acids (FA), the precursors of milk PUFAs, are high in the leaves of plants during the growing stage in the wet season [14]. Different feeding strategies that improve sheep milk quality have been tested under controlled conditions and reviewed recently; these considered zero grazing vs. full grazing crops vs. different fat supplementation [10,12,13,15].
In the Mediterranean basin, sheep milk production is seasonal (with herbage availability in late autumn-winter and early spring and a shortage in summer), with a short lambing period (November-December). In this context, a few constraints occur both for farmers (low pasture forage production in winter) and cheese factories (milk production mainly concentrated in winter and early spring). On the other hand, in late spring, grazing animals appear more economically sustainable (total dry matter intake from fresh herbage covers about 85% of the animals' requirements) with a lowering of feeding costs. This implies that several concerns exist when implementing the authentication process for pasture-based feed in sheep dairy products [16]. Recent papers underlined, at the experimental level, that several tools could be implemented to increase the authentication methods of dairy products, such as milk alkanes and fatty acid profile detection [17,18]. In the dairy sheep farming system, data on the influence of permanent pasture on milk quality traits are scanty, in particular regarding carotenoids, colour and phenol contents as reported in a recent review [13], at both controlled conditions and at farm level. To our knowledge, the effect of supplementation on milk components in sheep grazing permanent pasture are reported only partially, because the authors reported only the effect on the macro components and fatty acid profile [6], or the macro components, vitamins, and volatile organic compounds [1] and Valdivielso et al. cited by [1].
Nowadays, the availability of a large number of milk chemical components that can be detected, and the experimental factors analyzed, give some concern about the best statistical approach to implement "on farm" the best managerial practice to improve the milk quality traits. The evaluation of the simultaneous variation in milk composition due to different variables could help us to increase knowledge of the factors affecting milk's nutritional and sensorial properties; the use of multivariate statistical analysis allows one to capture the covariance structure of complex patterns of variables and create a simpler interpretation of the original multivariate system [19]. Recently, several papers were published on the effect of livestock systems on sheep milk fatty acid profile collected at the farm level, [20,21] considering only partially the farm managerial factors and often without any references to the grassland's botanical composition detected during the survey.
The aim of this paper is to investigate, in extensive dairy sheep farms, the main relationships between a few structural (stocking rate) and managerial (animal diet composition, pasture botanical composition, plant phenological stage, total supplementation) factors, and the milk's macro-composition (fat, protein, casein, somatic cell count, lactose) and micro-composition (FA profile, vitamins, colour, phenols and other antioxidant parameters). The study aims to assess a data methodological approach in a case study at the commercial farm level, to aggregate a large number of variables detected in sheep milk composition into fewer latent structures, namely the principal components (PC). These latent structures are subsequently related to some managerial factors to point out the milk biomarkers that could be useful to test the link between dairy products and the territory of origin. This approach is a first step in providing a small dataset at the farm level, and later needs to be validated with a large database that takes into account the effect of years. Our basal hypotheses to test at farm level are: (i) does the grassland's botanical composition and plant phenological stage interact with feeding supplementation to affect the milk's macro and micro chemical composition? (ii) Is it possible to synthetize into a small group of factors a large number of variables (about 100) related to the milk's macro-and micro-composition?

Dairy Farming System and Pasture Botanical Characteristics
The study was performed in 2019, in 11 dairy sheep (Sarda breed) farms located in a hilly central area of Sardinia (historic region of Marghine) (300-700 m a.s.l.). The Sarda breed is a native medium-sized sheep (45 kg of live body weight), mainly raised on a semi-extensive system based on permanent crop grassland. The average milk production is 220 L in 180 days of lactation. In each farm, milk samples were collected during 4 sampling days (SD), in January, March, May and July 2019 when all of the sheep had 58, 98, 138 and 178 days in milk (DIM), respectively. The SD were carried out every 40 days to better link the milk yield and composition to the evolution of the botanical composition phenological stage and the fresh herbage availability during the lactation period of the dairy sheep in the Mediterranean basin. In total, 42 milk samples were collected, as on two occasions the samples were lost. At every SD, structural and managerial data were also recorded: Utilised Agricultural Area (UAA), flock size (number of sheep), and feeding management, such as hay and concentrate supplementation, access time to grazing per day, and the forage conservation systems. Animal diet composition was indicated as the intake of (i) grazed herbage in natural pasture (HeI); (ii) hay (HI); (iii) concentrate (CI); and (iv) standing hay. Pasture herbage intake (Hel) was estimated by the difference between the potential intake capacity of the animals at the different SD and the encumbrance provided by the other feedstuffs fed to the sheep [22]. The potential intake capacity was estimated according to the herd's characteristics (milk yield and composition, and animal body weight) as detailed in the small ruminant nutrition system (SRNS) [23]. A complementary study was carried out to evaluate the botanical composition of the pasture grazed by the sheep on the same days as milk sample collection. In each SD, the herbage mass was determined by cutting the sward at 3 cm above ground level on a sampling area of 0.5 m 2 . The pasture botanical composition (BC) was detected by partitioning the herbage samples into Legumes, Grasses, Forbs species. The plant phenological growth stages (PPS) were detected at each SD: (i) Growing stage (GW, from 0 to 5 (Biologische Bundesanstalt, Bundessortenamt und CHemische Industrie) BBCH scale; Refs. [24,25] characterized by plants in active growth and with vegetative activity; (ii) flowering stage (FW, from 5 to 6 BBCH; Refs. [24,25] characterized by plants with inflorescence or flower buds visible until the end of flowering; (iii) maturity or senescence stage (MS, from 8 to 9 BBCH; Ref. [25] characterized by no vegetative activity, leaf fall and with most of the plants dead or dormant. The fresh forage was dried in an oven-dryer at 65 • C for three days to evaluate the proportion (%) of each partitioned group in terms of dry matter (DM).

Milk, Herbage and Feedstuff Chemical Characteristics
The ewes were milked twice a day (in the evening and morning) using automatic milking machines. In each SD, 1 bulk milk sample (1.5 L) deriving from the daily milking was taken from the tank of each farmhouse. Five sub samples were obtained from each farm: one for macro-composition, one for vitamins, cholesterol and the degree of antioxidant protection (DAP), one for the phenols profile, one for ferric reducing antioxidant power (FRAP), one for gallic acid equivalent detection (GAE), and one for the fatty acid (FA) profile. The fat, protein, casein, somatic cell counts (SCC) and lactose contents were determined by Milkoscan FT+, Foss Electric, (Hillerød, Denmark). Milk color was measured by Chroma Meter CR400/410 (Konica Minolta Sensing Japan) with a spectral sensitivity characteristic of the sensor and the spectral distribution of the xenon lamp. Three readings per sample were performed. The white calibration was performed with a white calibration plate where Y = 93.5, x = 0.3114, and y = 0.3190. The following parameters were detected: L * (light-ness; on a scale from 0 to 100, where 0 = black and 100 = white), a * (where −a * has a green color and +a * has a red color), and b * (where −b * has a blue color and +b * has a yellow color).
For the milk fatty acid profile detection, milk fat was obtained in agreement with [26]. Fatty acid methyl esters from milk fat (0.05 g of milk fat sample) were obtained according to the international standard (FIL-IDF, 1999), using base-catalyzed trans-methylation. Separation of FAMEs was performed [26]. To achieve a compromise between an acceptable separation of the cis and trans isomers (C16:1, C18:1 C18:2 and C18:3) and the time of analysis, one of the oven temperature programs proposed by [27] was used. Individual FAMEs were identified by comparison with a standard mixture of 38 FA methyl esters (Matreya Inc., Pleasant Gap, PA, USA). C18:1 isomers standards (C18:1t9, C18:2 t9t12, C18:1t11, Matreya Inc., Pleasant Gap, PA, USA) and CLA standards (CLAc9t11; CLAt10c12; CLAc9c11; CLAt9t11, Matreya Inc., Pleasant Gap, PA, USA) and published isomeric profiles [28] were used to identify the C18:1 isomers, C18:2 non conjugate isomers and the CLA isomers. Calibration curves with internal standards were performed to quantify fatty acid methyl esters (FAMEs) that are calculated as g FAME/g of milk fat and then expressed as g FAME/100 g of total identified FAMEs. The internal standard Me-C5:0 was used to quantify from Me-C4:0 to Me-C7:0, the internal standard Me-C9:0 was used to quantify from Me-C8:0 to Me-C10:0, the internal standard Me-C13:0 was used to quantify medium chain FAMEs from Me-C11:0 to Me-C17:0 and the internal standard Me-C19:0 to quantify long chain FAMEs from Me-C18:0 to Me-C26:0. The concentration of each internal standard, added to the milk fat sample prior to transmethylation, was 0.005 g in 0.050 g of milk fat.
Vitamin A, E and cholesterol were analyzed according to [29]. Briefly, aliquots (2 mL) of milk were digested with 2 mL of KOH (60% aqueous solution, w/v), 2 mL of 95% ethanol, 1 mL of NaCl (1% aqueous solution, w/v), and 5 mL of an ethanolic solution of pyrogallol (6%, w/v) added as an antioxidant. After digestion in a water bath at 70 • C, the suspension was cooled for 30 min, and 5 mL of an NaCl solution (1%, w/v) was added to prevent emulsification. The suspension was then extracted with 10 mL of n-hexane/ethyl acetate (9:1, v/v). The lower aqueous layer was extracted 3 more times, with 5 mL of n-hexane/ethyl acetate (9:1, v/v) each time. The pooled organic layers were evaporated with a rotary evaporator at 30 • C, and the dry sample was dissolved in 3 mL of methanol for HPLC. A sample volume of 20 µL was injected in HPLC equipment, previously filtered using a 0.20 µm PTFE filter. All determinations were carried out in duplicate. The degree of antioxidant protection (DAP), which represent a synthetic index, was calculated as the ratio between the amount of antioxidant element (e.g., α-tocopherol) and the element to be protected against oxidation (cholesterol) according to [30].
The original FRAP assay [31] protocol was modified to make the assay suitable for raw milk according to [32]. The FRAP reagent was prepared by mixing 300 mM sodium acetate buffer (pH 3.6) and 10 mM TPTZ (2,4,6-Tris(2-pyridyl)-s-triazine) in a 40 mM HCl and 20 mM iron (III) chloride solution in a volume ratio of 10:1:1, and the mixture was warmed to 37 • C in a water bath before use. Afterward, 900 µL of the prepared FRAP reagent was mixed with 30 µL of a diluted sample, and an absorbance at 593 nm was recorded after 5 min of incubation at 37 • C. A standard curve was constructed using iron sulfate heptahydrate (FeSO 4 ·7H 2 O), and the data were expressed as millimoles of FeSO 4 equivalents per kilogram of cheese. All determinations were carried out in duplicate. Gallic acid equivalents (GAE) were analyzed according to [33] using the Folin-Ciocalteu colorimetric method. Briefly, 100 µL of the diluted extract was mixed with 500 µL of a 0.2 N Folin-Ciocalteu reagent. Afterward, 400 µL of a sodium carbonate solution (7.5% aqueous solution, w/v) was added to the reaction mixture. The absorbance readings were taken at 765 nm after incubation at room temperature for 1 h. Gallic acid was used as a reference standard, and the results were expressed in milligrams of gallic acid equivalents (GAE)/liter of milk.
At the same time as the milk sampling, randomized samples of the herbage from grazed pasture (samples of the pasture on offer) and supplied feedstuffs were sampled in each farm and stored at −20 • (1999). The NDF (neutral detergent fiber assayed with a heat stable amylase and expressed exclusive of residual ash), ADF (acid detergent fiber, expressed exclusive of residual ash) and ADL (acid detergent lignin) were determined [37] using an Ankom 220 Fiber Analyzer equipment (Ankom Technology, New York, NY, USA). The powdered material was used for the extract preparation of the total phenolics [38]. The total phenolics, were analyzed using Folin-Ciocalteau reagent and expressed as a tannic acid equivalent (mg TAE g −1 DM) according to procedures previously described [38]. Briefly, the samples (0.5 g each) were incubated in a solution of acidified methanol and then filtered. The filtrate was added to a solution including Folin-Ciocalteau reagent and sodium carbonate. After two hours, the samples were assayed for total polyphenols by a spectrophotometer (765 nm) using tannic acid as a standard. To measure non-tannic polyphenols the same procedure was used but the filtrate was pretreated with methylcellulose and ammonium sulphate before adding the Folin-Ciocalteau reagent. Tannic phenolics were measured as the difference between total and non-tannic polyphenols.
A principal component analysis (proc princomp performed by SAS) was used to reduce the original number of variables to a new set of fewer variables (principal components, PC). The number of PCs retained was defined according to the amount of explained variance and by adopting Kaiser's Criterion (eigenvalue > than 1). The PCs loadings were interpreted as the coefficients of the linear combination with the original variables. Principal component scores (n = 42) were analyzed using the proc mixed model for repeated measurement (SAS) to investigate the structural and managerial studied factors affecting significant macro and micro milk composition. The means and significance of the factors were detected at p < 0.05 using Tukey test. In order to stress the effects of the different phenological stages, the data of GW and FW stages were considered in the same group as the "vegetative phase" and compared to the MS stage. In the same way, in the pasture botanical composition (BC), the occurrence of forbs and legumes species, less present in pasture, were considered in the same group (forbs) to compare with the grasses species group. The model adopted for repeated measurement analysis was:  (1999). The NDF (neutral detergent fiber assayed with a heat stable amylase and expressed exclusive of residual ash), ADF (acid detergent fiber, expressed exclusive of residual ash) and ADL (acid detergent lignin) were determined [37] using an Ankom 220 Fiber Analyzer equipment (Ankom Technology, New York, NY, USA). The powdered material was used for the extract preparation of the total phenolics [38]. The total phenolics, were analyzed using Folin-Ciocalteau reagent and expressed as a tannic acid equivalent (mg TAE g −1 DM) according to procedures previously described [38]. Briefly, the samples (0.5 g each) were incubated in a solution of acidified methanol and then filtered. The filtrate was added to a solution including Folin-Ciocalteau reagent and sodium carbonate. After two hours, the samples were assayed for total polyphenols by a spectrophotometer (765 nm) using tannic acid as a standard. To measure non-tannic polyphenols the same procedure was used but the filtrate was pretreated with methylcellulose and ammonium sulphate before adding the Folin-Ciocalteau reagent. Tannic phenolics were measured as the difference between total and non-tannic polyphenols.
A principal component analysis (proc princomp performed by SAS) was used to reduce the original number of variables to a new set of fewer variables (principal components, PC). The number of PCs retained was defined according to the amount of explained variance and by adopting Kaiser's Criterion (eigenvalue > than 1). The PCs loadings were interpreted as the coefficients of the linear combination with the original variables. Principal component scores (n = 42) were analyzed using the proc mixed model for repeated measurement (SAS) to investigate the structural and managerial studied factors affecting significant macro and micro milk composition. The means and significance of the factors were detected at p < 0.05 using Tukey test. In order to stress the effects of the different phenological stages, the data of GW and FW stages were considered in the same group as the "vegetative phase" and compared to the MS stage. In the same way, in the pasture botanical composition (BC), the occurrence of forbs and legumes species, less present in pasture, were considered in the same group (forbs) to compare with the grasses species group. The model adopted for repeated measurement analysis was: where: Yijhtm = is the studied variable; µ = is the mean; ijhtmo (1) where: Yijhtm = is the studied variable; µ = is the mean; HeI = is the fixed effect of the herbage intake (HeI, i = high, medium and zero); PPSj = is the fixed effect of the phenological stage (PPS, j = GW/FW or MS stage); BCh = is the fixed effect of plant botanical composition (BC = H = grasses or forbs); SDt = is the repeated fixed effect of the sampling date (SD, t = D1, D2, D3 and D4); Fm = is the random effect of farm (F, m = 1 to 11) . The NDF (neutral detergent fiber assayed with a heat d exclusive of residual ash), ADF (acid detergent fiber, al ash) and ADL (acid detergent lignin) were determined r Analyzer equipment (Ankom Technology, New York, NY, l was used for the extract preparation of the total phenolics analyzed using Folin-Ciocalteau reagent and expressed as AE g −1 DM) according to procedures previously described g each) were incubated in a solution of acidified methanol was added to a solution including Folin-Ciocalteau reagent two hours, the samples were assayed for total polyphenols m) using tannic acid as a standard. To measure non-tannic ure was used but the filtrate was pretreated with methyllphate before adding the Folin-Ciocalteau reagent. Tannic he difference between total and non-tannic polyphenols.
is, the data on FA composition were processed to calculate rt chain fatty acids (SCFA C4:-C11:0); medium chain fatty FA (C18:0-C26:0); odd and branched chain FA (OBCFA = 14:0 + isoC15:0 + anteisoC15:0 + C15:0 + isoC16:0 + isoC17:0 nounsaturated fatty acids (MUFA all FA with single double re than one double bond); saturated fatty acids (SFA all FA turated fatty acids (UFA all FA with one or more double nols were also calculated as flavonoids (cyanidin + luteolin), ) and non-flavonoids (ferulic acid + tyrosol + sesamin). ion (Fat, protein, casein, Lactose and SCC), 68 FA (expressed , 1 total gallocatechin equivalent, FRAP Vitamin A and E, (L *, a * and b *) were analyzed with PCA SAS PRINCOMP Cary, NC, USA). analysis (proc princomp performed by SAS) was used to of variables to a new set of fewer variables (principal r of PCs retained was defined according to the amount of dopting Kaiser's Criterion (eigenvalue > than 1). The PCs the coefficients of the linear combination with the original t scores (n = 42) were analyzed using the proc mixed model AS) to investigate the structural and managerial studied macro and micro milk composition. The means and e detected at p < 0.05 using Tukey test. In order to stress the gical stages, the data of GW and FW stages were considered etative phase" and compared to the MS stage. In the same l composition (BC), the occurrence of forbs and legumes re, were considered in the same group (forbs) to compare p. The model adopted for repeated measurement analysis = µ + HeI + PPSj + BCh + SDt + Fm + Ɛijhtmo le; ijhtmo = is the error term.
Total supplementation and stocking rates were not considered in the model since any significance with the PC scores was detected. Finally, Person's index correlation (r) was studied between the significant factors pointed out by repeated measurement analysis and the original milk composition variables with the higher loadings. Person's index correlation (r) between the milk's macro and micro chemical components that showed the highest loadings were also investigated.

Dairy Farming System and Feeding Managements
On average, the UAA was around 68 ha, and the flock size was 204 head of Sarda sheep ( Table 1). The total supplementation amount (hay and concentrate) supplied by the sheepherders is a complement of pasture herbage diet contribution and it depends on pasture herbage availability, milk yield, and the animals' physiological stage in terms of DIM as also reported by [8]. The MS plant stage led to a decrement of pasture palatability and nutritive value of the total daily diet intake as a consequence of the reduction in the leaves to steam ratio that led to an increase of the ADF as reported by [14]. The forage (herbage + hay)/concentrate ratio in the animal diet changed throughout the whole survey period, with the lowest values in January (70% of forage in the total diet) and the highest in May (more than 80% of forage in the total diet) when the herbage availability and sheep grazing activity were maximum. To outline, all of the sheepherders followed the same feeding strategy, using a higher quantity of hay and concentrate during the winter than in spring in response to the lack of herbage in the permanent pasture, and the higher amount of available pasture herbage or standing hay in spring and summer, respectively. This was also reported by [40] in similar extensive dairy sheep farms located in the Basque Country. The estimated pasture herbage intake and total supplementation ranged between 0 and 1777, and 100 and 1850 g DM/day/head, respectively, as reported in Table 1.

Grassland Botanical and Chemical Composition
The main representative forage species in pasture were grasses composing, on average, 65% pasture DM based, followed by forbs (15% DM) and legumes (10% DM), with different patterns of evolution from January to May. The frequency of the botanical families was similar to the results found by [7] in the same study area. The proportion of legumes and forbs increased, whereas grass species decreased in the pasture from January to May. During January and March, the plants were characterized by GW stages in agreement with their growing degree days linked to environmental drivers, such as temperature and rainfall [41]. In May, most plants were in the FW stage, whereas in July, most plants were in the MS stage without vegetative activity. The ether extract decreased from January to July as has also been reported previously [14]. The total phenols increased from January to May in parallel with the legumes and forbs pattern in the grassland, in agreement with [42]; a drop (30%) was observed in July compared to March (MS and GW, respectively) according to [38]. The pasture CP content significantly (p < 0.05) decreased from March (GW) to May (FW), passing from 16.69 to 13.16 CP% of DM mainly due to the natural drop of CP content in the grasses component; in addition, the ADF content increased mainly from March (GW) to May (FW) by 23% (p < 0.05). Overall, the observed reduction in terms of EE, CP and quality of fiber from vegetative to reproductive stage changes is probably more linked to the high value of the leaves to stems ratio during GW than the FW or MS stage, as also reported by [14]. The evolution pattern of the grasslands' botanical composition and their phenological stage affect the chemical composition of pasture with different magnitude. In fact, it is known that in the Mediterranean basin, when plants turn from the vegetative to the reproductive stage, the worse effects on CP drop and NDF increase is for grasses compared to legumes, as also reported by Cabiddu et al. [10].

Milk Chemical Composition and Relationship with Structural and Managerial Factors
In previous years, different attempts were carried out to explore the datasets from experimental studies using the multivariate approach, to evaluate the association between feeding management and milk fatty acid profiles [19][20][21]. In this study, we consider the "on farm" level additional factors that are not common in the bibliography, such as the natural pasture's botanical composition and the estimated animal intake of herbage, hay and concentrate. However, the total number of observations is not very large (n = 42) due to the difficulty in collecting the data for such a wide number of variables.
Results from PCA showed that 16 out of 101 PC were found to be able to explain about 92% of the total variance of the system. The variance explained ranged from 37.00% for PC1 to about 1.00% for PC16. Only the first five components were detailed in this paper that explained around 70% of the total variance ( Table 2). The five components (PCs) were chosen considering a threshold of 0.13 (absolute value) of eigenvalues.   Table 1.
The average values of the milk FA content confirmed, partially, the previous reports on extensive dairy sheep farm systems located in a few areas of Sardinia with different forage crops and a permanent pasture-based diet [11,43], whereas the milk total n-3 FAs (Table 1) appeared lower when compared to the milk of sheep grazing on crops [44].
SFA, MUFA and PUFA accounted for almost 68%, 25% and 7% of the total FAs, respectively. Palmitic, oleic, and linoleic acid are the main SFA, MUFA, and PUFA, respectively, in agreement with results from [43] in sheep, and [19] in cow's milk.
The PC1, (variance explained 37%) presented the highest positive loadings for lactose, de novo FA synthesis (SCFA) and biohydrogenated FA intermediates, whereas OBCFA had negative loadings values (Table 2). Figure 1 represents the plot of loadings between PC1 and PC2 and shows a high correlation between the milk yield (MY) and lactose (loadings were 0.118 and 0.136, respectively), in agreement with [45], which underlined the role of milk lactose as an indicator of energy balance intake in line with the distribution of oleic acid (C18:1 c9) located on the opposite side of lactose and MY. Therefore, as expected, animals in a positive energy balance could increase MCFA and SCFA in milk as also found by [46,47], according with the high availability of acetate and 3-hydroxybutyrate, which are the precursors of lipid de novo synthesis at the mammary level.
1. The average values of the milk FA content confirmed, partially, the previous reports on extensive dairy sheep farm systems located in a few areas of Sardinia with different forage crops and a permanent pasture-based diet [11,43], whereas the milk total n-3 FAs (Table  1) appeared lower when compared to the milk of sheep grazing on crops [44].
SFA, MUFA and PUFA accounted for almost 68%, 25% and 7% of the total FAs, respectively. Palmitic, oleic, and linoleic acid are the main SFA, MUFA, and PUFA, respectively, in agreement with results from [43] in sheep, and [19] in cow's milk.
The PC1, (variance explained 37%) presented the highest positive loadings for lactose, de novo FA synthesis (SCFA) and biohydrogenated FA intermediates, whereas OBCFA had negative loadings values (Table 2). Figure 1 represents the plot of loadings between PC1 and PC2 and shows a high correlation between the milk yield (MY) and lactose (loadings were 0.118 and 0.136, respectively), in agreement with [45], which underlined the role of milk lactose as an indicator of energy balance intake in line with the distribution of oleic acid (C18:1 c9) located on the opposite side of lactose and MY. Therefore, as expected, animals in a positive energy balance could increase MCFA and SCFA in milk as also found by [46,47], according with the high availability of acetate and 3-hydroxybutyrate, which are the precursors of lipid de novo synthesis at the mammary level.   These results partially agree with [48,49], who found a strong link between herbage intake and cow's milk and lamb's meat fatty acid composition, respectively. Loadings detected for SCC between the PC1-PC2 showed low values, and this agrees partially with These results partially agree with [48,49], who found a strong link between herbage intake and cow's milk and lamb's meat fatty acid composition, respectively. Loadings detected for SCC between the PC1-PC2 showed low values, and this agrees partially with [50], who observed that a pasture-based diet did not affect milk SCC. The high correlation between SCFA and lactose (0.78; p < 0.01), SCFA and MY (r = 0.68; p < 0.01), and lactose and MY (r = 0.76; p < 0.01) pointed out the effect of a high energy intake, associated with the high availability of acetate and 3-hydroxybutyrate, both contribute to increase milk SCFA by de novo synthesis [51], in agreement with [19,21,45], and in line with the previous studies in cows [19], buffalo [52], and in sheep [21]. The significant effect (p < 0.02) of herbage intake (HeI) on the PC1 scores (Table 3) probably partially explains the above results, and the position in Figure 1 of lactose, MY and de novo FA change, along the PC1 axis; in particular, milk total SCFAs are positively correlated to fresh herbage intake levels (r = 0.73; p < 0.01). A significant effect of plant PS on PC1 scores was detected too (Table 3), that could be due to an increase of the herbage ADF/NDF ratio (from vegetative to reproductive plant stage) with an implication for rumen microbial activity and, consequently, on milk OBCFA, as also reported by [39] in agreement with the positive correlation (r = 0.86; p < 0.01) between OBCFA milk content and plant phenological stage. Plant PS is a key driver of animal metabolism, effecting changes on ruminal and milk components, in particular for FA, as also reported by [53]. A negative correlation was found between milk SCFA content and PS (r = −0.82; p < 0.01), with a higher level seen when sheep were grazing plants at the GW stage, compared to the FW stage, when fiber digestibility decreased [54], and with a detrimental effect on FA precursor availability for de novo synthesis at the mammary level. The change in the diet's botanical composition (from grass to legumes and forbs) also affected the PC1 score (p < 0.01), probably as a consequence of plant secondary metabolites (PSM), which occur in legumes and forbs plants species more than in grass species [55]. As is well known, PSM in the diet could affect ruminal microbial activity with several effects on the milk's FA profile and other micro-components, likewise phenols. In this study, we found a low negative value relationship between milk LCFA and the occurrence of grass in pasture (r = −0.33; p < 0.03) with no clear explanation, probably due to their low loading value for PC1. In addition, no effect on milk phenol profiles were detected when fresh pasture contribution in the diet increased. Increased legumes and forbs are also well represented in Figure 1, where phenols are mainly located near the origin of the axis.
The PC2 (variance explained 16.1%) had a negative loadings values for SFA and MCFA, whereas positive loadings values were detected for DAP, LCFA and UFA, and MUFA. No significant effect was detected between the PC2 scores and herbage intake, plant phenological stage or plant botanical composition (Table 4). However, the positive correlation found between fresh herbage intake and milk linolenic acid content (r = 0.42; p < 0.01) partially explains the high and opposite loadings values of UFA and SFA, thanks to the plant's PS (GW and FW) with a pivotal effect on the fatty acid precursor content in the herbage [14,53]. In addition, a good correlation was found between herbage intake and milk C18:2 t11c15 content (r = 0.60; p < 0.01) and between linolenic acid and C18:2 t11c15 milk content (r = 0.72; p < 0.01); it is well known that an increase in herbage intake leads to an increase in linolenic acid (ALA) and C18:2t11c15 in milk, as also reported by [56]. The BH intermediate (in particular, trans C18:1 isomers) with positive loadings were all related to the availability in the diet of linolenic acid (ALA), in agreement with [57], even if linoleic acid (LNA) and ALA showed low loadings values for PC2. A positive correlation between total trans C18:1 trans isomers (r = 0.70; p < 0.01) and milk ALA content was detected due to a higher intake of daily fresh pasture herbage, which increases ALA intake, as also reported by [20,57] and is associated with an increase of total C18:1 trans [56]. In spite of any significant effect of hay or concentrate supplementation on PC scores, we found a positive correlation between hay intake and C15:0i, (r = 0.52; p < 0.01) as a probable consequence of the association effect of the concentrate. A high level of diet concentrate could lead to a detrimental effect on ruminal activity and OBCFA milk content, as also reported by [39]. In this case, it is very difficult to explain how the concentrate influences the milk component since no relationship was found between the forage to concentrate ratio and the OBCFA milk content. For instance, the C18:2 t9c12 appears positively correlated with the concentrate intake (r = 0.61; p < 0.01), and with the forage to concentrate ratio (r = 0. 32; p < 0.05) despite its low loadings values in both PC1 and PC2 (0.07 and 0.09, respectively).
The PC3 (variance explained 7.3%, Table 2) was characterized by high positive loadings values for ferulates non-flavonoids, total phenols, and CLA c9t11. A positive significant correlation was found between herbage intake and ferulates (r = 0.37; p < 0.05) and the CLA c9 t11 (r = 0.55; p < 0.01) in milk, respectively. This agrees with the positive correlation found between milk non-flavonoids, and ferulate with a n-3/n-6 ratio (r = +0.51 and 0.39; p < 0.05, respectively). Ferulates are the phenylpropanoid compounds that occur, in particular, on grass cell walls and in small amounts in legume cell walls as a main component of lignin [58]. These components are mainly responsible for ruminal fiber degradation, even if during the vegetative stage they are more degraded than in the plant at the reproductive stage. Between 10-50% of ingested lignin, following ruminal digestion, could be absorbed, as reported by [59]. These phenolic compounds, released at the ruminal level [60], showed a toxic microbial effect that could be the reason for the protective effect against the ruminal biohydrogenation of PUFA, which led to an increase in the n-3/n-6 ratio in the milk. This is also in line with the relationship found between ferulic acid and ALA (r = 0.41; p < 0.01) since linolenic acid represents the main fatty acid in fresh herbage. A negative correlation was found between CLA c9t11 and the plant's phenological stage, when passing from the growing and flowering to the maturity stage, as it probably resulted in a decrease in the precursor [56].
The PC4 (variance explained 4.9%) is characterized by a positive loading for milk macro-composition (fat protein and casein) and C18:1 t10. Negative loadings were detected for ALA, EPA and DHA (see Table 2). A ruminants concentrate-based diet, with a low level of fiber, mostly exposes them to increases in trans 18:1 isomers and, in particular, C18:1 t10 [61]. In our study, no relationship was found between concentrate intake and total trans C18:1 fatty acid contents in milk, probably because of the high F/C ratio detected through the whole lactation period. Moreover, in our case, the level of rapidly fermentable carbohydrates (RFCHO) (starch + sugar) are probably low through the whole period (below 40% DM intake), seriously affecting the level of biohydrogenated intermediates, in agreement with [62].
The PC5 (variance explained 4.7%) is characterized by high positive loadings for the a * and L * colours parameters, while the total flavonoids (sum of cyanidin and luteolin) showed negative loadings. As expected, the L *, a * and b * parameters (Table 1) are lower (in particular, the b * values) than the values found in grazing cows [13], in line with [63], probably because, in sheep, there is a higher conversion rate of β-carotene to Vit. A than in cows [63]. Considering previously results found in grazing (4 h/day) sheep on natural pasture [64], we found high values of b * (+21%), whereas the L* and a* values were lower (−8% and −65%, respectively). The occurrence in the permanent pasture of plant secondary metabolites as phenols could also be responsible for the low value of a * and the decrease in brightness (low value of L *). A slight positive correlation was found between milk fat and the colour parameters (r = 0.57, 0.54 and 0.46 for b *, a * and L *, respectively, p < 0.01), since milk carotenoids (lutein) are mainly located in the fat, as also reported by [65]. The permanent pasture-based diet could affect the intake of plant secondary metabolites [63] with subsequent interference upon ruminal metabolism. Color parameters are influenced by several molecules, such as carotenoids and phenols [65]. As expected, a little negative correlation was found between Vit. A and the b * parameter (r = −0.35 p < 0.05). Considering the phenols profile, only sesamin appeared to be positively correlated with b* (r = +0.34; p < 0.02); sesamin belongs to the lignans component and is characterized by being able to be absorbed at 550 nm near the yellow region, as reported by [66]. As such, it is probable that it is not only carotenoids that could affect the b * parameter in sheep milk. Even if difficult, considering the results from the recent study, we assume that our results partially agree with [67], who found effects on the cheese colour parameter when the authors compared pasture vs. indoor feeding, although the % of legumes and forbs detected in the pasture was high. In addition, more herbage in the diet is positively correlated with n-3 FAs in the milk's content. Moreover, it appears that flavonoids correlate negatively with total n-3 FAs (r = 0.30; p < 0.05) and positively with total n−6 FAs (r = 0.29; p < 0.05).

Conclusions
This study focuses at the farm level in an attempt to better evaluate the relationship between a few farm management choices and milk quality in grazing dairy sheep reared on permanent grassland. Moreover, it provides a tool by which to transfer research knowledge to farmers (even if the numbers of observations were restricted compared to the total variables considered) and provides new suggestions for research models. The fresh herbage intake and plant phenological stage appear to be the main drivers for OBCFA and SCFA and lactose content. Considering the other micro components in milk, this study pointed out the negative relationship between total supplementation and milk phenols content, in particular ferulates and non-flavonoids. The relationship between milk colour parameter (b *) and sesamin, and between the milk phenols and colour parameters could be considered in an authentication procedure to link animal products with the territorial brand, although this needs to be confirmed in experimental conditions. Our basal hypotheses are only partially explained by these results. Regarding the first hypothesis, the dataset exploration needs to be increased, taking into account the effect of years and to validate the results. For the second hypothesis, the methodological approach evaluated is a promising challenge to study the effect of farming practices on several traits of milk sheep reared on permanent grassland. In conclusion, even if these preliminary results are derived from commercial farms, they could help to establish future threshold values of biomarkers in milk sourced from grazing dairy sheep with natural, permanent pasture-based diets. These results need to be confirmed at the farm level using a larger database.