**Survival of Selected Pathogenic Bacteria during PDO Pecorino Romano Cheese Ripening**

**Giacomo Lai 1, Rita Melillo 2, Massimo Pes 1, Margherita Addis 1, Antonio Fadda <sup>2</sup> and Antonio Pirisi 1,\***


Received: 3 November 2020; Accepted: 2 December 2020; Published: 7 December 2020

**Abstract:** This study was conducted to assess, for the first time, the survival of the pathogenic bacteria *Listeria monocytogenes*, *Salmonella* spp., *Escherichia coli* O157:H7, and *Staphylococcus aureus* during the ripening of protected designation of origin (PDO) Pecorino Romano cheese. A total of twenty-four cheese-making trials (twelve from raw milk and twelve from thermized milk) were performed under the protocol specified by PDO requirements. Sheep cheese milk was first inoculated before processing with approximately 10<sup>6</sup> colony-forming unit (CFU) mL−<sup>1</sup> of each considered pathogen and the experiment was repeated six times for each selected pathogen. Cheese composition and pathogens count were then evaluated in inoculated raw milk, thermized milk, and cheese after 1, 90, and 150 days of ripening. pH, moisture, water activity, and salt content of cheese were within the range of the commercial PDO Pecorino Romano cheese. All the cheeses made from raw and thermized milk were microbiologically safe after 90 days and 1 day from their production, respectively. In conclusion, when Pecorino Romano cheese is produced under PDO specifications, from raw or thermized milk, a combination of factors including the speed and extent of curd acidification in the first phase of the production, together with an intense salting and a long ripening time, preclude the possibility of growth and survival of *L. monocytogenes*, *Salmonella* spp., and *E. coli* O157:H7. Only *S. aureus* can be still detectable at such low levels that it does not pose a risk to consumers.

**Keywords:** cheese safety; foodborne pathogens; sheep milk; *Listeria monocytogenes*; *Salmonella* spp.; *Escherichia coli* O157:H7; *Staphylococcus aureus*; raw milk; thermization

#### **1. Introduction**

Consumer demand for unpasteurized milk cheeses is constantly increasing because of their more intense flavor and varied aroma than those of pasteurized milk cheeses [1–3]. However, especially when made from unpasteurized milk, cheese can hold a risk for the consumer, because of the possible presence of some pathogenic bacteria such as *Listeria monocytogenes*, *Salmonella* spp., *Staphylococcus aureus*, *Campylobacter* spp., *Brucella* spp., and pathogenic *Escherichia coli* [3–5]. Based on a data collection of dairy products-associated outbreaks in the United States from 1993 to 2006, Langer et al. [6] determined that outbreaks attributed to the consumption of unpasteurized dairy products were approximately 150 times more frequent, based on the unit of consumption than those related to pasteurized milk or pasteurized milk products. Likewise, Costard et al. [7] reported that in the United States, from 2009 to 2014, the consumption of unpasteurized dairy products caused 840 times more disease and 45 times more hospitalizations than that of pasteurized products.

Among the pathogenic agents, *L. monocytogenes*, *Salmonella* spp., Shiga toxin-producing *E. coli* (STEC) as the serotype O157:H7 and enterotoxin-producing *S. aureus* are the most involved in foodborne outbreaks related to the consumption of raw milk cheese in industrialized countries [3,5]. These foodborne pathogens usually cause disease with acute symptoms restricted to the gastrointestinal tract such as diarrhea, abdominal cramps, nausea, vomiting, and limited in time and severity [8]. However, in some cases, they can cause serious diseases such as hemolytic uremic syndrome and thrombotic thrombocytopenic purpura associated to *E. coli* O157:H7 or meningitis and septicemia caused by *L. monocytogenes*, with a significant mortality rate in vulnerable groups such as infants, elderly, and immunocompromised adults [9,10].

These pathogenic bacteria can originate in raw milk by direct excretion from animals infected udder, by fecal contamination or from the farm environment more generally [5,11]. The foodborne pathogens can reach cheeses via contaminated milk or through the dairy plants environment and the processing equipment [4]. A recent survey shows how *L. monocytogenes* was first detected in a newly established cheese-making facility just nine months after the starting of the production [12]. Moreover, some pathogens like *L. monocytogenes* and *S. aureus* can colonize abiotic surfaces by forming biofilms, that make the bacteria immune to the action of antimicrobial agents [13,14]. Thus, they can persist for a long time in the manufacturing environment, where these microorganisms could be a potential cause of cross-contamination of the dairy products. Workers can also be an important source of contamination as a result of improper handling and some of them could be asymptomatic carriers of *S. aureus* [15,16].

The growth and survival of microbial pathogens during cheese-making is hindered by several factors such as milk heat treatments, curd cooking, rate of curd acidification by starter cultures, final product pH, salt addition, and competition with the native microflora present in milk [17]. Despite these hurdles, foodborne pathogens can express an adaptive response to different sublethal stresses, essential for their survival in harsher environments. Moreover, pathogenic bacteria adapted to a sublethal stress may exhibit cross-protection with enhanced resistance to different stresses [18–20]. Several challenge studies have shown that some foodborne pathogens are able to survive during the manufacturing and ripening of different kinds of cheese made from raw milk [9,21–23]. Some authors especially report that *L. monocytogenes*, *Salmonella* spp., and *E. coli* O157:H7 remained detectable after selective enrichment even for more than 200 days of ripening [24–26].

However, no studies have investigated the behavior of foodborne pathogens during the production and ripening of Pecorino Romano cheese. Pecorino Romano is an Italian protected designation of origin (PDO) semi-cooked hard cheese, which must be made exclusively from raw or thermized whole sheep milk, according to the PDO specifications [27]. PDO Pecorino Romano cheese is the most popular ovine cheese produced in Italy and it has a very important role from the economic point of view. Indeed, Pecorino Romano is one of the most exported Italian cheeses in the world [28]. The United States, with an export quota of around 13,000 tons in 2019, is the leading export destination, where this cheese is mainly used as an ingredient in the food industry [29].

Typically Pecorino Romano cheese has about 32% of moisture, a water activity (aw) value of around 0.85, and a pH value between 5.07 and 5.31. This cheese usually has a high salt content that ranges from 4.5% to 8.3%. The minimum ripening period is 150 days for the table cheese, while grating cheese requires 240 days [28].

Although unpasteurized milk cheeses are commonly consumed in a large number of countries, some of them, the United States primarily, have doubts about the safety of these products. For instance, the Food and Drug Administration (FDA) is frequently evaluating whether to propose more restrictive requirements on the sale of unpasteurized milk cheeses, like Pecorino Romano cheese [30,31]. Despite the healthiness of this product, also indirectly attested from the lack of foodborne infection or intoxication episodes bound to the consumption of this cheese, experimental studies are important to investigate the fate of pathogenic bacteria during Pecorino Romano cheese manufacturing and ripening.

The main objective of this work was to investigate the survival of *L. monocytogenes*, *Salmonella* spp., *E. coli* O157:H7, and *S. aureus*, during the ripening of PDO Pecorino Romano cheese made from raw or thermized milk.

#### **2. Materials and Methods**

#### *2.1. Experimental Design*

Experimental productions of Pecorino Romano cheese were carried out according to PDO specifications, using cheese-making facilities of Agris Sardegna (Olmedo, Italy). Although PDO Pecorino Romano cheese is now exclusively produced from thermized milk, PDO specifications do not exclude the use of raw milk. Therefore, we performed two series of experimental cheese-makings, from raw milk (RM) and thermized milk (TM). Each series consisted of twelve cheese batches, three replicates for each selected pathogen. On the same day two experimental cheese-makings were performed starting from a single batch of raw whole sheep milk inoculated with a given pathogenic microorganism, the first one from raw milk and the second one after thermization of the same sheep milk.

#### *2.2. Bacterial Strains and Inoculum Preparation*

Seven strains of *L. monocytogenes*, two of *Salmonella* spp., three of *E. coli* O157:H7, and five of *S. aureus* were used (Table 1). We have chosen both reference and wild strains, the latter were isolated from milk and dairy products. Each inoculum was prepared from a culture containing different strains (reference and wild strains) of the same species.


**Table 1.** Pathogen strains used in Pecorino Romano cheese-making trials.

<sup>a</sup> Strain designation provided by collection. <sup>b</sup> Collection: ATCC, American Type Culture Collection, USA; IZSLER, Istituto Zooprofilattico Sperimentale della Lombardia e dell'Emilia Romagna, Italy; IZS, Istituto Zooprofilattico Sperimentale della Sardegna, Italy. <sup>1</sup> Reference strains. <sup>2</sup> Wild strains isolated from milk and dairy products.

Inoculum preparation included several steps. Briefly, all strains were conserved in Microbank beads (Biolife, Milan, Italy) at −18 ◦C. For each experiment, the strains were previously reactivated in brain heart infusion broth (Biolife, Milan, Italy). The pre-inoculation preparation included serial passages on trypticase soy broth (TSB, Oxoid, Thermo Fisher Scientific, Basingstoke, UK). Two milliliters of each strain culture were inoculated on 35 mL of TSB and incubated overnight at 37 ◦C under continuous stirring. Then, for each species (*L. monocytogenes*, *Salmonella* spp., *E. coli* O157:H7 and *S. aureus*), subcultures have been combined in equal quantity and centrifuged at 6000 rpm for 10 min. Subsequently, the supernatant was discarded and the pellet was resuspended in physiological saline solution (NaCl 0.85%, pH 7). The optical density at 600 nm (Perkin Elmer Lambda 25 UV/VIS Spectrophotometer, Perkin Elmer, Waltham, MA, USA) was determined and counts were confirmed by serial decimal dilution and inoculation in the following selective agar media: Agar Listeria according to Ottaviani and Agosti (ALOA, Biolife, Milan, Italy) plates for *L. monocytogenes*; Baird Parker plates with Rabbit Plasma Fibrinogen supplement (RPF, Biolife, Milan, Italy) plates for *S. aureus*; cefixime tellurite sorbitol MacConkey agar plates (CT-SMAC, Biolife, Milan, Italy) for *E. coli* O157:H7; xylose lysine deoxycholate agar plate (XLD, Microbiol, Uta, Italy) for *Salmonella* spp. The agar plates were incubated at 37 ◦C for 24 h.

#### *2.3. Cheese Production and Sampling*

Each time, 700 L of raw whole bulk sheep milk collected from the "Bonassai" experimental farm of Agris Sardegna, were used. The milk was previously inoculated, under gentle stirring, with a multi-strain bacterial suspension of the same pathogenic microorganism to get a final concentration of approximately 10<sup>6</sup> CFU mL−<sup>1</sup> and then was split into two 350 L aliquots, for RM and TM cheese-makings. In TM trials, milk was subsequently heated to 65 ◦C, without resting at the set temperature and quickly cooled down to coagulation temperature (38 ◦C). Then, the manufacturing process (RM and TM) followed the same procedure reported in Figure 1. Two vats were used alternately, one for RM trial and one for TM trial, in order to eliminate any vat variation. A thermophilic freeze-dried starter culture (FD-DVS CO-02, CHR Hansen, Hoersholm, Denmark) including *Streptococcus thermophilus* and *Lactobacillus delbrueckii* subsp. *bulgaricus* was added at a concentration of 106 CFU mL<sup>−</sup>1.

Three cheese wheels of approximately 25 kg at 1 day were obtained from each cheese-making, for a total of 36 RM and 36 TM cheeses wheels. After molding, cheeses were subjected to drainage in hot room at 36 ◦C until reaching pH 5.20–5.30 and then at 20 ◦C up to 18–24 h. Dry salting and ripening were conducted in controlled conditions (10–12 ◦C and 78–85% relative humidity). In particular, cheeses were dry salted after 48 h (first application) from the manufacture. Later, during the first 90 days of ripening, cheeses were salted three more times, respectively after 12, 26, and 56 days. After 90 days, at the end of the salting, cheeses were washed, dried, and aged for two more months, for a total of 5 months of ripening. The cheese was sampled after 1, 90, and 150 days for physico-chemical analysis while inoculated raw milk, thermized milk, and cheese were sampled to enumerate pathogenic bacteria.

#### *2.4. Physico-Chemical Analysis*

All physico-chemical analysis were performed in duplicate. Samples of curd and cheese were analyzed for pH (pH-meter Basic 20+, Crison Instruments S.A., Alella, Spain). The following parameters were determined for the cheese samples: moisture (ISO 5534:2004) [32]; sodium chloride, determined by potentiometric titration with AgNO3 (ISO 5943:2006) [33] (automatic titrator, model DL55, Mettler-Toledo GmbH, Schwerzenbach, Switzerland); water activity (aw), determined at 25 ◦C (Aw Sprint instrument, Axair Ltd., Novasina Division, Lachen, Switzerland).

#### *2.5. Microbiological Analysis*

Aliquots of 25 mL or g for qualitative detection and 10 mL or g for quantitative detection were taken from each sample and used to prepare a suspension with appropriate diluents. The following parameters have been then researched:

(a) *Listeria monocytogenes.* Detection and enumeration were performed according to ISO 11290-1:2017 and ISO 11290-2:2017, respectively [34,35]. For detection, samples were mixed with Fraser broth base (Oxoid, Thermo Fisher Scientific, Basingstoke, UK), homogenized 90 s in a Stomacher Lab Blender 400 (International PBI S.p.A., Milan, Italy), and incubated for 24 h at 30 ◦C (primary enrichment). Subsequently, 100 μL of primary enrichment were transferred to 10 mL of Fraser broth supplemented by Fraser selective supplement (Oxoid, Thermo Fisher Scientific, Basingstoke, UK), which were incubated for 24 h at 37 ◦C (secondary enrichment). From primary and secondary enrichments, aliquots of 100 μL were streaked onto selective differential medium plates: Agar Listeria according to Ottaviani and Agosti (ALOA, Biolife, Milan, Italy) and Listeria selective agar (Oxford formulation, Oxoid, Thermo Fisher Scientific, Basingstoke, UK) and incubated for up to 24–48 h at 37 ◦C. All isolates with typical *L. monocytogenes* characteristics were subjected to morphological and biochemical proofs as confirmatory tests.

For the enumeration of *L. monocytogenes* ten-fold serial dilutions were made and aliquots of 100 μL were plated in duplicate on the surface of ALOA agar plates which were incubated for 48 h at 37 ◦C. Presumptive *L. monocytogenes* colonies were counted after confirmatory tests.

	- (1) Enrichment of the test portion homogenized in modified tryptone soya broth containing novobiocin (mTSB + N, Biolife, Milan, Italy) with incubation at 41.5 ◦C for 6 h and subsequently for a further 12 h to 18 h.
	- (2) Separation and concentration of microorganisms by means of immunomagnetic particles coated with antibodies to *E. coli* O157:H7.
	- (3) Isolation by subculture of the immunomagnetic particles with adhering bacteria onto cefixime tellurite sorbitol MacConkey agar (CT-SMAC, Biolife, Milan, Italy) and sorbitol MacConkey agar (SMAC, Biolife, Milan, Italy) incubated at 37 ◦C for 24 h.
	- (4) Confirmation of typical colonies.

*E. coli* O157:H7 count was determined by ten-fold serial dilutions and direct plating (100 μL in duplicate) on CT-SMAC agar plates incubated at 37 ◦C for 24 h. Typical *E. coli* O157:H7 colonies were subjected to confirmatory tests and were then enumerated.

(d) *Salmonella* spp. The detection was performed according to ISO method 6579-1:2017 [38]. The method required the following successive stages. A pre-enrichment in buffered peptone water (BPW, Microbiol, Uta, Italy) at 37 ◦C for 24 h. A selective enrichment in Rappaport-Vassiliadis with soy broth (RVS, Oxoid, Basingstoke, UK) and Müller-Kauffmann tetrathionate-novobiocin broth (MKTTn, Microbiol, Uta, Italy) for 24 h at 41.5 and 37 ◦C, respectively. Aliquots of the selective broths were streaked onto two selective isolation agar media, xylose lysine deoxycholate agar (XLD, Microbiol, Uta, Italy) and Salmonella detection and identification agar (SMID, BioMérieux, Marcy L'Etoile, France). The agar plates were incubated at 37 ◦C for 24 h. Confirmation of suspect colonies was carried out by biochemical and serological testing.

For the enumeration of *Salmonella* spp. ten-fold serial dilutions were prepared and aliquots of 100 μL were double plated on XLD agar plates which were incubated at 37 ◦C for 24 h. Presumptive *Salmonella* spp. colonies were subjected to confirmatory tests and were then counted.

**Figure 1.** Experimental flow diagram of protected designation of origin (PDO) Pecorino Romano cheese production. RH = relative humidity.

#### *2.6. Statistical Analysis*

Microbial loads were log transformed and expressed as means and standard deviation (SD). The means and SD of physico-chemical parameters were determined from twelve cheese-makings for each experimental series (RM and TM). Analysis of variance was carried out using Minitab statistical package release 16 (Minitab Inc., State College, PA, USA). The general linear model procedure was used to verify the effects of the two studied factors "milk heat treatment" (2 levels) and "ripening time" (3 levels) on the microbial loads and physico-chemical parameters, as far as their interaction is concerned. The comparison between means was performed using Tukey's significant difference test (*p* < 0.05). Data were also analyzed by the Pearson correlation to measure the degree of the linear relationship.

#### **3. Results and Discussion**

#### *3.1. Physico-Chemical Properties of Cheese*

Figure 2 and Table 2 show the pH values in RM and TM cheese-making trials at different time. The pH values differed significantly between RM and TM (*p* < 0.05) in various step of cheese manufacturing process. As shown in Figure 2, after molding was completed, RM curd had a higher pH than TM curd (6.4 vs. 6.2, respectively; *p* < 0.05). This difference suggests a slight delay in the starter culture acidification activity in the early stage of RM cheese-making. This time lag could be due to the bacteriostatic activity of milk proteins such as lactoferrin, lactoperoxidase, and immunoglobulins since it normally decreases with milk heat treatments [39]. We cannot exclude that the competition of the starter culture with the native microflora of raw milk may also be involved in this delay of the acidification process. During drainage in hot room at 36 ◦C, a difference in pH values between RM and TM cheese-making was kept, however, this gap became less evident after two hours of drainage. Moreover, the lower SD in TM acidification curve shown in Figure 2 suggests a more regular and repeatable trend of the acidification process during TM cheese-making compared to RM one.

#### **Sampling point**

**Figure 2.** Acidification profiles of Pecorino Romano cheese produced from raw (RM), and thermized milk (TM). Twelve replicates for each cheese-making technology. Error bars indicate standard deviations. Different letters at the same time point indicate significant differences (*p* < 0.05).


**Table 2.** Physico-chemical parameters of Pecorino Romano cheese at different ripening times, obtained from raw (RM) and thermized milk (TM). Values are given as means ± standard deviation.

Values in the same line with different superscript letters differ significantly (Tukey's test, *p* < 0.05). DM = dry matter; Aw = water activity; RM = cheese from raw milk; TM = cheese from thermized milk; T = thermal treatment; R = ripening time. NS, *p* > 0.05; \* *p* < 0.05; \*\* *p* < 0.01; \*\*\* *p* < 0.001.

As shown in Table 2, the pH values in RM and TM cheeses were significantly affected by the milk thermization treatment (*p* < 0.01) and ripening time (*p* < 0.01), furthermore the interaction between these factors was significant (*p* < 0.01). Despite the slight delay in the acidification process, after 1 day, RM and TM cheeses were not statistically different for pH values (5.3 and 5.21, respectively). On the contrary, after 90 and 150 days of ripening, RM and TM cheeses differed significantly for pH values (*p* < 0.01). In particular, pH in RM cheeses remained substantially unchanged during ripening (around 5.3), while pH in TM cheeses showed a slight decrease, reaching an average pH value of 5.0 at 150 days. The evolution of pH during cheese ripening is due to several factors involved in the metabolism of the main components of the cheese: hydrolysis of the residual lactose, which involves a decrease in pH, and degradation of proteins that, instead, leads to an increase in pH value [40,41]. Pecorino Romano cheese is characterized by modest proteolysis due to the intense salting that limits the activity of the proteolytic enzymes. On the other hands, lipolysis is more accentuated mainly because of the use of exogenous lipases contained in the lamb paste rennet, used for milk coagulation. These enzymes catalyze the biochemical process of triglyceride hydrolysis resulting in the release of short-chain fatty acids that can contribute to a decrease in the pH during ripening [28,42]. Finally, it is important to point out that the pH values found both for RM and TM cheeses are in the range of variation found in commercial PDO Pecorino Romano cheese [28].

In Table 2 we report also the remaining physico-chemical properties of RM and TM cheese, at different time of ripening. The milk thermization treatment did not significantly affect the moisture content, which instead was significantly influenced by ripening time and its interaction with heat treatment (respectively, *p* < 0.001 and *p* < 0.05). No significant differences were found between RM and TM cheese at each observed ripening time. The major changes in moisture occur in the first 90 days of ripening in both TM and RM cheeses, while the decrease from 90 up to 150 days was statistically significant only in TM cheese (*p* < 0.05). Moisture content was approximately 42% in cheese after 1 day and 31–32% in cheese after 90 and 150 days of ripening, according to Addis et al. [28,42].

As reported for moisture content, also the aw values were significantly affected by the ripening time and its interaction with heat treatment (*p* < 0.001) (Table 2). The aw, significantly differed in RM and TM 1 day cheeses, and dropped significantly from approximately 0.98–0.97 (TM and RM, respectively) up to 0.89–0.88 at the end of the salting process (90 days), then slightly decreased tendentially (*p* > 0.05) at 150 days. The aw values at the end of ripening (150 days) were higher than those reported by Addis et al. [28]. However, these authors referred to longer aged Pecorino Romano cheese (7–8 months).

At the end of the salting process (90 days), the NaCl on dry matter (DM) content was significantly lower in RM cheeses (6.3%) compared to TM cheeses (7.2%). Salt content slightly increased, up to the end of ripening (6.4% and 7.3%, respectively in RM and TM cheeses). This difference in NaCl content between RM and TM cheeses could be due to the intrinsic variability of the dry salting technique. However, NaCl values are comparable to those of commercial Pecorino Romano cheese, which is characterized by high and variable salt content [28].

These results showed that the experimental cheeses met the typical requirements of PDO Pecorino Romano cheese. Moreover, the different cheese-making technology (RM and TM) seems to affect only in part the characteristics of the product at the end of the ripening period, which substantially has the same peculiarities of the Pecorino Romano cheese available on the market.

#### *3.2. Pathogenic Bacteria Counts in RM Cheese-Making Trials*

The pathogens counts in raw milk are given in Table 3. Pathogenic bacteria levels in raw milk pointed out the accuracy of the technique used for the inoculum preparation. Indeed, all the counts were around 6 log10 CFU mL−1, as established from the experimental protocol. Table 4 shows the pathogens counts in cheese at different ripening stages. A reduction in bacterial loads was observed for all pathogens in 1-day old cheese. *L. monocytogenes* and *Salmonella* spp. showed the major reduction, with a decrease of about 4 log10. The number of viable *E. coli* O157:H7 decreased by about 2.5 log10, while *S. aureus* counts showed a reduction by about 0.5 log10. The behavior of the pathogenic bacteria was likely in part attributed to the number of viable cells lost with the whey separation and partly to the combined effect of various hurdles occurring during the first phase of Pecorino Romano cheese manufacturing, such as competition with the starter culture, acidification, and curd cooking [17,43]. Conversely, other studies showed an increase in pathogens counts during the early stages of cheese-making in different unpasteurized milk cheeses. Bellio et al. [21] reported a growth of 2 log cycle in *E. coli* O157:H7 levels during the first 24 h of PDO Fontina cheese production, despite a curd cooking phase at 48 ◦C. Similar findings were also given by other authors in Cacioricotta, Cheddar, and Gouda cheeses [9,24]. Chatelard-Chauvin et al. [26] found that *L. monocytogenes* increased by about 3.5 log10 in Cantal cheese at 24 h. An increase in *Salmonella* spp. populations over 1 log10 was observed in Gouda and Cheddar cheeses in the first 24 h of manufacturing [25,44]. In our study, *S. aureus* fell slightly, however, some authors observed a growth in *S. aureus* counts in the early stages of cheese production [45,46]. These authors considered that the increase in pathogenic bacteria counts during the early stages of cheese-making could be due both to an actual bacterial growth and to the entrapment of pathogens cells within the curd, during curd contraction and whey separation. In our study, the reduction in pathogens population levels during the early stages of cheese manufacturing might be related not only to the aforementioned hurdles to bacteria growth and survival but also to a protracted curd breaking (about 6 min) into extremely fine grains (3–5 mm). This is a typical feature of Pecorino Romano cheese-making technology which could lead to a higher level of pathogens viable cells lost in the cheese whey.

All tested pathogens, despite a reduction, were still present in 1-day old cheese in relevant bacterial loads. However, it must be pointed out that in the present study, high inoculum levels were used, which are difficult to find in practice, to simulate the worst-case scenario of contamination. Overall, our results are in agreement with other studies that indicate the ability of the pathogenic bacteria to overcome the obstacles of the first phase of cheese-making, first of all, acidification, probably owing to an acid tolerance response (ATR) widely reported in the literature [18–20].

**Table 3.** Results of pathogens counts (log10 CFU mL<sup>−</sup>1) in sheep milk before and after thermization. Values are given as means ± standard deviation. *Escherichia coli* O157:H7: ATCC 43984, 47, 719. *Salmonella spp*.: Typhimurium ATCC 6994, Enteritidis 670. *Listeria monocytogenes*: ATCC 15313, ATCC 19114, ATCC 9525, ATCC 153/3, 2, 90, V7. *Staphylococcus aureus*: ATCC 14458, ATCC 25923, 401, 466, 64494.


Values in the same line with different superscript letters differ significantly (Tukey's test, *p* < 0.05). \*\*\* *p* < 0.001.

**Table 4.** Results of pathogens counts (log10 CFU g<sup>−</sup>1) in Pecorino Romano cheese at different ripening times, obtained from raw (RM), and thermized milk (TM). Values are given as means ± standard deviation. *Escherichia coli* O157:H7: ATCC 43984, 47, 719. *Salmonella spp*.: Typhimurium ATCC 6994, Enteritidis 670. *Listeria monocytogenes*: ATCC 15313, ATCC 19114, ATCC 9525, ATCC 153/3, 2, 90, V7. *Staphylococcus aureus*: ATCC 14458, ATCC 25923, 401, 466, 64494.


Values with different superscript letters differ significantly (Tukey's test, *p* < 0.05). RM = cheese from raw milk; TM = cheese from thermized milk; T = thermal treatment; R = ripening time. \*\*\* *p* < 0.001.

Table 4 also shows the pathogens counts at the end of salting (90 days) and the end of the ripening period (150 days). At the experimental conditions, all inoculated pathogens except *S. aureus*, were not detectable after 90 days, even after selective enrichment. This result was confirmed in 150-days old cheese. *S. aureus* counts were below the detection limit for direct plating enumeration method (<1 log10 CFU g−1), both in cheese at 90 and 150 days of ripening. In the early phase of cheese-making (24 h), the pathogenic bacteria were subjected to initial stress because of the curd cooking and fast acidification. Afterwards, in the first 90 days of ripening, a significant correlation was found between microbial loads and a*w*, moisture and NaCl content. In particular, all the studied pathogens exhibited a strong positive correlation with moisture (E. coli O157:H7, *R*<sup>2</sup> = 0.987, *p* < 0.001; L. monocytogenes, *R*<sup>2</sup> = 0.982, *p* < 0.001; *Salmonella* spp., *R*<sup>2</sup> = 0.987, *p* < 0.001; *S. aureus*, *R*<sup>2</sup> = 0.993, *p* < 0.001) and aw (*E. coli* O157:H7, *R*<sup>2</sup> = 0.971, *p* < 0.001; L. monocytogenes, *R*<sup>2</sup> = 0.982, *p* < 0.001; *Salmonella* spp., *R*<sup>2</sup> = 0.971, *p* < 0.001; *S. aureus*, *R*<sup>2</sup> = 0.948, *p* < 0.001). Conversely, a strong negative correlation was found between microbial count and salt content (*E. coli* O157:H7, *R*<sup>2</sup> = <sup>−</sup>0.998, *<sup>p</sup>* <sup>&</sup>lt; 0.001; L. monocytogenes, *R*<sup>2</sup> = <sup>−</sup>0.997, *<sup>p</sup>* <sup>&</sup>lt; 0.001; *Salmonella* spp., *R*<sup>2</sup> = <sup>−</sup>0.987, *p* < 0.001; *S. aureus*, *<sup>R</sup>*<sup>2</sup> <sup>=</sup> <sup>−</sup>0.996, *<sup>p</sup>* <sup>&</sup>lt; 0.001).

Several studies show the ability of pathogenic bacteria to survive in cheeses even beyond 90 days of ripening. Nevertheless, none of these had similar physico-chemical properties to those of PDO Pecorino Romano cheese (pH ≤ 5.4; moisture ~32%; aw < 0.90; NaCl 4–8%), together with a long ripening period. Ioanna et al. [9] reported that *E. coli* O157:H7, when inoculated in milk at 2 log10 CFU mL<sup>−</sup>1, survived in 90 days-old Cacioricotta cheese with a bacterial load greater than 4 log10 CFU g−1. Similar findings were observed in Fontina cheese at 80 days of ripening [21]. In Cheddar and Gouda cheeses, *E. coli* O157:H7 was detectable after an enrichment procedure for up to around 300 days, with an initial raw milk inoculation of about 20 CFU mL−1, which was far lower than that used in the present study [24]. *E. coli* O157:H7 appears to be one of the pathogens having the greatest ability to overcome environmental stresses and it is also infectious at very low doses (5–50 viable cells) [47]. Therefore, proving its complete inactivation is relevant. Viable cells of *Salmonella* spp. were detectable in Gouda cheese for more than 200 days [25]. However, in Feta and Tulum cheese, with a salt content close to that of Pecorino Romano cheese, *Salmonella* spp. was no more found after 20 and 90 days, respectively [48,49]. *L. monocytogenes* was detectable for more than 250 days in Cantal cheese [26]. In agreement with our study, Wusimanjiang et al. [50] observed that *L. monocytogenes*, when inoculated in milk at 3 log10 CFU mL<sup>−</sup>1, reached undetectable levels after 28–35 days of ripening in a cheese with a salt content similar to that of Pecorino Romano cheese (4–10% on DM basis).

According to Regulation (EC) No 2073/2005, Pecorino Romano cheese, already at 90 days, belongs to ready-to-eat foods that do not support the growth of *L. monocytogenes*, since its aw value is less than 0.92. For these products, the Regulation establishes a tolerance for *L. monocytogenes* up

to 100 CFU g−<sup>1</sup> [51]. The cheeses we obtained were *Listeria*-free at 90 days, despite the high initial inoculum in cheese milk.

*S. aureus*, unlike the other pathogens, was still present in high population levels in 1-day old cheese. Since *S. aureus* is halotolerant, it could find the proper environment to survive or even grow during or after the salting process, as some authors have already reported [52,53]. However, according to the growth limits of *S. aureus* stated by Valero et al. [54], the conditions of ripening and the characteristics of Pecorino Romano cheese: ripening temperature (10–12 ◦C), constant low pH (~5.3), and decreasing of aw from 0.97 to 0.89 in the first 90 days of ripening, were enough severe to prevent proliferation and also to reduce the survival of *S. aureus* as shown by the results (<1 log10 CFU g<sup>−</sup>1). These findings were in agreement with those observed by Pexara et al. [45] in Feta cheese with a salt content close to that of Pecorino Romano cheese.

Unlike the other bacteria tested in this study, the pathogenicity of *S. aureus* is due to the possible production of different enterotoxins in contaminated food. A level above 105 CFU mL−<sup>1</sup> or g−<sup>1</sup> of *S. aureus* can produce an infective enterotoxins concentration of about 1 μg [55]. Therefore, Regulation (EC) No 2073/2005 sets the limit in cheeses for coagulase-positive staphylococci (including *S. aureus*) at 10<sup>5</sup> CFU g−<sup>1</sup> [51]. Above this limit, the presence of staphylococcal enterotoxins is expected, and there is a need to perform assays for enterotoxins detection, which must be absent in 25 g. However, staphylococcal enterotoxins production seems to be strongly limited in cheeses by the antagonistic effect of lactic acid bacteria. Microbial competition together with unfavorable environmental conditions, such as low pH, limits the growth of *S. aureus* and strongly downregulate virulence genes expression [52,53,56]. In particular, a rapid acidification in the first 6 h of cheese-making, as reported in our study (Figure 2), seems to be critical for *S. aureus* growth and enterotoxins production [46]. Furthermore, in our work, we simulated a worst-case contamination scenarios by employing unrealistically high starting *S. aureus* levels and despite this, the cheese-making conditions prevented the proliferation of *S. aureus*. Future studies will be carried out to examine the presence of staphylococcal enterotoxins in Pecorino Romano cheese produced from milk contaminated with *S. aureus*.

#### *3.3. Pathogenic Bacteria Counts in TM Cheese-Making Trials*

In Table 3 are reported the pathogens counts in milk (before and after thermization). Thermization is a sub-pasteurization heat treatment usually performed at 57–68 ◦C for 10–20 s to reduce the number of spoilage bacteria, with minimum collateral heat damage to milk components and thus provide a suitable environment for the growth of the added starter cultures [57]. Unlike the pasteurization, the thermization process is not able to ensure the complete inactivation of pathogenic microorganisms but, especially if they are present in high initial levels, it can only reduce their number, as confirmed by the pathogens counts in thermized milk shown in Table 3. In our study, we used a batch thermization at 65 ◦C, without resting at the set temperature. All the pathogens showed a reduction of around 3 log10 after thermization, starting from a bacterial count of about 6 log10 CFU mL−<sup>1</sup> (Table 3). These rates of thermal inactivation are close to those reported by other authors, although a comparison with the time-temperature conditions reported in our work is not easy [58,59]. These studies highlight the heterogeneity in thermal resistance not only in different microbial groups but also between different strains of the same species. In particular, different strains of *E. coli* and *S. aureus* showed varying degrees of heat resistance, while *L. monocytogenes* and *Salmonella* spp. strains displayed a more uniform thermal susceptibility [59]. Hence, it is important to underline that a multi-strain inoculum should be used to represent the behavior of a given pathogen, as done in our work. The data reported in Table 3 show that counts levels in thermized milk did not differ significantly among the four pathogenic microorganisms, which indicate similar thermal tolerance under the conditions applied in the present study. However, we cannot exclude that some of the strains used in our study were more susceptible toward heat than others.

After 1 day from the production, *L. monocytogenes*, *E. coli* O157:H7, and *Salmonella* spp. dropped below the detection limit and were no longer detectable up to 150 days, even after an enrichment procedure, as shown in Table 4. *S. aureus* counts were also below the detection limit for the direct plating enumeration method (<1 log10 CFU g−1) from 1 day up to the end of the ripening period. As already reported for RM trials, it was not possible to assess exactly the amount of the pathogens decrease related to the bacteria loss through the whey separation compared to the pathogens counts reduction due to the several environmental hurdles during the first stage of cheese-making. The latter may have a greater impact on TM trials, especially as a consequence of milk thermization, a second less severe thermal stress (curd cooking), and the faster and more intense acidification. In fact, after molding, the curd had a lower average pH value in TM trials (6.24) compared to RM trials (6.41). This difference was kept during the drainage process (Figure 2). Besides, as reported in Table 4, the microbial counts in RM and TM cheeses were significantly affected by milk thermization treatment (*p* < 0.001), as well as ripening time (*p* < 0.001). Furthermore, the interaction between these factors was also significant for pathogens counts (*p* < 0.001).

The introduction of the thermization in Pecorino Romano cheese-making is a further obstacle to the possible growth and survival of pathogens. In addition to other hurdles (curd cooking and acidification), it makes the cheese microbiologically safe before the salting process. In the case of RM cheeses, it is instead necessary to wait until the end of the salting process (90 days) to consider safety of the product. In any case, this happens 60 days before the mandatory minimum period for marketing PDO Pecorino Romano cheese (150 days).

Hence, our study shows that the production process of Pecorino Romano cheese both from raw and thermized milk is enough restrictive for the survival of the pathogens we have examined, except for *S. aureus* which may be still present in such low levels that it does not represent a hazard to the consumers. Currently, thermization is practiced by all producers of PDO Pecorino Romano cheese since this treatment has the advantage in standardizing the microbiological properties of the cheese milk by reducing the total microbial count and containing the undesirable bacteria, to obtain a constant quality of the product. This is a key element for a product aimed at an international market, such as Pecorino Romano cheese.

#### **4. Conclusions**

Cheeses made from unpasteurized milk are known to pose a health risk to the consumers. According to PDO specifications, Pecorino Romano cheese can be produced from raw or thermized whole sheep milk. In this study, Pecorino Romano cheese was produced from unpasteurized (raw and thermized) milk inoculated with *L. monocytogenes*, *Salmonella* spp., *E. coli* O157:H7, or *S. aureus*, in high initial loads (106 CFU mL−1). The results demonstrated that the cheese-making process may ensure the microbiological safety of this cheese after the minimum ripening period (150 days) established by the PDO specifications. In particular, the obtained cheese was free from microbiological hazard after 90 days and 1 day, when manufactured from raw and thermized milk, respectively. These results could encourage cheese manufacturers to diversify Pecorino Romano cheese by using raw milk, which is currently almost not used, for its production, to provide the market with a product having a more intense flavor than that made from thermized milk. This is the first preliminary work to evaluate the survival of selected pathogenic bacteria during Pecorino Romano cheese ripening. Further studies are necessary to investigate the behavior of the pathogens in the early stages of manufacturing in order to understand more in depth on how their survival is hampered by different factors. The eventual production of staphylococcal enterotoxins will also be tested in subsequent investigations.

**Author Contributions:** Conceptualization, A.P. and A.F.; methodology, A.P., M.A. and A.F.; validation, G.L., R.M., M.P., M.A., A.F., and A.P.; formal analysis, G.L., R.M., and M.P.; investigation, G.L., R.M., and M.P.; resources, A.P., M.A., and A.F.; writing—original draft preparation, G.L., R.M., M.A., and M.P.; writing—review and editing, A.P. and A.F.; supervision A.P. and A.F. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

*Dairy* **2020**, *1*

**Acknowledgments:** The technical staff of both the Technology and Chemistry sectors of Agris Sardegna is gratefully acknowledged.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

### *Article* **Functional Odd- and Branched-Chain Fatty Acid in Sheep and Goat Milk and Cheeses**

**Anna Nudda, Fabio Correddu \*, Alberto Cesarani, Giuseppe Pulina and Gianni Battacone**

Dipartimento di Agraria, Sezione di Scienze Zootecniche, University of Sassari, viale Italia, 39, 07100 Sassari, Italy; anudda@uniss.it (A.N.); acesarani@uniss.it (A.C.); gpulina@uniss.it (G.P.); battacon@uniss.it (G.B.) **\*** Correspondence: fcorreddu@uniss.it; Tel.: +39-079-229-308

**Abstract:** The inverse association between the groups of odd-chain (OCFA) and branched-chain (BCFA) and the development of diseases in humans have generated interest in the scientific community. In experiment 1, the extent of the passage of odd- and branched-chain fatty acids (OBCFA) from milk fat to fresh cheese fat was studied in sheep and goats. Milk collected in two milk processing plants in west Sardinia (Italy) was sampled every 2 weeks during spring (March, April and May). In addition, a survey was carried out to evaluate the seasonal variation of the OBCFA concentrations in sheep and goats' cheeses during all lactation period from January to June. Furthermore, to assess the main differences among the sheep and goat cheese, principal component analysis (PCA) was applied to cheese fatty acids (FA) profile. Concentrations of OBCFA in fresh cheese fat of both species were strongly related to the FA content in the unprocessed raw milk. The average contents of OBCFA were 4.12 and 4.13 mg/100 mg of FA in sheep milk and cheese, respectively, and 3.12 and 3.17 mg/100 mg of FA in goat milk and cheese, respectively. The OBCFA concentration did no differed between milk and cheese in any species. The content of OBCFA was significantly higher in sheep than goats' dairy products. The OBCFA composition of the cheese was markedly affected by the period of sampling in both species: odd and branched FA concentrations increased from March to June. The seasonal changes of OBCFA in dairy products were likely connected to variations in the quality of the diet. The PCA confirmed the higher nutritional quality of sheep cheese for beneficial FA, including OBCFA compared to the goat one, and the importance of the period of sampling in the definition of the fatty acids profile.

**Keywords:** sheep and goat milk; cheese; odd and branched chain fatty acids

#### **1. Introduction**

The fatty acid (FA) composition of dairy products has assumed considerable interest in consumers from a nutritional and healthy point of view. Indeed, specific FA of dairy products can affect human health and can have an important role in the prevention of metabolic diseases. This increasing attention is also demonstrated by recent studies that investigated the feasibility of improving milk FA profile in sheep and goats though breeding schemes [1–3]. A significant amount of attention has been focused in the last decade on polyunsaturated fatty acid of the omega3 family (PUFA n-3) and conjugated linoleic acid (CLA) concentrations of dairy products. Moreover, the groups of odd-chain fatty acids (OCFA) and branched-chain fatty acids (BCFA), long neglected due to their low incidence on the total amount of FA, have sparked interest in the scientific community due to an inverse relationship with the development of human diseases.

In particular, the group of BCFA comprises mainly saturated fatty acids characterized by the presence of one or more methyl groups in the iso or anteiso position. Such molecules represent the main FAs in some microorganisms (e.g., Bacilli and Lactobacilli), and they can be also observed in mammal tissues. In laboratory animals these FA evidenced antiinflammatory properties [4], reduced the incidence of necrotizing enterocolitis and altered the ecology of the gastrointestinal microorganisms in a neonatal rat model [5].

**Citation:** Nudda, A.; Correddu, F.; Cesarani, A.; Pulina, G.; Battacone, G. Functional Odd- and Branched-Chain Fatty Acid in Sheep and Goat Milk and Cheeses. *Dairy* **2021**, *2*, 79–89. https://doi.org/10.3390/dairy2010008

Received: 9 December 2020 Accepted: 27 January 2021 Published: 5 February 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

The OCFA, which include pentadecanoic (C15:0) and heptadecanoic (C17:0) acids and their isomeric forms, are mostly derived from ruminal bacteria cell wall and then partitioned firstly to organs and tissues and, therefore, detectable in ruminant-derived foods. Thus, dairy foods or ruminant fats represent the major dietary source of OCFA for humans.

Moreover, the milk concentrations of C15:0 and C17:0 are considered biomarkers of rumen microbial fermentation and microbial de novo lipogenesis [6]. In addition, the mammary gland plays a role in their synthesis by using propionate [7], and the subcutaneous adipose tissue and, therefore, after their mobilization, can be incorporated into milk fat [8].

Humans are not able to synthetize C15:0 and C17:0 so they could be considered valuable markers of dairy fat intake [9,10]. The role of these FA in human health has recently been reinforced in view of new biological and nutritional observations. The C15:0 and C17:0 have been inversely associated with cardiovascular disease [11–14] and incidence of type 2 diabetes [15–18]. The direct role of C15:0 in attenuating pro-inflammatory state, cytotoxicity in human cell lines, anemia, and dyslipidemia lowering glucose and cholesterol in vivo has been recently evidenced [19].

Among the different factors influencing the OBCFA concentration in milk and dairy products, animal diet is the main one. Indeed, the diet of animals can modulate microbial growth in the rumen, de novo microbial FA synthesis and the uptake of blood circulating FA by the mammary gland. Variations in milk OBCFA have been observed in response to lipid supplementation [20–22], change in forage to concentrate ratio [23] and use of by-products rich in bioactive compounds [24]. Changes in the OBCFA concentration in milk has been also observed with lactation stage [25] and with energy balance of animals [26], probably as a consequence of the fat mobilization from adipose tissue and the extent of mammary uptake for milk fat synthesis [27].

In Mediterranean areas, almost all sheep milk is processed into cheese, and milk fatty acid profile has important effects on cheese fat quality under nutritional point of view. The total transfer from milk into cheese of fatty acids of nutritional interest, as PUFA n-3 and CLA, has been reported in sheep [28], but effects of cheese-making technology has been also evidenced [29]. However, the mechanisms influencing the relationship between OBCFA concentrations in unprocessed milk and the derived dairy products are still unclear. The presence of FA in the FA composition of different microorganism can be of significant importance for the dairy industry where different microbial cultures, as starters, are widely used for making different products. The consequence could be the change of the native OCFA and BCFA milk concentrations after processing in cheese or other dairy products.

The aims of this study were (a) to evaluate the extent of OBCFA transfer from ovine and caprine milk to cheese (Experiment 1); (b) to describe the seasonal variation of OBCFA in sheep and goat cheeses (Experiment 2). Both were sequentially investigated with two independent surveys.

#### **2. Materials and Methods**

In Experiment 1, bulk tank milk of Sarda breed sheep and Sardinian goats and cheeses were collected from processing plants in North Sardinia for sheep (n. 2) and in west Sardinia for goats (n. 2). Samples were taken every two weeks during spring (March, April and May) for a total of 12 milk samples per species. The transfer of OBCFA from milk to fresh cheeses after 24 h from processing was evaluated.

In Experiment 2, mild sheep and goat cheeses, with a ripening time of about 20–30 days, were sampled monthly from January to July from nine sheep cheese-making plants located in different areas of Sardinia (Guspini, Nurri, Marrubiu, Macomer, Dorgali, Oliena, Onifai, Thiesi, Berchidda), five plants of which produced also goat cheese. The samples collection started in January, which corresponds to the begin of lactation both in sheep and goats reared in Sardinia's breeding system. The goats cheese samples were from four lots/month from January to April and three different lots/month from May to July for a total of 25 goat cheese samples. The sheep cheese samples were Pecorino Sardo PDO type from

8 lots/month from January to May, and from four different lots/month in June and July for a total of 48 sheep cheese samples.

The determination of FA concentrations in milk and cheese samples was carried out by a gas chromatographic method, as previously described [20]. Briefly, after the lipid extraction, a base-catalyzed trans-esterification FIL-IDF standard procedure [30] was used to prepare FA methyl esters (FAME). Identification of individual FAME was allowed by comparing their retention time with that of a series of analytical standards. Those included the Supelco 37 component FAME MIX (Supelco, Bellefonte, PA, USA), the GLC-110 MIX (Matreya Inc. Pleasant Gap, PA, USA) and some individual BCFA (Matreya Inc. Pleasant Gap, PA, USA). Identification of OBCFA was also supported by the consultation of previous studies [31,32]. The FA concentration was reported as mg/100 mg of total FAME.

Using R software (R Core Team, 2020), general linear model (GLM) was applied to data on OBCFA and others main nutritional FA (R Core Team, 2020). The considered fixed effects were species (sheep and goat) and product type (milk and cheese) for experiment 1, whereas species and month for Experiment 2. Significative differences were declared at *p* < 0.05.

Principal component analysis (PCA) was applied on sixty-three FA profile of sheep and goats' cheeses using *prcomp* R function. A total of 42 FA was considered and, therefore, 42 principal components were extracted. Before PCA, all considered FA were scaled to ensure unit variance.

#### **3. Results and Discussion**

Means (±SE) of fat concentration in milk and cheese of sheep and goat, respectively were 6.16% (±0.26) and 5.42% (±0.03), 28.81% (±0.65) and 21.17% (±0.54).

#### *3.1. Experiment 1: Fatty Acid Transfer from Milk to Cheese*

FA composition in milk and cheese for sheep and goats is reported in Table 1. In total, 10 OBCFA were identified in both types of milk (Table 1), including 3 OCFA (C13:0, C15:0 and C17:0), 4 isoBCFA (isoC14:0, isoC15:0, isoC16:0, and isoC17:0) and 3 anteisoBCFA (anteisoC13:0, anteisoC15:0 and anteisoC17:0).

Among them, C15:0 and C17:0 were the most abundant fatty acids, in agreement with other research on sheep [33], goats and cows [34,35]; these OCFA accounted for 29% and 16% of the total concentration of OBCFA, respectively. In our study, C13:0 represents only 2.3%, value that agrees with the observations in dairy cows in which it represents 2% of OCFA [34] or it has not been found [36], but lower than values previously reported in sheep's milk, where accounted for 6.3% of total OCFA [33].

Concentrations of individual OBCFA in the fat of fresh cheeses did not differ from that of raw unprocessed milk (*p* > 0.05) in both species: it means that overall, OBCFA content of milk is not significantly modified by cheese making processes, both in sheep and goats.

In sheep milk, the content of total OBCFA is higher than that observed in goat milk, due to a higher content of isoBCFA form of C15:0, C16:0 and C17:0, and anteiso C15:0 and C17:0 (*p* < 0.05). No differences among species have been observed only for isoC13:0 and isoC14:0. The differences between species could be related to their different feeding regimen, as OBCFA in milk may change according to the bacterial population present in the rumen and, particularly, because of the differences in their relative abundance. Differences in OBCFA concentrations between these the species can occur also because of different extents of de novo synthesis of these FA in the mammary gland, differential efficiency of intestinal absorption, and differential storage of some OBCFA in adipose tissue, followed by their release during periods of fat mobilization.

The C17:1c9 was not detected in sheep and it was found in negligible amount in goat milk and fresh cheese. In the mammary gland, the delta-9 desaturase activity can promote the metabolization of some odd-chain fatty acids; only the conversion of C17:0 to C17:1cis9 has been found of quantitative importance in cows [37]. It has also been reported with a mean concentration of 0.2% in milk of Assaf ewes raised in controlled and intensive farming system [33]. Therefore, the lack of C17:1c9 in Sarda dairy ewes might be related to specific environmental condition related to the typical extensive breeding system of Sarda dairy ewes and Sardinian goats.

The total content of OBCFA in sheep milk is higher in this experiment than that reported for cows [36] and sheep (3.19 g/100 g of total FA) [33]. Differences could be related to different feeding strategies at farm level, as grazing pasture, or vegetable lipid supplementation, forage to concentrate ratio or other dietary factors, that deserve to be further investigated.

**Table 1.** Odd and branched chain fatty acid in milk and cheese of sheep and goat species.


<sup>1</sup> SEM = standard error of means; <sup>2</sup> <sup>M</sup>*vs*C = Milk vs. Cheese; <sup>3</sup> Fatty acids: expressed as mg/100 mg of FAME; \*\* *<sup>p</sup>* ≤ 0.01, \* *<sup>p</sup>* < 0.05, ns not significant *p* ≥ 0.05.

#### *3.2. Experiment 2: Seasonal Variation of Branched and Odd Chain FA*

Seasonal variation of individual and group of branched and odd chain fatty acid has been investigated by collecting cheese samples throughout all lactation period of sheep and goats, that in Sardinia, as in all Mediterranean environmental conditions, occurs from January (start of lactation) to July (when animals are dry off).

The temporal evolution of these FA in cheese showed an increase of total OCFA and BCFA in both species (*p* < 0.05) across all the sampling period (Figure 1a,b); the highest contents of these FA were observed in the cheese produced with milk in advanced lactation, which corresponds also to the hot season. The increase of OBCFA across lactation has been previously reported in dairy cows [34].

Among individual OCFA, the content of odd C15:0 during lactation followed the same pattern of the total OCFA, being the main representative of this group of FA (Figure 2a), whereas C17:0 showed constant concentration during the early and mid-lactation, then increased at the end of lactation, in June and July for both species (Figure 2b). As far as the changes of individual BCFA during lactation are concerned, both anteisoC15:0 and anteisoC17:0 (Figure 2c,d) followed a similar pattern observed for total BCFA.

**Figure 1.** Temporal evolution of (**a**) odd chain fatty acids (FA) and (**b**) branched chain FA from January to July in cheeses from sheep (light grey line) and goats (dark grey line) milk.

**Figure 2.** Temporal evolution of (**a**) C15:0, (**b**) C17:0, (**c**) anteisoC15:0) and (**d**) anteisoC17:0 from January to July in sheep (light grey line) and in goats (dark grey line) cheeses.

Changes in the animals' diet during the period under investigation, and in particular the progressive modification of pasture grazed by animals, can be the most probable explanation for the observed FA patterns. In fact, the typical lactation curve of the Sarda sheep breed follows the natural availability, both in quantitative and qualitative terms, of pastures [38]. As lactation progress, there is a worsening of nutritional quality of pasture characterized by increase of fiber and a reduction of protein and FA contents, especially alpha-linolenic acid (C18:3n3), which concentrations decrease in mature grass [39,40]. This is supported by evolution of C18:3in cheeses samples as lactation progress (data not reported), that confirmed the pattern of Sarda ewes previously observed [28].

Different mechanisms could be hypothesized to explain this pattern. One mechanism could be related to the decrease of polyunsaturated fatty acids (PUFA) in the pasture that could have reduced the toxic effect of unsaturated lipids on microbial growth [41,42], especially of cellulolytic bacteria [43]. This could increase the rumen microbial abundance, followed by higher rumen outflow of fatty acid of microbial origin. A second possible explanation could be the positive association between OBCFA and dietary fiber as with progress of lactation there is an increase of neutral detergent fiber (NDF), especially acid detergent lignin (ADL) in pasture. This is supported by previously observation in dairy goats where dietary NDF was found the most important factor of variation in lipid composition of bacteria; and in dairy cows where the proportion of odd- and branched-chain FA increased and those of even-chain saturated FA decreased with increasing forage [44].

The contents of almost all OBCFA in goat cheese differed from that of sheep, confirming the findings of the experiment 1. It is noteworthy that the goat's milk processed into cheese in this Mediterranean area comes from animals raised in extensive system, characterized by natural pastures rich also in shrubs and essences of the Mediterranean maquis, which contains tannins. These different pasture conditions among species could be a possible explanation of the different concentrations in OBCFA observed between sheep and goats [45].

In order to assess the overall differences among sheep and goat cheese in terms of FA composition, a multivariate statistical analysis was carried out on the detailed FA profile. In particular, PCA was used because it is a useful instrument to reduce the complexity of a multivariate space, such as that of FA of cheese, and has demonstrated to be able to separate samples (e.g., of different origin, dietary treatment, season of production, lactation stage), based exclusively on the FA composition [46–48]. An important application of PCA could be, as an example, the authentication of dairy products according to their lipid profile [49].

The first five principal components (PCs) explained the 80% of the total variance, with the two first PCs accounting for more than 50% (PC1 and PC2, 36% and 22%, respectively). The plot of the scores of the PC1 and PC2 (Figure 3a) allows the clustering of cheeses according to the season of milk production (PC1) or to their species of origin (PC2). Cheeses produced during winter–early spring season had positive scores for PC1, whereas those produced during late-spring–summer season had negative scores. Regarding the PC2, goat and sheep cheese had positive and negative scores, respectively. In the same plot, can be also observed that the Fa profile of sheep cheese was characterized by a larger variability compared to goat one, across the seasons of production.

The loadings of the first two PCs are plotted in Figure 3b. PC loadings can be used to describe and explain the main FA differences among the clusters of cheese observed in Figure 3a, i.e., the plot of the scores (cheese produced in different season or from milk of different species).

The original variables (i.e., single fatty acids) exhibiting the highest loading values for the 2 PCs (both negatives and positives) can be used to describe the main differences, in terms of FA, among samples belonging to the cheese of the two species and, also, produced in different seasons.

PC1 had high positive loadings for short chain FA (C4:0, C6:0, C8:0, C10:0), alpha linolenic acid, vaccenic and rumenic acid. High negative loadings for PC1 were observed for different long chain saturated FA (C20:0, C22:0, C24:0) and for the C18:1cis-9. This pattern was quite consistent with the milk (and cheese) FA variation during seasons, due to the worsening of pasture quality: in particular, the decrease of alpha-linolenic, vaccenic and rumenic acids in dairy products is correlated to the decrease in pasture of alpha-linolenic acid content from winter to summer. In addition, PC1 loadings can be also related to the lactation progress (high short chain FA in early lactation and increase in C18:1 cis-9 in mid and late lactation, in response to the energy requirement). This double effect of pasture quality and animal lactation stage on the FA composition of cheese (particularly in sheep cheese) was recently reported by [48]. As aforementioned, it is likely to be a consequence of the typical seasonality of lambing of Sarda dairy sheep.

PC2 had larger negative loadings for rumenic acid, C4:0 and some OBCFA, in particular, anteisoC15:0, C15:0, and some isomer of C18:1, including vaccenic acid. Higher positive loadings of PC2 were observed for C18:0 and C10:0. This pattern suggested a better FA composition of sheep cheese, from a nutritional point of view, considering that FA with negative loadings for PC2 (including some OBCFA), and, therefore, mostly related to sheep cheese, have been associated to beneficial aspect on human health, that were the subject of this work. Due to its nutritional importance, the rumenic acid (most abundant and important CLA isomer) should be mentioned. This FA has been received particular attention due to the important beneficial healthy effects shown in vitro, animal and human as well [50,51]. The higher concentration in sheep milk (and cheese) compared to that of goat has previously been reported [45].

**Figure 2.** *Cont.*

**Figure 3.** Scores (**a**) and loadings (**b**) plots of the two first principal components explaining 36 and 22% of the total variance, respectively. In Figure 3a, circle identifies goat cheese, triangle identifies sheep cheese; green symbols identified cheese produced in winter early spring, red symbols identified cheese produced in early spring–summer.

#### **4. Conclusions**

Results of the present survey showed that odd and branched chain fatty acids concentrations in the fat of fresh cheese is related to their content in both sheep and goats' milk. The seasonal evolution evidenced that OBCFA concentrations in milk fat increase with advancing lactation, probably due to variation in feeding technique typical of extensive system of Mediterranean area characterized by animal grazing on pasture. As lactation progress, the pasture availability and quality gets worse, leading to a reduction in the content of PUFA and protein and to an increase in fiber content. These results showed that cheese could be an important source of OBCFA, but the concentration could vary according to the lactation stage. The use of PCA on the detailed FA profile of cheese confirmed the higher nutritional quality of sheep cheese for beneficial FA including OBCFA compared to the goat one, and the importance of the period of sampling in the definition of the fatty acids profile.

**Author Contributions:** Individual contributions: A.N.: Conceptualization, Statistics and Writing. F.C.: Methodology and Analysis. A.C.: Software and PCA analysis. G.P.: Supervision and Writing. G.B.: Writing Original Draft, Preparation. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by "Fondo di Ateneo per la ricerca 2019 of University of Sassari".

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available on request from the corresponding author.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


### *Article* **LC-QTOF/MS Untargeted Metabolomics of Sheep Milk under Cocoa Husks Enriched Diet**

**Cristina Manis 1, Paola Scano 1, Anna Nudda 2, Silvia Carta 2, Giuseppe Pulina <sup>2</sup> and Pierluigi Caboni 1,\***


**\*** Correspondence: caboni@unica.it; Tel.: +39-070-6758617

**Abstract:** The aim of this work was to evaluate, by an untargeted metabolomics approach, changes of milk metabolites induced by the replacement of soybean hulls with cocoa husks in the ewes' diet. Animals were fed with a soybean diet integrated with 50 or 100 g/d of cacao husks. Milk samples were analyzed by an ultra high performance liquid chromatograph coupled to a time of flight mass spectrometer (UHPLC-QTOF-MS) platform. Multivariate statistical analysis showed that the time of sampling profoundly affected metabolite levels, while differences between treatments were evident at the fourth week of sampling. Cocoa husks seem to induce level changes of milk metabolites implicated in the thyroid hormone metabolism and ubiquinol-10 biosynthesis.

**Keywords:** mass spectrometry; animals management; thyroid hormone metabolism; ubiquinol-10 biosynthesis

#### **1. Introduction**

In the Mediterranean basin, where dairy sheep and goats breeding is widespread, the use of agricultural by-products in the diet of small ruminants is an ancient practice. Today the great availability of agro-industrial by-products produced worldwide opens the door to other products to be tested. Noteworthy, the use of co-feeds can contribute to reducing the impact associated with ruminant's management. Several positive aspects were observed when agro-industrial by-products (e.g., grape, olive, tomato, citrus pulp and myrtle residues) were included in the diet of small dairy ruminants. In particular, it has been reported that the diet implementation with these by-products has beneficial effects on ruminal metabolism, animal health, and quality of derived products [1–3].

Cocoa husks represent the part of cocoa pod left over and are the principal by-products derived from *Theobroma cacao* L., representing an important crop for many tropical developing countries. Cocoa husks are obtained after the removal of the cocoa beans, the principal commercial product of cocoa processing [3]. Husks represent 70–75% of the fruit weight [4] and, besides limited local use, are considered as waste [5]. Cocoa pod husks along with other industrial by-products such as cocoa bean shells, cocoa bean meal and cocoa germ can be used as feed especially in the developing countries.

Purine alkaloids, including theobromine, theophylline and caffeine, and other methylxanthines, play a substantial role in pharmacology and food chemistry. A limited number of plant species accumulates purine alkaloids, such as caffeine and theobromine, which are synthesized from xanthosine, a catabolite of purine nucleotides. The most widely distributed methylxanthine in the plant kingdom is the antioxidant caffeine, which accumulates in leaves and seeds of tea (*Camellia sinensis*), coffee (*Coffea arabica*) and a limited number of other species. Considerable amounts of theobromine are stored in the seeds of cacao (*Theobroma cacao*) [6]. Moreover, from a pharmacokinetic point of view, theobromine can be easily adsorbed, with a low protein binding affinity and distributed into body

**Citation:** Manis, C.; Scano, P.; Nudda, A.; Carta, S.; Pulina, G.; Caboni, P. LC-QTOF/MS Untargeted Metabolomics of Sheep Milk under Cocoa Husks Enriched Diet. *Dairy* **2021**, *2*, 112–121. https://doi.org/ 10.3390/dairy2010011

Received: 15 December 2020 Accepted: 26 January 2021 Published: 22 February 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

tissues with no specific tissue accumulation [7]. Despite the positive health benefits of consuming foods containing purine alkaloids [8] one of the limits in the use of cocoa husks as supplementary feed for animals is the fact that in certain amount, theobromine, can be toxic [9]. The EU Scientific Opinion of the Panel on Contaminants in the Food Chain reported that when exposed to theobromine dairy cows showed reduced milk yield, increase in fat levels, and adverse effects such as hyperexcitability, sweating, increase in respiration and heart rates [10]. Reported levels of theobromine in cocoa husk meal are in the range of 1.5–4.0 g per kg of feed material and UE regulations set for ruminants the maximum levels of theobromine in feed materials at 300 mg/kg. On the other side, cocoa husks on a dry matter basis (DM) are a good source of neutral detergent fiber (NFC 38– 44%) and high variable concentrations of non-fiber carbohydrate (NFC 17.5–47% on DM), crude proteins (CP 2.1–9.1% on DM), and lipids (ethyl esters 0.6–4.7% on DM); phenolic compounds account for 4.6–6.9% of gallic acid/100 g DM. However, the high content of acid detergent lignin (19.6% on DM) derived from cocoa pericarp affects the digestibility of this by-product, limiting its further use in animal feeding [11].

Metabolomics is a largely used branch of omics sciences used to explore the composition, relative levels, dynamics, and interactions amongst metabolites in an organism or biological system in response to various stimuli, diet or treatments. In the field of dairy, metabolomics has been used for example in the authentication of beef production systems [12], in the milk to investigate metabolite differences of healthy, subclinical, and clinical mastitis cows [13] or sheep milk grazed on different grazing systems [14], while in cheese samples used for discriminate cheese produced from raw or thermized milk [15].

To our knowledge, experimental data on theobromine carry over in milk are not available, as well as the metabolic effects of cocoa husks diet supplementation. In a previous work, milk and blood parameters of ewes supplemented with different doses of cocoa husks were evaluated [5]. In this experiment DMI and milk yield were not affected; however, milk protein content increased with cocoa husk supplementation. Furthermore, the cocoa husks supplementation lowered milk somatic cell count (SCC), suggesting a beneficial role of these by-products on mammary health status. Both systemic and local health conditions were good as evidenced by hematological parameters and electrophoretic profile of serum protein fractions [5]. Both systemic and local health conditions were good as evidenced by hematological parameters and electrophoretic profile of serum protein fractions [5]. Taking into consideration these results, we decided to investigate the potential beneficial nutritional outcomes by exploring, for the first time, milk metabolite profile changes of healthy ewes fed with a diet supplemented with 100 g/d per head soybean hulls that were replaced with 50 and 100 g/d per head of cocoa husks, respectively. In order to pursue these objectives, we used a liquid chromatography mass spectrometry metabolomics platform and multivariate statistical analysis.

#### **2. Materials and Methods**

The experiment was approved by the Ethics Committee of the University of Sassari (no. 54584/2018).

*Animals and diets.* As previously described [5], twenty-four Sarda dairy ewes in mid lactation stage (120 days in milking), at third parity, with an average body weight of 42.5 kg were divided into three homogenous groups. Each experimental group received a total mixed ration (2.58 kg/d per head) as the basal diet and received (a) a supplement of 100 g/d per head of soybean hulls (CH0 group); (b) half of the soybean hulls supplement was replaced with 50 g/d per head of cocoa husks (CH50 group); and (c) soybean hulls was totally replaced with 100 g/d per head of cocoa husks (CH100 group). Theobromine concentration was 130 ± 12 and 253 ± 27 mg/kg of DM for the CH50 and CH100 diets, respectively. The CH was administered individually during milking, mixed with beet pulp to increase palatability. The experiment lasted 8 weeks (from May to June 2018), with the last 4 weeks of sampling collection. Individual milk samples were collected weekly at morning milking.

*Chemicals.* Sigma Aldrich (Milan, Italy) reagents were used for the analysis. Bi-distilled water was obtained with a MilliQ purification system (Millipore, Milan, Italy).

*Sample preparation for UHPLC-QTOF/MS analysis.* Individual milk samples from the morning milking were stored into sterile plastic Falcon tubes at −80 ◦C before analysis. Prior to UHPLC-QTOF/MS analysis, the 107 individual sheep milk samples were thawed and thoroughly vortex mixed. A total of 150 μL of milk samples were transferred to Eppendorf tube containing 10 uL of the internal mixture of standards (Splash, Lipidomics, Sigma Aldrich, Milan, Italy) and added with 525 μL of methanol and 525 μL of MTBE. Samples were then mixed by vortexing for 1 min. Next, samples were centrifuged at 4000 rpm for 15 min. Subsequently, a second round of centrifuge was performed at 12000 rpm for 5 min. Before transferring to autosampler vials, the supernatant was filtered through a 0.22 μm MS nylon syringe filter.

*UHPLC-QTOF/MS analysis.* After extraction with MTBE/methanol, the supernatant of milk samples was analyzed with a 6560 Q-TOF (Agilent Technologies, Palo Alto, CA, USA) coupled with an Agilent 1290 Infinity II LC system. An aliquot of 1.0 μL from each sample was injected in a Kinetex 5 μm EVO C18 100 A, 150 mm × 2.1 μm column (Agilent Technologies, Palo Alto, CA, USA). The column was maintained at 50 ◦C at a flow rate of 0.5 mL/min. The mobile phases consisted of (A) methanol:water (90:10, *v*/*v*) with ammonium acetate (10 mM) and (B) acetonitrile:methanol:2-propanol (20:30:50, *v*/*v*) with ammonium acetate (10 mM). The chromatographic separation was obtained using the following gradient—0 min 70% B kept for 1 min; 1–3.5 min 86% B; 3.5–10 min 86% B; 10.1–17 min 100% B; 17.1–19 min 70% B. The analytical setup used was equipped with an Agilent jet stream technology source which was operated in both positive and negative ion modes with the following parameters: gas temperature,200 ◦C; gas flow (nitrogen) 10 L/min; nebulizer gas (nitrogen), 50 psig; sheath gas temperature, 300 ◦C; sheath gas flow, 12 L/min; capillary voltage 3500 V for positive and 3000 V for negative; nozzle voltage 0 V; fragmentor 150 V; skimmer 65 V, octapole RF 7550 V; mass range, 40–1700 *m/z*; capillary voltage, 3.5 kV; collision energy 20 eV in positive and 25 eV in negative mode. An Agilent MassHunter software was used for instrument control (revision B.09.00).

*Multivariate statistical data analysis (MVA).* The QTOF-MS data were submitted to the web platform XCMS [16]. The workflow allows feature detection, retention time correction, alignment, annotation, statistical analysis, and data visualization. Retention time, *m/z* values, and intensities of each feature obtained from XCMS were then submitted to multivariate statistical analysis (MVA), as implemented in SIMCA-P+ software (version 14.1, Umetrics, Umeå, Sweden). Prior to MVA, QTOF-MS features were mean centered and scaled to unit variance column-wise. Principal component analysis (PCA) was performed to investigate sample distributions, deviating features and prevailing trends. This technique reduced the dimensionality of data set retaining most of the variance. Results are shown as scatter plots of scores and loadings in the first principal components, reporting sample and variable displacement in the hyperplane, respectively. The partial least squaresdiscriminant analysis (PLS-DA) and its orthogonal variant (OPLS-DA) were performed for classification of samples and identification of the most discriminant variables. The quality of the models was evaluated based on the cumulative parameters R2Y and Q2Y, the latter estimated by cross validation, and tested for overfitting using a y-table permutation test (*n* = 400). The variable importance in projection (VIP) score summarizes the contribution of each variable to the model. The VIP scores in the predictive component were analyzed and only those metabolites having VIP values > 1 were considered as discriminant between the classes [17].

#### **3. Results and Discussion**

This untargeted metabolomics approach on dairy ewes was conceived to provide a molecular and biological insights on milk metabolite changes induced by a diet integrated with 50 or 100 g/d of cacao husks. Moreover, we designed this cost-efficient metabolomics study to measure predictive metabolite biomarkers to be used in exploring feed suitability and efficiency.

For these purposes, a total of 107 sheep milk samples were analyzed with a UHPLC-QTOF/MS platform, (Figure S1 reported the LC-MS chromatograms). After analysis, raw mass spectrometry data acquired were then uploaded to the XCMS platform, which generated an average of 7700 and 3300 features for the positive and negative modes, respectively. The XCMS outputs were submitted to multivariate statistical analysis using the software SIMCA.

A PCA was performed to explore similarities of milk sample metabolic profiles of the three diets over the experimental sampling period. Results indicated that metabolic similarities of milk samples were mainly due to the time of sampling, as clearly visible in the score plots shown in Figure 1 for the positive and negative ionization modes.

**Figure 1.** Principal component analysis (PCA) score plots of milk samples. Samples are colored by week of treatment (first, blue; second, red; third, green; and fourth, light blue circles). Metabolite data collected in the liquid chromatography-mass spectrometry (LCMS) (**a**) positive and (**b**) negative ionization modes.

Information collected for the two ionization modes yielded similar results (i.e., samples collected in the two first weeks of treatment clustered together as a result of similar metabolite profile, while milk samples from the third and the last week of experiment can be easily differentiated along the second principal component (PC2) in the positive and negative ionization modes). Interestingly, the same pattern was observed for milk yield, as shown in Figure S2, whereby the first two weeks gave comparable average values for all the group samples, followed by a drastic milk yield drop after the second week of sampling (from the 13 to 26 of June). A decrease of milk yield and lactose content and an increase in fat and protein content, due to a milk concentration effect, during the late-spring/summer seasons were already reported as a natural evolution of the lactation curve in dairy Sarda sheep [18]. Probably these milk qualitative and quantitative changes due to lactation seasonality have caused an overall modification of the metabolic profiles captured by the PCA.

Attempts to differentiate samples based on the diets (CH0, CH50, and CH100) were carried out by performing three-classes PLS-DA for each time point for positive and negative modes. The multivariate statistical analysis failed to discriminate samples for the first three weeks of sampling (data not shown) whereas at the fourth week of treatment, samples were properly classified, with high statistical confidence (Figure S3), indicating that the effects of the cocoa husks supplementation on milk metabolic profile were evident only in the last time point. Subsequently, with the purpose of finding discriminant metabolites, pairwise OPLS-DA were carried out comparing the three groups of milk samples collected in the last week of treatment (Figure 2). In Tables 1 and 2 we reported the list of discriminant metabolites, their biological significance, and regulation in sheep milk samples.

**Figure 2.** Score plots of pairwise orthogonal partial least squares-discriminant analysis (OPLS-DA) of milk samples data in liquid chromatography-mass spectrometry positive ((**a**) CH0 vs. CH100 and (**b**) CH50 vs. CH0), and negative ((**c**) CH0 vs. CH100 and (**d**) CH50 vs. CH0) ionization modes. CH0 = green circles; CH50 = red boxes; CH100 = light blue triangles.


**Table 1.** Orthogonal partial least squares-discriminant analysis (OPLS-DA) main discriminant metabolites as variable importance for prediction (VIP > 1) detected in positive ionization modes.

**Table 2.** Orthogonal partial least squares-discriminant analysis (OPLS-DA) main discriminant metabolites (VIP>1) detected in negative ionization modes.


From the analysis of discriminant metabolites, as summarized in Tables 1 and 2, we found that the addition of cocoa husks in the dairy lactating sheep diet affected levels of different metabolites which are involved in different metabolic pathways, including ubiquinol synthesis and thyroid hormone metabolism. In particular, we found that in milk samples of lactating ewes the thyroid hormone L-thyroxine (T4) and tetraiodothyroacetate, or tetrac (T4A), were down regulated for both CH50 and CH100 groups, thus suggesting a modification of the thyroid hormone metabolism. Levels of T4 and T4A are shown in Figure 3, indicating that, in the first two weeks of treatment, group samples showed similar average values. Levels where higher in the third week and in the last week CH50 and CH100 groups showed lower values when compared to CH0.

**Figure 3.** Box plots of (**a**) thyroxine (T4) and (**b**) tetraiodothyroacetate (T4A) levels in milk samples for the CH0, CH50, and CH100 groups in the last four weeks of treatments (week 1 = blue, 2 = red, 3 = green, 4 = cyan).

The thyroid hormones T4 and 3,5,3 -triiodo-L-thyronine (T3) are regulated at many metabolic levels. In addition to deiodination (thyroid hormone metabolism I pathway), other routes of thyroid hormone metabolism have been described and metabolically active metabolites have been identified. In the thyroid hormone metabolism pathway (via degradation), small amounts of T4 may also be oxidatively deaminated at the alanine side chain to yield T4A. The oxidative deamination due the kidney L-amino acid oxidase seems the most probable explanation [19,20]. The acetic acid analogues of thyroid hormones are formed in the liver and other peripheral tissue by a two-stage process involving transamination followed by oxidative decarboxylation [21]. T4A may be conjugated to form ether glucuronides or ester glucuronides of these metabolites may form by conjugation of the phenolic hydroxyl group, or carboxyl group respectively [20].

Thyroid hormones in milk should derive from circulating blood, then undergoing both desiodative and non-desiodative metabolism. A reduction in milk T4 could; therefore, be derived from a reduction in circulating T4 levels or from the reduction in T4 level available at the mammary gland, due to the effect of interferents on the transporters of thyroid hormones, or, finally, by increased desiodation of T4. The reduction of milk T4A level seems the logical direct consequence of the reduction of T4, deriving the T4A directly from T4. Thyroid hormones play a relatively important role in lactational processes. In small ruminants and, in particular, in free grazing animals which are subjected to seasonality, endogenous (breed, age, gender, physiological state) and environmental factors (climate, season, with a primary role of nutrition) are able to impact on thyroid activity and hormone blood levels [22]. During heat stress, an increase of blood T4 concentration and a decrease of milk production were reported in sheep [23,24]. This trend was evident also in our milk samples, which were collected in June, where, during the experimental period, in the CH0 group the milk T4 levels inversely follow the milk production trend (Figures S2 and S4). Differently from CH0, in the last week of sampling, milk T4 levels dropped in the CH50 and CH100 groups (Figure S4), thus indicating a goitrogenic effect of cocoa husks on milk T4. In the rat, Wolff and Varrone (1969) reported that methyl-xanthines (i.e., theobromine, theophylline, and caffeine) showed a goitrogenic effect and are; therefore, to be considered weak anti-thyroid agents if present in the diet at concentrations varying between 5–20 microM [25]. Further studies will be carried out to better evaluate this effect in sheep, measuring the volume of the thyroid glands and the concentration of T4 and TSH in ewes' blood after consumption of cocoa husks. At the milk level, which only contains traces of iodothyronines, which; therefore, do not represent a significant nutrient of the milk itself, as opposed to iodine, any reduction in concentration should not be particularly relevant for its nutritional value.

Notably, phenylalanine was found up-regulated in CH50 and CH100 milk samples. Considering that phenylalanine is the precursor of tyrosine and the thyroid hormone T4 is a tetraiodinated derivative of tyrosine [22], this fact confirms a dysregulation of thyroid hormones metabolism.

The discriminant metabolites: 2-methoxy-6-(all-trans-octaprenyl)phenol, 6-methoxy-3 methyl-2-all-trans-decaprenyl-1,4-benzoquinol, and 3,4-dihydroxy-5-all-trans-decaprenyl benzoate, involved in the ubiquinol-10 biosynthesis, were found down regulated in milk samples of the CH50 and CH100 groups with respect to CH0 (Table 1). Coenzyme Q (CoQ10), also known as ubiquinone, is an essential component of the mitochondrial respiratory chain and participates in aerobic cell respiration. The reduced ubiquinol form of CoQ10 is a chain-breaking antioxidant, decreasing oxidative damage caused by lipid peroxidation within mitochondria. The cocoa husks feeding regime, showed markedly lower levels of intermediate metabolites of the CoQ10 pathway compared to the controls. To the best of our knowledge, no studies of theobromine effects on CoQ10 pathway have been widely reported, but interference of phenolic compounds has been found in liver of rats [26]. Specifically, the authors found that pentagalloylglucose, which is a component of tannic acid, was able to interfere with the CoQ10 pathway: At low concentration inhibiting the electron transport system, while at high concentration impairing the structural integrity of the mitochondrial membrane. Consistently, CoQ10 was found to affect the pharmacokinetic parameters of theophylline [27], an alkaloid of the xanthine family.

Remarkably, in the negative ionization mode we found different discriminant metabolites involved in the purine metabolism (Table 1). Inosine is a purine nucleoside formed by deamination of adenosine and, in the context of purine metabolism, is further transformed to hypoxanthine by the enzyme phosphate alpha-D-ribosyltransferase and then to xanthine and to uric acid. The fact that inosine was found upregulated in milk samples of sheep consuming cocoa husks may be associated to the inhibition of the enzyme phosphate alpha-D-ribosyltransferase, in order to reduce the increase of xanthines and uric acid derived from the metabolism of theobromine.

Our experimental findings are in accordance with scientific opinion of EFSA [10] on theobromine as an undesired substance in animal feed. In the EFSA report, theobromine was reported as showing moderate acute toxicity on dog and rodents, while in horses, which are particularly susceptible to theobromine, the liver and thyroid were affected, and pigs showed growth retardation, diarrhea and lethargy. Particular attention should be paid on the levels of theobromine and other alkaloids present in the industrial by-products in the preparing of the dairy sheep's diet.

#### **4. Conclusions**

In this paper we show evidence, with an untargeted metabolomics approach, that sheep milk metabolite profiles are profoundly altered on the basis of an eight-week diet containing 50 and 100 g/d of cocoa husks. In particular these effects were more pronounced in the eighth week, at the end of the treatment. Moreover, this approach allowed us to highlight changes in multiple metabolic pathways. Importantly, the thyroid hormone Lthyroxine (T4) and tetraiodothyroacetate (T4A) were found down regulated for both CH50 and CH100 groups, thus suggesting a modification of the thyroid hormone metabolism. Furthermore, different ubiquinol-10 biosynthesis metabolites were found down regulated in milk samples of the CH50 and CH100. Biochemical mechanisms related to these changes remain unclear, while thyroid hormones levels should be further investigated in the ewes' blood. In conclusion, taking into account these results, particular attention should be paid when using cocoa by-products in small ruminant feed management.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/2624-862 X/2/1/11/s1: Figure S1. Milk sample overlaid chromatograms of CH0, CH50, and CH100 groups acquired in positive (a) and negative (b) ionization mode. Figure S2. Temporal evolution of milk yield in the control (CH0) and in the two cocoa husks supplemented groups (CH50 and CH100). Figure S3. Partial least square discriminant analysis (PLS-DA) score plots of milk samples of the

experimental groups (CH0, CH50, and CH100) collected in the last week of treatment. Lipid data collected in the LCMS positive (top, *R*2Y = 0.62, Q2Y = 0.27) and negative (bottom, *R*2Y = 0.98, Q2Y = 0.63) ionization mode. Figure S4. Temporal evolution of thyroxine (T4) milk levels in the control (CH0) and in the two cocoa husks supplemented groups (CH50 and CH100).

**Author Contributions:** Sample analysis, statistical analysis, wrote original draft, C.M. and P.C.; statistical analysis, wrote original draft, P.S.; conceptualization, wrote original draft, A.N. and G.P.; animal experiment and milk sample collection, S.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** The research was partly funded by Mignini-Petrini company.

**Institutional Review Board Statement:** The study was conducted according to the guidelines of the Declaration of Helsinki, and approved by the Ethics Committee of the University of Sassari (no. 54584/2018).

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available on request from the corresponding author.

**Acknowledgments:** The authors are grateful to Stefano Mariotti and Elio Acquas for helpful discussions.

**Conflicts of Interest:** The authors declare that the sponsor had no role in the study design and interpretation of the results.

#### **References**


## *Article* **Quality Control in Fiore Sardo PDO Cheese: Detection of Heat Treatment Application and Production Chain by MRI Relaxometry and Image Analysis**

**Roberto Anedda \*, Riccardo Melis and Elena Curti**

Porto Conte Ricerche s.r.l., S.P. 55 Porto Conte-Capo Caccia, Km 8.400 Loc. Tramariglio, 07041 Alghero, SS, Italy; melis@portocontericerche.it (R.M.); curti@portocontericerche.it (E.C.) **\*** Correspondence: anedda@portocontericerche.it; Tel.: +39-079-998-578

**Abstract:** Fiore Sardo (FS), a traditional Italian cheese, is present in the market as a heterogeneous variety of products. The use of heat-treated (HT) milk is forbidden by the official production protocol, but no official analytical method able to detect heat application is yet available. Here, a combined magnetic resonance imaging (MRI) relaxometry and image analysis approach to recognize FS made from raw milk is presented. Artisanal FS cheeses were produced from raw milk (RC) by five shepherds in accordance with the official protocol. They were compared to HT-milk counterparts (HTC). Additionally, industrially manufactured commercial FS cheeses (I) were also purchased and compared to RC and HTC. Relaxometry data of FS indicated the presence of two water populations; the ratio of characteristic relaxation time constant T2 and area fraction (Score, Ṩ) of the fastest relaxing population was used to compare RC, HTC and I samples. RC from HTC were successfully discriminated, the latter exhibiting lower Ṩ (enhanced protein hydration). I cheeses exhibited the lowest Ṩ values, sometimes comparable to HTC. Since visual appearance of RC and HTC is appreciably different, an image analysis deep learning approach using MRI and photographic pictures was adopted to discriminate the two productions, with promising percentages (>93%).

**Keywords:** magnetic resonance imaging (MRI); image analysis; Fiore Sardo; cheese; microstructure; dairy chemistry; thermised milk; raw milk; protected designation of origin

#### **1. Introduction**

Protected Designation of Origin (PDO) describes a group of agricultural/food products originating from a specific geographical area, for which production, processing and preparation must take place according to a recognised know-how [1].

In the European Union, 185 cheeses are registered as PDO products; this number rises to 231 if cheeses registered as Protected Geographical Indications (PGI) are considered [2]. A first classification of these products can be made based on milk processing. A total of 91 cheeses out of 231 (PDO + PGI) compulsorily require the use of raw milk; in other cases, the use of both raw and pasteurized/heat-treated milk is allowed. The lack of a specific requirement allows producers to freely decide how to process milk, leading to cheeses characterized by very different attributes [3].

The effects of milk heat treatment on cheese features have been extensively investigated and explained. From a microbiological and sensory point of view, the heat treatment depletes pathogens along with the characteristic microbiota, resulting in a more uniform final product that exhibits a less characteristic flavour profile compared to the raw milk counterpart [4]. From a technological point of view, though the application of a heat treatment guarantees improved hygienic conditions during processing and a better safety of the final product, it can result in an impaired coagulation and syneresis and in an alteration of cheese rheology and texture [5,6].

**Citation:** Anedda, R.; Melis, R.; Curti, E. Quality Control in Fiore Sardo PDO Cheese: Detection of Heat Treatment Application and Production Chain by MRI Relaxometry and Image Analysis. *Dairy* **2021**, *2*, 270–287. https:// doi.org/10.3390/dairy2020023

Academic Editors: Pierluigi Caboni and Paola Scano

Received: 4 March 2021 Accepted: 13 May 2021 Published: 26 May 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Fiore Sardo PDO (FS) cheese is the oldest cheese of Sardinia (Italy), being historically produced by shepherds in very small artisanal cheese factories. It must be obtained exclusively using raw whole milk from the Sarda Sheep breed. The official cheesemaking protocol still follows the ancient traditional process [7]. The peculiarities of the manufacturing process confer to Fiore Sardo a firm, crumbly and grainy texture, which easily breaks into flakes, and a savoury, piquant and smoky flavour [8–11]. Market offer is characterized by a quite heterogeneous variety of Fiore Sardo PDO, with very different quality attributes, with almost 90% of the product supplied by the dairy industry (Agris Sardegna, 2015).

In particular, the use of heat-treated milk for FS production has been associated with common industrial processing practices [12–15]. Previous studies concluded that this practice cannot be excluded [16] and that it is not unrealistic that some FS are manufactured using heat-treated milk, against the specifications indicated in the EC Regulation. From a scientific point of view, several chemical and physical characteristics of some productions were found to be reasonably explained by assuming the effects of thermal treatments on milk [16]. This is also supported by the fact that no official analytical methods able to discriminate between ewe milk cheese obtained from raw and thermised milk are available [17].

Some analytical techniques able to discriminate products made from raw milk and heat-treated milk have been proposed for product quality safeguards [17–21].

Analytical approaches based on nuclear magnetic resonance (NMR) relaxometry are particularly interesting. This technique was previously applied to describe the thermal denaturation of whey proteins and other milk components. Lambelet et al. [22], for example, demonstrated that irreversible thermal denaturation of whey proteins can be detected by low-field NMR. They observed that when NMR transverse relaxation rates were plotted versus heating temperature, a sigmoid curve was obtained, with relaxation rate values showing a steep increase from about 60 ◦C to about 80 ◦C, and a flex corresponding to milk thermisation/pasteurization temperatures usually adopted in dairy industrial processing (64 ◦C to 72 ◦C range). Later work by the same research group extended the results to other milk components and dyalized skim milk [23], suggesting the suitability of NMR relaxometry for measuring the thermal denaturation of milk proteins. Recent works, based on a high sample numerosity and exploiting an experimental setup more closely resembling real industrial processing conditions, demonstrated that the effect of heat treatment on milk is transferred to NMR relaxometry characteristics of both rennet curd and aged cheese. Raw and heat-treated milk [24] and derived curds [25,26] could be successfully discriminated by benchtop low-field NMR relaxometry. A similar discriminative ability of NMR relaxometry was also observed in ewe milk mature cheese by high-field magnetic resonance imaging (MRI) [19].

The aforementioned investigations demonstrated the sensitivity of NMR relaxation parameters to key cheesemaking processes (mainly thermal treatments of milk) in real dairy systems (whole milk, curds, mature cheese), suggesting the suitability of NMR relaxometry as a valid and useful tool in quality control activities in the dairy industry and for the safeguard of typical dairy productions.

In this study, heat treatments of milk were intentionally applied to produce FS cheese, and NMR relaxometry of raw and heat-treated FS were compared. Moreover, a comparison between artisanal and some industrial FS, the latter purchased from the market and visibly different from the typical FS PDO, was carried out for a more comprehensive evaluation of the effect of production chain on relaxometry features and aiming to optimize the exploitation of NMR relaxometry to assess the adherence to the official production protocol. MRI was used to compare heat-treated samples to their raw counterparts. A supplemental investigation was also carried out by exploiting image analysis protocols on a selected subset of MRI images and digital pictures in order to correlate process-induced changes with MRI image features and the visual appearance of samples.

#### **2. Materials and Methods**

Two sets of experiments were carried out, as summarized in Table S1. The first dataset (Dataset 1) was acquired by analysing samples made by different artisanal producers in different seasons. Additionally, a selection of industrial cheeses was purchased from the local market and subjected to the same analytical protocol used to characterize samples from Dataset 1.

#### *2.1. Production of Fiore Sardo (FS) Cheese*

#### 2.1.1. Dataset 1

Five shepherds (S1, S2, S3, S4, S5) who were artisanal manufacturers of Fiore Sardo were selected for cheese sample production. Each shepherd produced a total of 20 cheese wheels in two different seasons (March–April 2019, Season 1; January–February 2020, Season 2; 100 cheese wheels in total for the two seasons). Ten wheels out of 20 were produced from raw milk, while the remaining ten were produced from the same milk after thermisation. Milk thermisation was carried out in different industrial plants, using plate heat exchanger modules. Cheesemaking was performed according to the official cheesemaking protocol [7,9], but milk thermisation was used for selected samples. Aged cheese samples were collected after 105 and 180 days of ripening in each shepherd's cellar. One quarter portion of each cheese wheel was sampled at each of the abovementioned time points, for a total of 10 wheels per producer per ripening time, 5 of which were made from raw milk (RC) and 5 from heat-treated milk (HTC).

#### 2.1.2. Dataset 2

Industrial FS cheeses with similar ripening (6–8 months) were purchased from the local market in different seasons (January–November). Both cheeses manufactured by industrial dairy factories (11 samples) and from a maturer industry (2 samples) were purchased. Maturer industries buy Fiore Sardo from artisan producers and ripen these cheeses in their industrial cellars for several months, until they can be released to the market, i.e., after at least 105 days of ripening, according to the official cheesemaking protocol. Fiore Sardo PDO cheeses from the maturer industry are sold with the industrial brands, but milk is processed by shepherds and small artisans according to traditional processes [9]. A schematic representation of sample production and experimental plans (relative to both Dataset 1 and Dataset 2) is reported in Figure S1 (Supplementary Material).

#### *2.2. Magnetic Resonance Imaging (MRI) Analysis*

Each analytical sample of cheese subjected to MRI analysis represented the most central part of the cheese wheel, as schematically depicted in Figure S2.

In brief, each quarter was carefully cut in the portion corresponding to the wheel centre to obtain one cylindrical sample (diameter 2.3 cm, height 3–4 cm) to easily fit into a 50 mL Falcon tube. All samples were equilibrated to 22 ◦C for 1 h before MRI analysis. MRI experiments were performed using a Bruker Avance 300 MHz equipped with a microimaging Micro2.5 probe (Bruker Biospin, Karlsruhe, Germany) at room temperature (approximately 22–24 ◦C). A conventional Carr Purcell Meiboom Gill (CPMG) spin echo sequence was used with the following parameters: single 1 mm slice; echo time (TE) = 7.907 ms; repetition time (TR) = 3 s; echoes = 64; matrix = 128 × 128; number of excitations = 1–3; acquisition time = 25 min. To selectively analyse signals from water, the fat signal was suppressed, taking advantage of the different chemical shifts of fat and water at 7.05 T, with a 90◦ selective pulse applied before the conventional CPMG spin echo sequence at a frequency offset of 3.5 ppm with respect to water (preparatory fat saturation scheme provided by Bruker Topspin). All MR images were characterized by a signal to noise ratio higher than 100.

The AnalyzeNNLS software package was used to deconvolute the multi-exponential transversal relaxation time (T2) decay of water proton signals [27]. Bruker MRI image data were converted into multiecho image data (MEID) files (Matlab, R2020b, The Mathworks, Sherborn, MA). Regions of interest (ROIs) were manually drawn on the MEID file (five ROIs for each image) and automatically processed to obtain the geometric mean T2 (T2 relaxation time constant, corresponding to the maximum intensity of the peak, on a logarithmic scale) and the area fraction (obtained by dividing the sum of the T2 distribution within a desired region between T2 min and T2 max by the sum of all T2 distributions).

#### *2.3. Moisture Content*

Moisture content was measured using the rapid determination method for cheese (130 ◦C, 90 min; [28]).

#### *2.4. Statistical Analysis*

The effect of thermisation was evaluated by comparing MRI parameters (geometric mean T2 and area fractions) between RC and HTC by using univariate Student's *t*-test. An additional statistical comparison (one-way ANOVA, followed by Tukey post-hoc test) was performed to compare each shepherd's (i.e., artisanal) production to the industrial and maturer industry cheeses. Industrial cheeses (11 samples) were considered as a unique group of data, as were maturer cheeses (2 samples). Both tests were carried out with the MetaboAnalyst tool (https://www.metaboanalyst.ca, v.5, accessed 20 November 2020) [29], considering a significance threshold limit of *p* < 0. 05. Box plots were obtained using MetaboAnalyst.

#### *2.5. MRI Data Conversion*

T2 MRI Bruker Paravision raw 2dseq data from full Dataset 1 were first converted to Neuroimaging Informatics Technology Initiative (NIfTI, .nii) format using Bru2nii software (Bruker2NIfTI v1.0.20180303: by Matthew Brett, Andrew Janke, Mikaël Naveau, Chris Rorden, Windows 64-bit) with a resizing factor of 10. Then, the 64 slices from each NIfTI file were converted in JPEG format using a Matlab script first developed by Alex Laurence [30]. Only MRI images derived from the first MSME slice (TE = 7.907 ms) were used for deep transfer learning classification purposes.

#### *2.6. Photographic Acquisition and Processing of Cheese Paste Surfaces*

Snapshots of the surface of Fiore Sardo cheese from Dataset 1 produced in Season 1 at 105 days of ripening were taken using a Machine Vision System (MVS) assembled in our laboratory. Briefly, our MVS consisted of a portable plastified (PVC) white light-box (24 × 23 × 22 cm) equipped with a USB LED strip of white-light-emitting diode lamp (input 5 V, lumen 550, colour temperature 6500 K), assembled on a metallic support and directed upwards at a 90◦ angle with a cardboard box wall (Figure S3). Photos were taken without flash using a 25-megapixel-wide camera integrated in a Galaxy A50 smartphone. Cheese quarters were placed under a dark background and manually adjusted to be below the camera with a fixed distance of 25 cm. All operations were carried out in a dark room to minimize the effect of outside light. The acquired images (Joint Photographic Expert Group—JPEG; 24-bit colour depth, spatial resolution 3072 × 2048 pixels) were stored on a removable memory card. Images were finally processed by removing black background and manually cropped by using in-house-modified scripts available in the Matlab Image Processing Toolbox (Matlab, R2020b, The Mathworks, Sherborn, MA, USA).

Both digital pictures and MRI images were reordered in a concatenated folder structure according to the scheme illustrated in Figure S4. Each subfolder, related to the two ripening stages for each season, was separately compressed in zip format as input data for the deep learning classification described in the following section.

#### *2.7. Deep Transfer Learning Based Classification*

The proposed classification method relies on the deep transfer learning approach based on deep neural network (DNN) architectures. DNNs consist of a series of functions that take an input and return a predicted label. For this purpose, labeled input data were initially split into a training dataset (i.e., the set of data used to fit the model) and a validation dataset (i.e., the set of data used to provide an unbiased evaluation of a model fit on the training dataset while optimizing, or tuning, the model parameters). In brief, labeled input data pass through the DNN network using a loss function to evaluate how well the model performs at correctly identifying the true classes. The model optimizes for this loss function by computing the gradient of the loss function with respect to the model parameters. In parallel, the validation process optimizes the model parameters iteratively to minimize the loss and perform validation of the trained model [31]. Six pre-trained DNN frameworks were tested: ResNet [32], AlexNet [33], VggNet [34], SqueezeNet [35], DenseNet [36] and Inception-Net [37]. All models were implemented using the learning framework Pytorch, written in Python v. 3.7 and integrated into the open source GPU web server Google Colaboratory [38]. For the training step, all six DNN frameworks were trained under the following conditions: number of epoch = 1000; batch size = 8; learning rate = 0.001; momentum = 0.9. When the loss function converged and stabilized, the training step was automatically stopped and the training model saved. For the validation set, 20% of the training dataset data was used, and a 5-fold cross validation was performed to evaluate the performance of each model. Each model was trained 5 times on each dataset, and related performance parameters (training and validation accuracy and loss as well as training time) were computed. Finally, for all pre-trained models, the average values of the five replicates were calculated and listed in tables in order to compare classification performance.

#### **3. Results and Discussion**

#### *3.1. MRI Analysis of Dataset 1*

NMR relaxation of protons can occur only if stimulated by a fluctuating field of proper frequency that would induce the spin transition, in a non-spontaneous relaxation process that finally leads to equilibrium magnetization. Longitudinal relaxation, defined by the constant T1, is triggered by field fluctuations at the Larmor frequency. Contrarily, in solid or semi-solid systems such as cheese, the major contribution of transverse relaxation, described by T2, comes from the slow molecular dynamics of the studied liquid at zero Larmor frequency. For this reason, observed T1 values usually present a very broad distribution, while populations are sharper and generally better defined for T2. We chose the latter NMR relaxometry constant to describe the analysed cheeses.

Informative relaxometry features of the samples can be derived by analysing the multiexponential decay of the transverse magnetization from CPMG experiments at variable interpulse spacing [19,27,39]. Fitting of CPMG decays of MRI images of FS cheeses indicated a multiexponential behaviour (i.e., two T2 populations), deriving from the presence of two water proton pools (Figure 1).

Multiexponential decay of the transverse component of NMR magnetization arises from the combined effect of diffusive and chemical exchange phenomena. In this condition, water molecules in cheese could be described, at a first approximation, as either free water or hydration water molecules (i.e., water embedded in the protein matrix and strongly interacting with it). Of course, the term "free" could be misleading, since water dynamics in this fraction are in fact much slower than in bulk water. NMR relaxometry is able to clearly detect the presence of different water pools, since water molecules sample different domains within the time scale of the CPMG experiment (i.e., the interpulse spacing). When diffusion is slow (and/or when relaxation domains are large), water molecules in the free water domain and hydration water retain their characteristic relaxation times, and multiexponential distribution is observed (e.g., in aged cheese); when diffusion is fast (and/or heterogeneity of domains is small) (e.g., in protein solutions, milk and curd), relaxation of water in different compartments is averaged, and proton transverse relaxation is mainly monoexponential [24,39,40].

**Figure 1.** Representative T2 distribution profile of an aged FS cheese at 300 MHz. Two proton populations are observed, characterized by geometric mean T2 values (T21 and T22) and area fractions of each population (AF1 and AF2).

At high field, diffusive exchange of water and fast chemical exchange between water protons, labile protein protons, and exchangeable protons in small molecules (such as sugars, vitamins, etc.) play a significant role in characterizing the observed relaxation distributions (Figure 1).

It was shown that in skim milk chemical exchange is mainly modulated by whey proteins [40], but in whole milk, the role of fat globules should be taken into account, since fat influences NMR relaxometry results [24]. In addition, it is clear that in cheese fat, protons might reasonably have a role in defining the observed NMR relaxometry profiles [41]. However, in the following discussion, the effect of fat protons on the observed T2 distributions can be reasonably neglected, since fat signal has been saturated by a solvent suppression pulse in the MRI sequence, as previously explained [19].

For each pool, the mean area fractions (AF1 and AF2) and the corresponding T2 values (i.e., the maximum intensity of each peak, T21 and T22) were obtained (Table 1).

These data are in agreement with previous studies on aged ewe milk cheese and FS PDO [19], which also reported a first T2 population relaxing in the range of 7–12 ms and a second one at 45–53 ms. The main novelty with respect to the mentioned data consists of the fact that the cheeses analysed in the present work were all made by following the official cheesemaking protocol of FS PDO, except for the thermal treatment of milk. Furthermore, a larger amount of data were used here to preliminarily validate observed data trends; the present results were indeed acquired during a much larger lag time (two production seasons), and two ripening times (105 and 180 days) were compared. Moreover, products from five different artisanal producers plus industrial Fiore Sardo PDO cheeses were analysed.

*Dairy* **2021**, *2*

**Table 1.** T21 and T22 values (ms) and corresponding mean area fractions AF1 and AF2 (%) of artisanal FS cheeses produced in Season 1 and Season 2; asterisks (\*) in each row indicatesignificant differences (*t*-test, *p* < 0.05) between RC and HTC parameters (T21 and T22, AF1 and AF2) for each sample (M = Mean; SD = Standard Deviation).


Statistical analysis of MRI relaxometry values (Table 1) indicated significant differences in T21 and T22 values between RC and HTC cheeses. Evaluation of AF suggested the same consistent trend, with a significantly larger amount of protons in the faster relaxing population (AF1) in HTC cheeses than in RC. This is in quite good agreement with previous results on MRI relaxometry of heat-treated milk cheeses [19] and on the sensorial evaluation of Fiore Sardo cheese, suggesting that sensorial attributes and relaxometric features may correlate to each other and that they are both more appreciably affected by milk properties and cheesemaking processes than ripening [11].

Since for each production reported in Table 1 (S1–S5) the same milk was used, the manufacturing process was carried out by the same dairyman, and ripening time was exactly the same for both treatments (RC and HTC), the only factor that accounts for the observed relaxometry differences is milk thermisation.

From a previous investigation [19] and in accordance with pertinent literature [42–44], the fast relaxing proton pool (described by the pair T21; AF1) can be associated with a fraction of water molecules strongly interacting with protons and entrapped in the para-casein network (protein hydration water); the slower population (T22; AF2) can be associated with more mobile water molecules, which less strongly interact with the cheese protein network (which is sometimes referred to as bulk or "free" water).

Several different mechanisms may reasonably affect observed transverse relaxation rates. The first mechanism is related to an altered proportion between more mobile (free) water pools and water in the vicinity of the protein surface (hydration water), which is reasonably well described by AF1 [16,19,45]. This mechanism is in turn modulated by the diffusive exchange rate between water entrapped within the protein network of cheese, confined in cheese holes or pores and hydration water, and chemical exchange [39,46] between water protons and labile protein protons (mainly whey proteins [40]). The second mechanism involves internal water molecules in casein micelles, sometimes referred to as "bound" water molecules, which represent structural water molecules in the internal cavities, engaged in slowly modulated intermolecular dipole couplings with protein protons [46]. A significant effect of these "bound" water molecules on the observed proton transverse relaxation rate has been already suggested in dairy systems [40]. A third mechanism is explained by considering changes in the correlation time of macromolecules as a function of the ripening and processing methods [16,47]. Finally, relevant contribution could be expected from the mineral equilibrium of milk and its effect on micelles hydration. Equilibrium between colloidal (solid) and soluble calcium phosphate is in fact affected by temperature [48,49]; this has also been observed in Fiore Sardo cheeses manufactured from raw or heat-treated milk [6]. In particular, higher total calcium was found in Fiore Sardo produced from raw milk than from thermised milk, and it was hypothesized that mineral equilibrium in milk was shifted from soluble calcium to insoluble colloidal calcium as an effect of a high treatment temperature, leading to a higher retention of calcium in Fiore Sardo from raw milk [6]. Micellar calcium phosphate is a highly hydrated colloid, which greatly influences proton transverse relaxation. Ca and P association to casein and micelle hydration are strictly related phenomena, having a marked effect of NMR transverse relaxation of protons [40,48,49]. Heat treatment of milk was in fact found to cause precipitation of calcium phosphate and correspondingly affect proton T2 [48].

#### *3.2. MRI Analysis of Dataset 2*

Industrial Fiore Sardo samples were characterized by the presence of two water proton pools, in agreement with results on artisanal samples RC and HTC discussed above. Mean area fractions (AF1 and AF2) and T2 values (T21 and T22) are shown in Table 2.

Both relaxation time constants T2 and area values AF were consistent with the HTC artisanal cheeses analysed in Dataset 1, with the exception of samples from the maturer industry (I12, I13).

Both the increase of AF1 and the decrease of T21 are the result of a strong interaction between water protons and the para-casein protein network of cheese, i.e., of cheese protein hydration. From a certain point of view, AF1 more directly describes the level of protein hydration, being determined by the amount of protons in this pool. However, a decrease of the T21 value, being the result of the weighted average between the (long) T2 of protons in free water pools and that (very short) of labile protein protons, characterized by a different chemical shift and affected by the typically long reorientational times of biological macromolecules, can be indirectly influenced by the amount of protein hydrating water.

Aiming to find suitable and objective descriptors of the relaxometric changes induced in cheese by milk thermisation, we define a new parameter, namely the Score (Ṩ) factor, as a representative indicator of the state of cheese microstructure. The Score (Ṩ) is defined as follows:

$$\mathbf{S} = \mathbf{T}\_{21} / \text{AF1} \tag{1}$$

where T21 and AF1 are the relaxation time and the area fraction of population 1, respectively, of the MRI T2 distribution profile (Figure 1). Based on what was discussed above, the lower the Score, the more affected cheese microstructure is by milk thermisation. Correspondingly, the higher the Score, the less influenced are cheese proteins by any thermal effect on milk, i.e., milk can be considered as raw or untreated.

The calculated Ṩ values (Season 1 and Season 2) are presented in Figure 2.

Even at first glance, two considerations can be made when observing Figure 2: first, thermisation of milk in HTC cheeses (red symbols) leads to a reduction of Ṩ with respect to their raw RC counterparts (blue symbols); second, RC samples show much higher Ṩ variability than their HTC counterparts.



Cheeses from one producer (S5) showed subtle, though significant, NMR relaxometry changes upon milk thermisation (blue and red circles are superimposed in Figure 2). In fact, for producer S5, both AF1 and T21 values always showed only slightly significant changes when RC and HTC cheese samples were compared at 105 days of ripening and at 180 days in Season 1 (Table 1). It should be considered that, at the shortest ripening time allowed for Fiore Sardo commercialization (105 days of ripening), the effect of water redistribution and protein hydration might still be not sufficient in order to clearly discriminate RC from HTC. In these cases, NMR relaxometry analysis should be replicated after a longer ripening time to obtain clearer results (Figure 2b,d). If a discrimination between RC and HTC is still not present at six months of ageing, then very likely milk has been overheated to some extent,

and cheese features are not comparable to raw milk cheeses anymore. This result can be explained by assuming that, though thermisation was not intentionally performed in S5 production, raw milk was likely subjected to excessive heat before rennet addition. In fact, it is worth recalling here that, even at an artisanal level, milk is always treated to some extent during cheesemaking, since rennet is usually added when milk reaches, upon heating, a temperature in the range of approximately 35 and 38 ◦C, depending on the season, site, climate and usual practices of production. Among Fiore Sardo artisan producers, milk heating, according to traditional cheesemaking processes, is mostly carried out by directly warming up (often by means of an industrial stove gas cooker) a copper boiler containing raw milk. In some other cases, but still at artisanal level, more modern stainless steel multipurpose cheesemaking tanks are used to heat milk by means of water vapour in somewhat larger dairy plants. Only bigger industrial dairy plants adopt flow pasteurizers equipped with sophisticated plate heat exchanger modules, in order to perform more controlled high temperature–short time heat treatments. This consideration could explain the very subtle difference between the NMR relaxometry parameters of RC and HTC cheeses at 105 days and 180 days of ripening in S5 samples.

**Figure 2.** Score factors (Ṩ) of Fiore Sardo cheeses from Dataset 1: 105 days (**a**) and 180 days (**b**) of ripening, Season 1 (March–April 2019); 105 days (**c**) and 180 days (**d**) of ripening, Season 2 (January– February 2020). Blue circles represent samples made from raw milk (RC) and red circles represent cheeses made starting from heat-treated milk (HTC).

Notably, S5 and other shepherd's productions show superimposed error bars in Figure 2. This is mainly due to the large variability of RC and to similar technological considerations as already discussed for S5.

As explained above, the score value Ṩ represents and describes the status of water molecules more tightly interacting with the paracasein network of cheese. Low score values (i.e., higher AF1 and/or lower T21) indicate a high degree of protein hydration, arising from either ripening [16,45] or from the effect of heat treatments [19]. For cheeses having the same ripening time, the effect of the thermal treatment of milk on the following processing steps [5] dominates relaxation.

In fact, previous studies carried out at variable fields by fast field cycling NMR relaxometry (FFC-NMR) have demonstrated that the hydration of cheese proteins is affected by cheese ageing [45,47]. Further FFC-NMR studies have demonstrated that, for a given ageing time, cheese protein hydration can also be significantly affected by the different processing methods preceding cheese maturation [16]. As discussed above, two similar Fiore Sardo cheeses, which only differ in their heat treatment of milk, also differ in their NMR relaxometry characteristics. The reasons for the observed relaxometry differences can be found in cheese microstructure (such as fractal dimension, protein proton to water proton ratio, water mobility and protein hydration) [16].

The variability of Score values (Ṩ) within each group and among groups of cheeses is presented in the form of box plots in Figure 3. Standard deviations of RC of shepherds (blue boxplot) were generally larger than the other groups. Median values of RC are higher than all other groups, followed by maturer cheeses, HTC and industrial samples, with the latter showing the lowest median values. RC samples show considerably higher interquartile range than HTC and much larger whiskers. This could be reasonably ascribed to both the effect of the artisanal manufacturing process on cheese microstructure and to the microbiological complexity of the systems arising from the use of raw milk. In fact, while artisanal cheesemaking entails mainly manual practices, industrial production implies more standardized protocols and a higher degree of process automation. In addition, previous reports on raw and heat-treated milk cheeses demonstrated a larger variability of raw-milk cheeses than pasteurized-milk cheese counterparts, and it was suggested that the effect of heat on milk proteins, indigenous enzymes and microbial flora of milk may reduce variability, resulting in a more homogeneous peptide profile, sensory properties and NMR relaxation parameters for pasteurized-milk cheese [19,50,51].

**Figure 3.** Box plots of the Score values (Ṩ) of each group of cheeses (RC, HTC, industrial and maturer cheeses).

Differences in the manufacturing process are reflected in different physico-chemical and sensory features, resulting in the distinctiveness of each wheel from the same shepherd. While heating milk to 35–38 ◦C, necessary for a proper action of rennet, is usually not strictly standardized in terms of temperature control, industrial processing is carefully controlled. Mild thermal treatments of raw milk, typical of artisanal practices, allow for the

preservation of the microorganisms that represent a sort of fingerprint of each shepherd's production. The development of indigenous milk microflora and the related molecular changes (e.g., changes in pH, proteolysis and lipolysis during ripening, etc. [5]) occurring during the following processing steps, together with the exclusively manual manufacturing operations on curd, are reflected in the microstructural, molecular and sensory properties of the final cheese, making every wheel a unique product. That is why artisanal Fiore Sardo PDO cheeses are very different from each other, being strongly linked to the territory and to the specific production practices. It is interesting to note, in this regard, that milk and curd handling practices for making Fiore Sardo are passed down from generation to generation and often differ from each other in subtle, but likely relevant, details. On the other hand, when heat is applied to the same milk, though (almost) the same manufacturing practices are carried out by the same cheesemongers, part of the characteristic microflora is depleted and milk functionality is altered [51,52]. The Ṩ data distribution in HTC cheeses indeed indicated a reduced heterogeneity in comparison to RC. The box plot of industrial cheeses more clearly shows that the spreading of Ṩ values is considerably reduced and very narrow in comparison to both RC and HTC. This is mostly ascribable to industrial practices, where cheesemaking is strictly controlled and standardized in terms of both microbiological conditions and technological operations. It is worth noting that industrial Fiore Sardo samples were purchased from different industrial producers, so the observed low interquartile range should be considered quite representative of the standardization of the industrial production in the market and cannot be ascribed to a specific producer. Mean and standard deviations of Score values Ṩ for all groups of samples (RC, HTC, I

and maturer), together with their statistical significances, are presented in Table 3.


**Table 3.** Statistical comparison of Score (Ṩ) values (ANOVA, *p* < 0.05) of the four groups of Fiore Sardo samples, namely RC, HTC (produced in Season 1 and Season 2), maturer and industrial cheeses. Each shepherd production (RC and HTC) was individually compared to maturer and industrial cheeses; data are expressed as Mean ± Standard Deviation; means sharing the same superscript letters are not significantly different from each other; the letter "a" is assigned to the higher mean.

In general, ANOVA results suggested a consistent differentiation among the four groups of cheeses. Industrial Fiore Sardo showed the lowest Ṩ values in all comparisons, except for one case (i.e., S3 180 days, Season 1), which had Ṩ values higher than HTC only. RC often had significantly higher Ṩ values than other groups of cheeses, except for two cases (i.e., S4 and S5, 105 days, Season 1), when maturer cheeses showed higher Ṩ. This result can be explained by considering that maturer cheeses were purchased after approximately 6–8 months of ripening. In some other cases, Ṩ values were comparable for RC, HTC and maturer cheeses. Statistical comparison of HTC and industrial cheeses usually resulted in a significantly lower Ṩ for the latter, but they were comparable in five cases. Maturer samples were never comparable to industrial samples.

Based on the aforementioned statistical evaluations, some conclusions can be drawn on the manufacturing process of each group of samples and on the quality features of the resulting cheeses. RC and HTC were entirely produced following artisanal procedures. Maturer cheeses share the same artisanal processing but were ripened in industrial-like conditions, in refrigerated cellars with controlled humidity and temperature. Despite different ripening conditions, RC and maturer cheeses were described by very similar Ṩ values, clearly suggesting that technological transformations that precede ripening confer to cheeses a common relaxometric behaviour described by comparable Ṩ. RC and maturer cheeses were generally described by higher Ṩ than industrial cheeses, although ANOVA suggested that RC, HTC and maturer cheeses can still differ to some extent among each other in terms of the new quality-related parameter, Ṩ. Since the RC and HTC cheeses of each shepherd were technologically processed in the same way, the decrease in Ṩ observed between RC and HTC (Figure 3), supported by ANOVA comparisons (Table 3), can be reasonably and solely related to the application of a heat treatment of the milk. On the other hand, industrial samples, besides having a much narrower Ṩ variability (Figure 3) clearly exhibited considerably and significantly lower Ṩ values with respect to the other three groups of cheeses. These relaxometric behaviours could be reasonably associated with the features of industrial productions, for which milk and curd processing and handling are more mechanized and standardized, ensuring more uniform features in the final products, as previously reported [16]. Moreover, a contribution to the low Ṩ may arise from the application of a heat treatment to the milk that is consistent with common industrial processing practices, as realistically suggested [12–15]. Interestingly, the discrimination ability of Ṩ and its potential in identifying relevant quality features of Fiore Sardo cheeses does not seem to be affected by different seasons of lactation (i.e., to milk compositional quality), being instead more strictly related to milk processing. This consideration derives from the significant statistical differences among the four groups of samples, regardless the season (Figure 3, Table 3). As far as ripening conditions are concerned, the significant role of ripening cannot be excluded. However, the effect of ripening conditions seems to be less influential than the milk and curd processing practices preceding cheese ageing. This makes NMR relaxometry a suitable analytical tool for characterizing the effect of heat treatments of milk on quality features of Fiore Sardo cheese.

#### *3.3. Moisture Content*

The moisture content of all cheese samples are reported in Table S2.

Moisture content was not influenced by heat treatment, indeed showing no consistent trend in RC and HTC samples. This is in contrast to MRI parameters, which showed consistent changes according to heat treatments. Since MRI parameters T21 and AF1 could be considered representative of the portion of water in the samples, such changes could not be related to total moisture content [19]. In addition, for industrial samples, no noticeable correlation was found between MRI parameters and moisture content.

Water status can be described at different scales of investigation—at a macroscopic level by moisture content and at a molecular level by MRI—and can lead to different information about the analyzed matrix. Moisture content is representative of the extractable water molecules by drying in certain conditions and is among the most common parameters used to describe food quality. However, it should not be considered as an exhaustive parameter as it does not give any information on water proton pools in the system (if any) and their interactions with the protein matrix in cheese [53].

#### *3.4. Image Analysis of Fiore Sardo Pictures and MRI Images*

Representative pictures of the visual appearance of RC and HTC at the minimum allowed ripening time for Fiore Sardo PDO (i.e., 105 days) are shown in Figure 4.

**Figure 4.** Representative images of the visual appearance of RC, HTC and industrial samples: Producer S1, Season 2, 105 days of ripening, Fiore Sardo from raw milk (upper row) and from heat-treated milk (lower row) (**a**); details of RC (**b**) and HTC (**c**) samples ripened 105 days upon cutting; industrial Fiore Sardo purchased from a local market, with at least 105 days of ripening (**d**).

Even at first sight, RC, HTC and industrial cheeses appear very different from each other. The HTC samples surface appears more homogeneous, smoother and with fewer holes, while RC cheeses are more grainy, friable and show more and bigger holes. When cut with a knife, RC cheeses are prone to easily break and crumble, while HTC appear more creamy and can be cut very easily and without breaking into pieces.

The different visual appearances of HTC samples with respect to RC cheeses is reasonably ascribable to the modifications resulting from the heat treatment of milk, which is known to have, in Fiore Sardo PDO, multiple consequences on cheese composition, microstructure and texture, which in turn affect appearance and sensory properties in general [6,10,52].

In general, MRI images confirmed, by a non-invasive approach, what already appeared at a first visual inspection of samples. RC cheeses exhibited a coarser paste, often with multiple cracks and holes, while HTC samples exhibited a smoother and more homogenous paste (Figure S5).

Two kinds of scientific considerations can be made based on the MRI images acquired. On one hand, as discussed above, MRI images give relaxometric information that is able to provide useful insights on the state of biopolymers and on the mobility of the system at molecular level. On the other hand, important information can be derived from the texture of MRI images by performing a detailed image analysis.

In recent years, image analysis techniques rapidly emerged not only at the academic level but also as rapid, economic, consistent and objective inspection techniques in the agricultural and food industry [54,55]. In particular, for the dairy industry sector, modern computer vision and image analysis methods have been used to evaluate cheese external or internal quality attributes such as colour, shape, texture, amount and distribution of added ingredients or production defects [56–58].

In this study, we pave the way for a deep-learning-based method for the classification of Fiore Sardo cheese images. The method consists of using data collected both from MRI and from a computer vision system based on photographic images taken with a smartphone camera. Deep learning is a subfield of artificial intelligence (AI) concerned with specific algorithms for thinking simulations and data processing inspired by the neurons of the brain [31]. Deep-learning-based image classification has proven to be highly effective in image classifying, object detection and other computer vision related problems [59,60].

All pre-trained DNN models run well in the open web GPU\_Google colab service, with compatible and relatively short training times. For both MRI and photo datasets, the best classification result was achieved by the SqueezeNet model (Tables S3 and S4). In particular, for the T2 weighted MRI images, the SqueezeNet model exhibits a training accuracy from 93% to 98%, with the highest validation score (100%) (Table S3). As far as the photo dataset is concerned, SqueezeNet shows a training accuracy of 96%, with the same validation of 100% (Table S4).

The other pre-trained models (ResNet, AlexNet, VggNet, Densenet and Inceptionnet) also exhibit satisfactory performance results, ranging from a minimum of 72% to a maximum of 92%, with related validation ranging from 77% to 100% (Tables S3 and S4).

We propose here an IA approach to classify photos and T2 weighted MRI images of Fiore Sardo cheese made with raw or thermised milk, by applying a deep neural network (DNN), which includes several chained layers of processing between the input (photos or MRI images) and the output (classification/label of the input image).

As with most state-of-the-art IA methods, deep learning models usually require fitting thousands or millions of input data; thus, their application to those studies where the sample size is limited appears prohibitive. To address this challenge, the deep transfer learning approach offers a more suitable alternative. Transfer learning is the method by which the model uses the knowledge gained during the training of a relatively large dataset in a different but related problem (i.e., image classification and recognition). During the transfer learning process, only the selected classifier(s) is trained in the new network, while the features learned from the large dataset are transferred. Therefore, with the deep transfer learning approach, we do not need to retrain the entire network for a new dataset, thus allowing the training of deep learning models with relatively few data, reduced computational power and less training time [61].

With the main disadvantage of the present study being the relatively low number of available cheese samples, compared to the large data size usually needed for training deep learning frameworks, we chose to transfer the learning technique by testing six pre-trained DNN frameworks (ResNet, AlexNet, VggNet, SqueezeNet, DenseNet and Inception-Net), a widely used approach for addressing classification and image recognition challenges of MRI data.

In the end, all tested pre-trained DNNs yielded promising performances, giving an average classification accuracy of over 72% and over 93% in particular for the best performing SqueezeNet model (Tables S3 and S4). We did not perform further accuracy tests of the DNN pre-trained model in this study; however, we are proposing to implement and optimize this model in a dedicated study when more data are available.

#### **4. Conclusions**

Our data demonstrate that NMR relaxometry is able to assess molecular mobility changes induced in cheese by the most common industrial thermal treatments carried out on milk. In the present work, the experimental plan was optimized to exclude possible confounding factors. The same milk was used for producing both raw-milk or thermisedmilk cheese. Moreover, analyses were repeated in different seasons, and the products from different producers were compared. Results were compared with a selection of industrial Fiore Sardo cheeses purchased in the local market.

NMR relaxometry data of Fiore Sardo cheese indicated the presence of two water populations, described by characteristic T2 relaxation times and area fractions. In particular, the first and fastest relaxing population showed a good predictive value. In order to obtain a relaxometric descriptor of FS cheese microstructure and molecular mobility, a new qualityrelated parameter was introduced, namely the Score factor (Ṩ). The Ṩ describes cheese in relation to the hydration status of the casein network and the state of the water-fat-protein network in cheese. We were able to discriminate FS cheeses produced from raw (RC) and heat-treated milk (HTC), with the latter exhibiting a lower Ṩ, indicating an enhanced hydration of the casein network. Industrially manufactured FS cheeses had generally lower Ṩ values than RC and HTC. Samples from the maturer industry, which buys cheeses before ripening and ages cheeses in industrial-like refrigerated cellars, showed similar values to artisanal cheeses, suggesting that NMR relaxometry is more sensitive to changes occurring in the cheesemaking steps preceding ripening. A further discriminative approach was finally carried out by processing acquired MRI images and photos with an IA method based on pre-trained deep learning algorithms. Preliminary classification results suggested that this innovative approach is promising in providing suitable discrimination between RC and HTC Fiore Sardo cheese.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/10 .3390/dairy2020023/s1, Figure S1: sample production and sampling scheme for MRI analysis, Figure S2: MRI sample preparation scheme, Figure S3: Computer Vision System (CVS) for photographic analysis, Figure S4: Image Dataset Directory and File Structure used, Figure S5: Representative T2 weighted MRI images of Fiore Sardo cheeses, Table S1: Production of FS cheeses and datasets, Table S2: Moisture content data of cheeses, Table S3: Performance and classification results of pre-trained DL frameworks on MRI images, Table S4: Performance and classification results of pre-trained DL frameworks on photographic images.

**Author Contributions:** Conceptualization, R.A.; Methodology, R.A., E.C. and R.M.; Software, E.C. and R.M; Formal Analysis, R.A., E.C. and R.M; Investigation, R.A., E.C. and R.M.; Data Curation, E.C. and R.M.; Writing—Original Draft Preparation, R.A., E.C. and R.M.; Writing—Review and Editing, R.A., E.C. and R.M.; Visualization, E.C. and R.M.; Supervision, R.A.; Funding Acquisition, R.A.; Project Scientific Coordination, R.A. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Sardinia Regional Government by means of Sardegna Ricerche (art 9, L.R. n. 15, 17.11.2010).

**Data Availability Statement:** The data presented in this study are available on request from the corresponding author.

**Conflicts of Interest:** The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

#### **References**


### *Communication* **GC-MS Metabolomics and Antifungal Characteristics of Autochthonous** *Lactobacillus* **Strains**

**Paola Scano 1,\*, M. Barbara Pisano 2, Antonio Murgia 1, Sofia Cosentino <sup>2</sup> and Pierluigi Caboni <sup>1</sup>**


**\*** Correspondence: scano@unica.it; Tel.: +39-070-675-4391

**Abstract:** *Lactobacillus* strains with the potential of protecting fresh dairy products from spoilage were studied. Metabolism and antifungal activity of different *L. plantarum*, *L. brevis*, and *L. sakei* strains, isolated from Sardinian dairy and meat products, were assessed. The metabolite fingerprint of each strain was obtained by GC-MS and data submitted to multivariate statistical analysis. The discriminant analysis correctly classified samples to the *Lactobacillus* species and indicated that, with respect to the other species, *L. plantarum* had higher levels of organic acids, while *L. brevis* and *L. sakei* showed higher levels of sugars than *L. plantarum*. Partial Least Square (PLS) regression correlated the GC-MS metabolites to the antifungal activity (*p* < 0.05) of *Lactobacillus* strains and indicated that organic acids and oleamide are positively related with this ability. Some of the metabolites identified in this study have been reported to possess health promoting proprieties. These overall results suggest that the GC-MS-based metabolomic approach is a useful tool for the characterization of *Lactobacillus* strains as biopreservatives.

**Keywords:** *Lactobacillus*; antifungal activity; fresh cheese; biopreservation; GC-MS; metabolomics

#### **1. Introduction**

In the Mediterranean area, to revive the ovine dairy industry, the development of new fresh dairy products is increasing. Taking into account that the physico-chemical characteristics of these products greatly favor the growth of pathogenic and spoilage microbes, many research activities have been devoted to finding solutions to protect fresh dairy products from deterioration and to prolong their shelf life while preserving the organoleptic and nutritional properties. Furthermore, in recent years, the increased request for natural foods has challenged the food producers and researchers to find natural alternatives to synthetic preservatives. Among natural food preservation techniques, biopreservation has gained particular attention; this term refers to the use of microbial strains or metabolites able to inhibit the growth of spoilage and pathogenic microbes in foods and thus improving safety and extending their shelf-life [1]. Lactic acid bacteria (LAB), naturally present in many fermented and not fermented foods of animal and vegetable origin, are considered good candidates as bio-preservative, as they possess numerous functional properties related to the improvement of food safety and health benefits [2,3]. Among LAB, several strains belonging to the genus *Lactobacillus*, frequently present in milk and cheese as the prevalent nonstarter LAB, have shown distinct antimicrobial and antifungal activities [4]. Additionally, *Lactobacillus sakei* strains, characteristic of meat products, have shown antifungal activity [4]. The latter strains can be found in raw meat products stored under vacuum and refrigerated, as well as in fermented sausages, due to their metabolic properties and phenotypic features that are particularly well adapted to growth and survival under the conditions found during meat processing and storage [5,6].

**Citation:** Scano, P.; Pisano, M.B.; Murgia, A.; Cosentino, S.; Caboni, P. GC-MS Metabolomics and Antifungal Characteristics of Autochthonous *Lactobacillus* Strains. *Dairy* **2021**, *2*, 326–335. https://doi.org/10.3390/ dairy2030026

Academic Editor: Burim Ametaj

Received: 11 May 2021 Accepted: 18 June 2021 Published: 23 June 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

The antifungal activity of LAB may be related to the synergistic action of several compounds, e.g., organic acids (acetic acid, lactic acid, propionic acid, and phenyllactic acid), hydrogen peroxide, cyclic dipeptides, proteinaceous compounds (bacteriocins), and fatty acids [2,7], and it is recognized that the capacity to synthetize these compounds is a strain-linked trait [4]. In our previous work, we found that selected *Lactobacillus* strains of food origin could control mold contamination without altering the sensory characteristics in miniature Caciotta ovine cheese, produced in a laboratory scale [4]; one limitation of this study is that the compounds responsible for the antifungal activity were not studied. Metabolomics is one of the most valuable techniques for the study of metabolite profiles in food matrices [8], and gas chromatography coupled with mass spectrometry (GC-MS) is a well-suited analytical platform for the study of microbial metabolites [9,10]. By the combined approach of GC-MS and multivariate statistical data analysis it was possible to find links between metabolites of dairy products and their microbial profiles [11]. The aim of this work was to test the suitability of this approach for the identification of the metabolites involved in the bioprotective characteristics of autochthonous *Lactobacillus* strains isolated from Sardinian dairy and meat products and to seek correlations between metabolites (or group of metabolites) and antifungal activity. To this goal, we studied the polar and semi-polar low molecular weight metabolites present in the cultured broth of different *Lactobacillus* strains. Antifungal activity of these strains against mold species, commonly occurring in cheese spoilage, was also assessed and results correlated to the GC-MS data by multivariate regression analysis.

#### **2. Materials and Methods**

#### *2.1. Microorganisms and Cultivation Conditions*

A total of 9 autochthonous *Lactobacillus* strains (3 *Lactobacillus plantarum*, 2 *Lactobacillus brevis*, and 4 *Lactobacillus sakei*), deposited in the MBDS culture collection (www.mbds.it, accessed on 25 May 2021) were studied. They were isolated from ovine raw milk and artisanal cheeses (*L. plantarum* and *L. brevis*), and sausages (*L. sakei*) manufactured in Sardinia (Table 1), and they were classified on the basis of phenotypic characteristics and molecular methods through polymerase chain reaction amplification. *Lactobacillus* strains were stored at −20 ◦C in DeMan Rogosa Sharpe (MRS) broth (Microbiol, Cagliari, Italy) with 15% (*v*/*v*) glycerol and routinely cultivated on MRS agar plates for 48 h at 30 ◦C in microaerophylia. Each strain was subcultered twice in MRS broth prior to experimental use. Carbohydrate fermentation profiles of the *Lactobacillus* strains was assessed and reported in Table S1.


#### **Table 1.** List of the *Lactobacillus* strains and food origin.

#### *2.2. Extraction of Extracellular Samples*

Cell-free supernatants were obtained from approximately 15 mL of overnight cultures of the *Lactobacillus* strains in MRS broth at 30 ◦C. Then, the cells were separated by centrifugation at 7000× *g*, for 20 min at 4 ◦C and supernatant was taken and filtered through a cellulose acetate membrane filter (0.2 μm pore size) to eliminate residual bacterial cells. Cell-free supernatant was distributed in 1 mL aliquots, then dH2O (10 mL) and internal standard (0.2 μmol of 10 mmol solution of 2,3,3,3-d4 D,L-alanine) were added to each sample. Three replicate samples were prepared for each *Lactobacillus* strain cell-free supernatant. Un-inoculated MRS broth samples were collected as control. An aliquot of 100 μL of each sample was stored at −80 ◦C.

#### *2.3. In Vitro Antifungal Activity of Lactobacillus*

Antifungal activity of *Lactobacillus* strains against *Alternaria alternata* (UNICAM2), *Cladosporium herbarum* (UNICAM10), *Paecilomyces variotii* (UNICAM17), and *Penicillium chrysogenum* (ATCC 9179) indicator strains was assessed in vitro using the agar plate method previously described [4]. The inhibitory activity against the fungal strains *Aspergillus flavus* (ATCC 46283), *Fusarium oxysporum* (UNICAM12), and *Mucor recurvus* (UNICAM15), which develop faster than lactobacilli, was assessed using the dual-culture overlay assay [4]. For all assays, the antifungal activity of each strain was determined by measuring the width of the halo around the bacterial streaks, according to the following 4 steps in a semiquantitative scale: (1) ≥8 mm, (2) 5–7 mm, (3) 3–4 mm, and (4) <3 mm.

#### *2.4. Sample Preparation for GC-MS Analysis*

Samples were thawed on ice and subsequently 375 μL of a methanol and chloroform mixture (2/1, *v*/*v*) was added. After 1 h, 380 μL of chloroform and 90 μL of aqueous KCl 0.2 M were added. Samples were vortexed for 10 s and centrifuged at 15,294× *g* for 10 min at 4 ◦C. A volume of 200 μL of the aqueous fraction was taken and transferred into 1.5 mL sterile glass vials with 10 μL of a 80 mg/L of 2,2,3,3-d4-succinic acid solution. The aqueous fraction was then dried by a nitrogen stream, and derivatized with 40 μL of MSTFA (N-Trimethylsilyl-N-methyl trifluoroacetamide). Samples were kept for 30 min at 70 ◦C, before adding 600 μL of hexane and then vortexed.

#### *2.5. GC-MS Analysis*

Samples were analyzed with a Hewlett Packard 6850 Gas Chromatograph, 5973 mass selective detector, and 7683B series injector (Agilent Technologies, Palo Alto, CA, USA) with helium as carrier gas at a flow of 1.0 mL/min. One microliter of each sample was injected with 1 min of split flow delay and resolved on a 30 m × 0.25 mm × 0.25 μm DB-5MS column (Agilent Technologies, Palo Alto, CA, USA). Inlet, interface, and ion source temperatures were 250, 250 and 230 ◦C, respectively. Oven starting and final temperatures were set at 50 and 230 ◦C, respectively, with a rate of 5 ◦C/min for 36 min and then for 2 min at a constant temperature. Electron impact mass spectra were recorded from *m*/*z* 50 to 550 at 70 eV. Metabolite annotation was achieved by mass spectra comparison with analytical standards, in house library and the NIST14 database (National Institute of Standards and Technology, Gaithersburg, MD, USA).

#### *2.6. Multivariate Statistical Data Analysis*

As input for multivariate statistical data analysis (MVA), a 27 × 59 matrix, composed of samples × the intensities of the most abundant fragment ion for each metabolite, was constructed. MVA were performed using the software SIMCA-P+ (version 14.1, Umetrics, Umeå, Sweden). GC-MS variables were mean centered and scaled to unit variance. The principal component analysis (PCA) was carried out to study the sample distributions and presence of outliers. For classification of samples to the three classes of *Lactobacillus*, and to find discriminant metabolites, the supervised partial least squares-discriminant analysis (PLS-DA) was used. Importance of the variables in PLS-DA was assessed by considering the

variable importance in the projection (VIP) values. Results of PLS-DA are shown as scatter plot of scores and loadings in the first two principal components, and as a list of metabolites having VIP > 1 attributed to each class of samples on the basis of their coefficients. In the Coefficients Overview Plot the X-variable coefficients, together with their standard errors in cross-validation, for the Y-class are shown. The supervised partial least squares (PLS) regression was applied for the study of linear relationships between antifungal activity and metabolite profiles. Antifungal activity scores were organized in the Y matrix. Results are shown as PLS loading weights plot in first and second components that displays the relation between the X-variables (GC-MS data) and the Y-variables (antifungal activity). Metabolites positively correlated with antifungal high activity lie next to the Y-variable, metabolites negatively correlated lie in the opposite side of the plot. The quality of the PLS and PLS-DA models and the proper number of principal components were assessed on the basis of the cumulative parameters R2Y and Q2Y (prediction power calculated in cross-validation), and of their difference value, with a <0.50 threshold.

#### **3. Results**

#### *3.1. GC-MS Metabolite Profiles of Lactobacillus Species*

The hydrophilic extracts of samples were analyzed by GC-MS and 59 low molecular weight metabolites were annotated and reported in Table S2. Short-chain carboxylic acids (such as lactic acid, citric acid, succinic acid, and phenyllactic acid), amino acids (valine, leucine, isoleucine, threonine, and pyroglutamic acid), mono- and disaccharides, and polyols were annotated together with fatty acids and analogues, such as oleamide. Analysis of MRS broth chromatograms indicates that it was mainly composed by monoand disaccharides, citric acid, phosphate, and amino acids essential for growth (results not shown).

To assess differences among the three species of *Lactobacillus*, a PLS-DA was carried out. The discriminant analysis, with a high degree of confidence, well classified samples (score plot in Figure 1a), indicating that the 3 *Lactobacillus* species showed overall different metabolite profiles. By analysis of the loading plot (Figure 1b) and metabolite VIP scores reported in Table 2 (coefficients of metabolites are shown in Figure S1), the fermentation medium of *L plantarum* was found to contain more lactic acid, 3-hydroxy butyric acid (BHBA), 2-hydroxy isovaleric acid (AHVA), 2-hydroxy isocaproic acid (AHCA), succinic acid, 3-phenyllactic acid, and malic acid, together with lower levels of mono- and disaccharides than *L. brevis* and *L. sakei*. In contrast, *L. sakei* showed more sugars. When compared to the other two *Lactobacillus* species, *L. brevis* had more monosaccharides, AHVA, and 4-gamma-aminobutyric acid (GABA) and *L. sakei* showed higher levels of disaccharides. *L. plantarum* and *L. sakei* showed higher levels of pyroglutamic acid. The higher level of saccharides in *L. sakei* and *L. brevis* than in *L. plantarum* was confirmed by the measured carbohydrates fermentation profiles reported in Table S1, where the *L. plantarum* strains, compared to the other 2 species, showed the better fermentation activity.

Broth sugar consumption from *Lactobacillus* strains and production of organic acids are confirmed by Pearson correlation values between GC-MS data depicted as a heat map in Figure S2. In this map, it is clearly visible that organic acids were strongly positively correlated with themselves (r > 0.75) and negatively correlated (r < 0.75) with saccharides (the linear strong negative relationship r = −0.97 between glucose and lactic acid is clearly visible in the Figure S3).

**Figure 1.** PLS-DA of GC-MS data, R2Y = 0.93 and Q2Y = 0.88, over 2 validated components. (**a**) score plot (*Lb*, *Lp*, and *Ls = L. brevis*, *L. plantaris*, and *L. sakei*, respectively); (**b**) loading plot, red hexagons = organic acids; green circles = amino acids; light blue stars = saccharides and polyols; yellow circles = fatty acids and analogues; black diamond = loadings of *Lactobacillus* classes. Unknown metabolites are not displayed. Metabolites are abbreviated as in Table S2.


**Table 2.** PLS-DA VIP <sup>a</sup> scores of discriminant GC-MS metabolites.

<sup>a</sup> Variable importance in the projection, R2Y = 0.93 and Q2Y = 0.88 over 2 validated components; <sup>b</sup> *Lactobacillus* species with higher level of the metabolite, as indicated by the coefficients; <sup>c</sup> OA = organic acids, AA = amino acids and analogues, S = saccharides and polyols, FA = fatty acids and analogues.

#### *3.2. Metabolite Profiles of Lactobacillus Strains*

A visual analysis of chromatograms suggests that all the strains consumed sugars and produced lactic acid, though to a different extent, and synthetized oleamide. *L. plantarum* 1/14537 and *L. brevis* DSM 32516 consumed almost all the citric acid of broth origin (data not shown). To study the different production/consumption of metabolites by the nine *Lactobacillus* strains, metabolites indicated by the PLS-DA as the most discriminant between the three species, together with oleamide, were individually measured and results reported as column plot with means and standard deviations (Figure 2). Phenyllactic acid, AHCA, BHBA, AHVA, and malic acid were produced at a greater extent by *L. plantarum*. Pyroglutamic acid was mainly produced by all the *L. plantarum* strains and *L. sakei* S3/1. *L. brevis* M8/1 produced the greatest, by far, quantity of GABA, confirming that different *Lactobacillus* strains within the same species can have their own individual metabolism.

**Figure 2.** Relative abundance of metabolites in each strain: green *L. plantarum* (1 = 1/14537, 2 = 4/16898, 3 = C1/15); light blue *L. brevis* (4 = DSM 32516, 5 = M8/1 S4); yellow *L. sakei* (6 = S4, 7 = S3, 8 = S3/1, 9 = S5). Y-axis values are in A.U. and refer to the intensity of a selected *m*/*z* ion fragment. Abbreviation of metabolites as in Table S2.

#### *3.3. Antifungal Activity and Correlations with Metabolite Profiles*

Antifungal activity of *Lactobacillus* strains was also studied. Seven mold species were chosen for their widespread presence in cheese spoilage and for their ability to produce mycotoxins [4]. As shown in Figure 3, the strains exhibited a wide range of antifungal activity, dependent on both fungal species and *Lactobacillus* strain. All *L. plantarum* strains and *L. brevis* DSM 32516 showed the strongest activity with inhibition zones greater than 8 mm against all the fungal strains. The *L. sakei* strains were the less active, with S4 strain showing inhibition zones less than 3 mm against *M. recurvus*, *P. variotii*, and *P. chrysogenum*. *A. alternata* and *C. herbarum* growth was strongly affected (inhibition zone higher than 8 mm) by all the *Lactobacillus* strains except for *L. sakei* S5 which had an inhibition zone lower than 8 mm.

**Figure 3.** Antifungal activity of *Lactobacillus* strains: inhibition zones against *A. alternata*, *C. herbarum*, *M. recurvus*, *P. variotii*, *F. oxysporum*, *A. flavus* ATCC 46283, *P. chrysogenum* ATCC 9179.

These results are in agreement with a previous study [4] where single *L. plantarum* C1/15 strain or multiple *Lactobacillus* strains (*L. plantarum* 4/16898, 1/14537 and *L. brevis* DSM 32516) were able to significantly inhibit the growth of *P. chrysogenum* and *A. flavus* when used as adjunct cultures in the production of miniature Caciotta cheese. With respect to the single *L. plantarum* strain (C1/15), the multiple *Lactobacillus* strain combination resulted more active against the fungal strains tested confirming the synergistic action of different antimicrobial compounds highlighted in the current study [4]. To correlate the metabolite profile of strains to the antifungal activity a PLS regression of GC-MS data was carried out considering the extent of the antifungal activity as Y variable. The resulting loading plot in the first two principal components with superimposed scores of antifungal activity is shown in Figure 4.

In this plot, metabolites that lie next to the antifungal activity have positive relationship with this latter, metabolites in the opposite side of the plot, passing from the origin are inversely related with the activity. On the basis of these rules, we can note that organic acids and arabitol were directly involved in the activity against *P. chrysogenum*, *P. variotii*, and *M. recurvus*. Oleamide was correlated with *A. flavus* and *F. oxysporum*.

**Figure 4.** PLS variable loadings plot of GC-MS data (X-variables) superimposed with antifungal activity (Y-variables). R2Y = 0.89 and Q2Y = 0.65 over 5 validated components. Black triangles = Y coefficients for antifungal activity; red hexagons = organic acids; green circles = amino acids; light blue stars = saccharides and polyols; yellow circles = fatty acids and analogues; black diamond = loadings of *Lactobacillus* classes. Unknown metabolites are not displayed. Metabolites are abbreviated as in Table S2.

#### **4. Discussion**

The inherent ability of *Lactobacillus* to produce compounds antagonistic to food microbial contaminants as a part of their natural defense mechanism can provide the dairy food industry with natural and healthier alternatives to chemical preservatives [15], especially for fresh products. The antifungal activity of *Lactobacillus* has been related to the synergistic effect of several compounds, and part of these compounds shows nutraceutical properties.

Compared to the other species here studied, *L. plantarum* strains showed more organic acids, among which 3-phenyllactic acid, which shows a broad spectrum of antimicrobial properties, is effective against bacteria and fungi including yeast [16]. In sourdough bread started with *L. plantarum* 21B, the onset of growth of *Aspergillus niger* was delayed 7 days, with respect to bread started with a *L. brevis* strain not producing phenyllactic acid [17]. However, it has been reported that phenyllactic acid is active against fungal species only at mg mL−<sup>1</sup> concentrations, suggesting an overall antifungal effect in association with other factors produced by LAB [17]. 3-Phenyllactic acid is a by-product of phenylalanine metabolism in *Lactobacillus*, where phenylalanine is transaminated to phenylpyruvic acid and further reduced to 3-phenyllactic acid by 2-hydroxy dehydrogenase [15]. GABA is another product of amino acid metabolism, biosynthesized by microorganisms through decarboxylation of glutamate by glutamate decarboxylase [18]. GABA is an important bioactive compound and its production by various microorganisms has been actively explored [18]. Higher levels of GABA were detected solely in the *L. brevis* M8/1 strain. Besides lactic acid, a number of 2-hydroxy organic acids, that show antimicrobial properties [10], are generally produced by *Lactobacillus*. Organic acids levels were found correlated with the antifungal activity against *P. chrysogenum*, *P. variotii*, and *M. recurves*. *L. plantarum* was found to produce hydroxy organic acids, some of which (AHVA and AHCA) have shown antifungal proprieties [10]. In particular, AHCA has been studied for its potential in inhibiting cell growth and biofilm formation of the pathogenetic *Aspergillus* and *Candida* species in topical use [19,20]. One of the most important virulence factors of *C. albicans* is the ability to form biofilms that protect this yeast against endogenous and exogenous inhibitory substances, in fact almost 65% of microbial infections in humans are biofilmrelated [19]. The effects of oral administration of lactobacilli against biofilm-associated infections have been studied [21]. All the strains synthetized oleamide, not present in the broth, and levels of this compounds have been correlated with the antifungal activity against *A. flavus* and *F. oxysporum*. Oleamide in an amide of oleic acid present in small amount in animal brains and can be produced by microorganisms [22]. In dairy products' fermentation processes, oleamide is synthesized from oleic acid, which is an abundant component in these food products, owing to lipase enzymatic amidation [23]. Epidemiological and clinical data have shown that fermented dairy products possess preventive effects against dementia, including Alzheimer's disease, and those effects have been linked to oleamide [23]. This compound has been identified as the agent responsible for reducing microglial inflammatory responses and neurotoxicity [23]. Oleamide is also an endogenous bioactive signaling molecule that acts in various cell types and could elicit different biological effects [23]. Among the wide range of functions, the most acknowledged one is its sleep-inducing effect [24]. It has been found that *L. plantarum* 1/14537 and *L. brevis* DSM 32516 consumed almost all the citric acid of broth origin. Citric acid fermentation by LAB could lead to acetoin and diacetyl production which are aromatic compounds. However, obligately heterofermentative lactobacilli such as *L. brevis* when present at high levels and facultatively heterofermentative lactobacilli such as *L. plantarum*, under certain conditions, by co-metabolization of citrate and different sugars may produce great amount of CO2, giving rise to holes, splits and blowing defects of aged cheeses [25]. *L. brevis* strains are also characterized by other biochemical properties such as production of CO2 from glucose and ammonia from arginine, as shown in Table S1, that are generally used as tool for lactobacilli classification [26]. Moreover, *L. brevis* could contribute to the degradation of arginine in cheese that can lead to ornithine accumulation by the arginine-urease or the arginine deaminase pathways [27]. On the other hand, ornithine and other free amino acids could be decarboxylated to putrescine and other biogenic amines that in high concentration could be responsible for toxic effect on human health.

In conclusion, results of this study provide reference data for interpreting the differences in the metabolite profiles of *Lactobacillus* strains isolated from different fermented food products and their correlation with antifungal activity in vitro. These autochthonous strains could be considered potential good candidates to be used in cheese manufacturing as bioprotective cultures and for in situ production of postbiotic substances. However, further studies in cheese model systems are warranted to evaluate in situ positive metabolic activities such as organic acids production, GABA and oleamide formation, as well as negative activities including biogenic amines production and the potential to cause blowing defects in the initial and later stages of ripening.

**Supplementary Materials:** The following are available online at https://www.mdpi.com/article/10 .3390/dairy2030026/s1, Figure S1: PLS-DA coefficient overview plot; Figure S2: Heat map of Pearson correlation values of GC-MS data matrix; Figure S3: Correlation plot of lactic acid (X-axis) vs. glucose (Y-axis) levels. Table S1: Carbohydrates fermentation profiles of Lactobacillus strains; Table S2: List of GC-MS metabolites.

**Author Contributions:** Conceptualization, P.C., M.B.P., S.C., and P.S.; Methodology, P.S., A.M., and M.B.P.; Software, P.S.; Validation, P.S. and M.B.P.; Formal Analysis, P.S., A.M., and M.B.P.; Writing— Original Draft Preparation, P.S., P.C. and M.B.P.; Writing—Review and Editing, P.S., P.C., S.C. and M.B.P. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not report data available.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Dairy* Editorial Office E-mail: dairy@mdpi.com www.mdpi.com/journal/dairy

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel: +41 61 683 77 34

www.mdpi.com

ISBN 978-3-0365-4536-3