With the introduction of HTS techniques, several studies addressed the microbiota characterization of the major PDO cheeses, with the aim to improve their quality, safety and commercial values 1,2,19. Scant information is available using HTS approach in traditional mountain cheese production linked to summer pastures biodiversity, including the Italian context12,19. For these reasons in this study, HTS and GC techniques were used for the first time as systemic approaches to characterize the biodiversity of bacterial communities, the volatilome, terpenes and fatty acid profiles in HR cheese.
HTS produced a total of 41,427,924 reads (avg. of 767,184 +/- 240,765 per sample), which, after filtering, resulted in 24,388,097 reads (avg. of 451,631 +/- 143,653 reads per sample). For computational reason, a random subset of 150,000 reads was extracted from each sample. After clustering (97% identity), a total of 550,606 OTUs were obtained, the most of which were singletons (i.e.: OTUs with only 1 supporting read across all samples) that were removed from further analysis, originating a high-confidence dataset of 20,272 OTUs (avg. of 1,360 +/- 463 OTUs per sample).
Proximate analysis
The % of moisture ranged between 34.30 and 37.36 g/100 g, the % of protein ranged between 23.94 and 25.29 g/100 g and the % of fat between 31.11 and 34.33 g/100 g. No significant differences were detected among producers. Some significant differences on the contrary were observed for moisture and fat across different grazing periods for the producers PP1, PP3 and PP4, with an increase of fat from July to September [Supplementary Files, Table S2].
Variability of the microbiota composition of Historic Rebel cheese samples
Historic Rebel cheese-associated microbiota was dominated by members of the Firmicutes phylum which accounted for about 99% of the entire abundance, with a minor presence of Proteobacteria. Carafa et al.12 and Dalmasso et al.6 noticed the same pattern in traditional mountain cheese produced with no starter strains in other alpine areas.
Generally, the overall bacterial composition of HR cheese showed a low diversity, with only 5 genera constituting the vast majority of the community: Streptococcus (mean relative abundance: 76.4%), Lactobacillus (mean rel. ab.: 10.4%), Lactococcus (mean rel. ab.: 6.0%), Leuconostoc (mean rel. ab.: 2.3%) and Pediococcus (mean rel. ab.: 2.1%) together represented an average of 97.2% (range among samples: 87.2% - 99.6%) of the total abundance. Only a minority of samples showed presence of Escherichia (7 samples, max. abundance 7.1%) and Enterococcus (11 samples, max. abundance 3.9%) [Supplementary Figure 1]. As a general trend, triplicated samples per each GP- PP were correlated, as evident from unweighted Unifrac distances analysis, with showed a significant lower distances among replicated samples of the same PP-GP than that from samples of different PPs-GPs (P-value of Mann-Whitney U-test: <0.001) [Supplementary Figure 2].
Analysis of the microbial community among GP-PPs, considering all time points and replicates together (n=9), revealed how the cheese from each producer had a peculiar bacterial composition. They differed both on biodiversity (Kruskal-Wallis P-value<0.01 on Shannon, observed species and Faith’s phylogenetic diversity metrics), with producer PP4 having the most biodiverse community and PP3 the least, and composition (P-value adonis test: 0.001, on both weighted and unweighted Unifrac distances)[Figure 1]; in particular, similar bacterial profiles were observed for PP2 and PP4, and for PP6 and PP3 (pairwise adonis test P-value >0.05), whereas PP5 and PP1 appeared to have a peculiar composition (pairwise adonis test P-value <0.05 in all comparisons). At genus level, HR cheese samples from PP2 were characterized by a somehow higher relative abundance of Lactococcus (15.3% on average) and Leuconostoc (5.4%); PP4 was characterized by a high rel. ab. of Lactococcus (15.4%) as well as the lowest presence (between producers) of Streptococcus (60.9%); PP5 microbiota, instead, showed the highest rel. ab. of Lactobacillus (19.1%) and a subdominant presence of Lactococcus (3.0%). On the other hand, PP6, PP3 and PP1 samples were all characterized by higher percentages of Streptococcus (>84.5%, whereas PP2, PP4 and PP5 were all <75%), associated to Pediococcus in PP6 (5.3%), Leuconostoc in PP3 (3.3%) and Lactobacillus in PP1 (9.8%) [Figure 2]. On the other hand, comparing GPs, we observed a partial difference in alpha-diversity (Kruskal-Wallis P-value<0.01 on Shannon, observed species and Faith’s phylogenetic diversity metrics), due to a lower biodiversity in GPII; however, no difference was visible in composition (P-value adonis test >0.05, on both weighted and unweighted Unifrac distances, Supplementary Figure 2).
Species-level analysis
Since Streptococcus, Lactobacillus, Lactococcus, Leuconostoc and Pediococcus constituted more than 90% of the total bacterial population, a focus on these genera has been performed (Supplementary Figure 3). Species-level characterization identified a number of interesting features of the HR cheese microbiota coming from the six Pastures areas-Producers (PP).
Among Streptococcus, S. thermophilus represented most of the total, accounting for > 98% of total rel. ab. of this genus. PP2 and PP4 samples only showed about 1% of other (unidentified) Streptococcus and a minimal presence (~0.1%) of S. uberis. The abundance of S. thermophilus could depend on different reasons. A higher presence of S. thermophilus was revealed in cooked cheeses where at the high cooking temperatures a selection of the strains relatively resistant to higher temperatures occurred, as in HR cheese. Moreover S. thermophilus is generally used as starter during cheese-making process. Similar results were found by Cremonesi et al.20 in milk from Rendena dairy cattle, a local Italian breed, where at same farming conditions, a greater quantity of S. thermophilus was found in the local breed milk than in the mainstream breed Italian Holstein, indicating that the animal genetics may have an influence on the composition of the milk microbiota. The abundance of S. thermophilus in HR cheese could be also due to the use of typical wooden equipment during the cheesemaking process, from milk collection to ripening, with wooden poles, bands and shelf, acting as a possible natural inoculate source of S. thermophilus strains of wooden vat origin, as demonstrated in other studies21,22, providing microbial diversity in dairy traditional products. The acidification of curd induced by milk acidification and coagulation during cheese-making process, probably due to the presence of indigenous S. thermophilus, could be of paramount importance for the HR cheese production, where commercial starter cultures are not allowed.
On the other hand, Pediococcus pentosaceus summed up to 98% rel. ab of this genus, apart for PP6 and PP5 samples, was due to P. acidilactici, 2,6% and 9,8%, respectively. Among Lactobacillus, the majority of species were unidentifiable, due to high similarities in their V3-V4 region of 16S rRNA gene, with some notable exceptions: PP5 showed 45.4% of L. delbrueckii and 2.3% of L. fermentum, whereas PP2, PP3 and PP4 were characterized by a somehow higher presence of L. coryniformis (rel. ab. of 6.4%, 12.3% and 20.0%, respectively); PP1 microbiota had a higher rel. ab. of L. buchneri (2.4% compared to an average of 1.1% in the other Alpine areas); finally, PP6 included about 2.2% of L. brevis (compared to an average of 0.8% in other PP). Lactococcus population was dominated by L. lactis with PP4, PP5 and PP2 showing the presence of other unidentifiable species (average rel. ab. 4.8%, 4.7% and 16.9%, respectively). Finally, among Leuconostoc, PP3 was composed for about 87% by L. mesenteroides, whereas PP1, PP5 and PP6 showed a consistent presence of L. paramesenteroides (average rel. ab. of 40.1%, 18.2% and 27.7%, respectively); PP2 and PP4, on the other hand, retained a majority of unidentifiable Leuconostoc (rel. ab. of 72.2% and 86.8%, respectively).
Volatilome
Cheese volatilome represents the final step of catabolism of sugars, proteins and lipids by endogenous and microbial enzymes. Along with the terpenoid content, volatiles characterize the richness of flavor of mountain cheeses, highly appreciated by consumers.
The main volatile compounds detected in the headspace of the HR cheese were: 8 alcohols, 7 free fatty acids, 4 ketones, 3 esters, and 1 aldehyde (Supplementary Files Table S3a-d). Other components present at trace levels were disregarded. Among all classes of compounds, ketones represented the majority of the volatilome, followed by fatty acids and alcohols. Ketones and volatile free fatty acids derive from hydrolysis and subsequent β-oxidation of milk lipids promoted by the intracellular enzymes after bacterial cell lysis. Moreover, ketones have a strong impact on cheese flavor, because of their high odor strength, described as ethereal and sweet with a piquant note. In our study, the major detected ketones were: 3-hydroxy-2-butanone (acetoin), butan-2-one and pentan-2-one. The first two derive from the catabolism of citrate, important flavor compounds in several cheese varieties23. In particular, acetoin is responsible for the buttery notes in cheeses24,25. The predominance of ketones in volatilome was also detected by Panseri et al.26 and Povolo et al.27 in Bitto, a PDO cheese produced in a neighboring area of the HR, in the same grazing period.
Fatty acids were the second most abundant compounds detected in HR volatilome. These carboxylic acids derive from the hydrolysis of triglycerides through the activity of endogenous and microbial enzymes, and are detectable in the headspace of the cheese only when the aliphatic chain is less than ten carbon number. They also contribute strongly to cheese flavor, like butanoic acid, characterized by a rancid cheese-like odor, which plays an important role in Camembert and Grana Padano cheeses flavor28. Fatty acids are precursors of other aroma compounds, such as methyl ketones, alcohols, lactones, aldehydes, and esters29. Branched-chain fatty acids like 3-methylbutanoic and 2-methylpropanoic acid, characterized by a “sweaty” and “ripened cheese” aroma, derive from isoleucine and leucine in the last steps of protein catabolism25 and were also found in PDO Bitto cheese26.
Ethyl butanoate, ethyl hexanoate and ethyl acetate were the most abundant esters isolated in the HR volatilome. They derive by the enzymatic esterification of free fatty acids, generally high in pasture-derived milk and cheese30. They play an important role in the formation of the fruity flavor in cheese and they are able to give this main fruity note to the flavor of some Italian cheeses31, being also the most potent odorants of Grana Padano, Ragusano and Cheddar24,32.
In our study, other important volatile compounds were the branched alcohols 3-Methylbutan-1-ol (medium odor strength, described as alcoholic, fuel), and Phenylethan-1-ol (medium odor strength, described as sweet almond, fresh), derived from the catabolism of amino acids leucine and phenylalanine, respectively25.
Regarding the volatile metabolome as a whole, the comparison among the volatilomes of the six PPs considered, showed some similarities (Figure 3). Some producers obtained cheeses characterized by a higher proteolysis, accompanied by a lower lipolysis like PP2, PP5 and PP6, (see the levels of 2-methyl propanoic acid and hexanoic acid respectively). The higher levels of ethanol from carbohydrates catabolism allowed the formation of ethyl esters from butanoic and hexanoic acid (from lipid catabolism), as for PP1, PP3 and PP4 that showed the higher levels of these volatiles. The odor strength of these esters is very high, and is described as fruity, sweet pineapple, very typical of mountain cheeses during grazing season.
The differences in volatilome seem unrelated to the valley (Supplementary Files Table S3a-d). In fact, cheeses from the same valley, do not show common characteristics; therefore, the single producer has a paramount role in defining the cheese trait, both qualitative and quantitative. Moreover, within the single PP and the single period among producers (Figure 4 and Supplementary Files Table S3a-d; respectively) different volatilomes were observed. Regarding the whole volatilome (Figure 4) PP2 producer showed the maximum overall content. The complexity of HR volatilome, like in other mountain cheeses, derives, mainly, by the usage of raw milk, and by the complexity of the microbiota peculiarity of that particular environment and of the whole cheese-making process, starting from feeding (grazing) till ripening in cellars.
Volatilome and microbiome association
As the microbial activity is of a paramount importance in producing volatiles, several studies were carried out investigating how the metabolism of different genera can contribute to cheese volatilome. At this regard, it must be stressed that the metabolic activity is also strain-dependent33. Lactococcus spp. and Lactobacillus spp. are characterized by low lipolytic activity34. Lactobacillus species produce a variety of metabolites, including alcohols, aldehydes, esters, ketones, and acids, which give rise to the characteristic flavor of fermented milk products35. Streptococcus species also contribute significantly to flavor formation in dairy products36. Lactobacillus and Streptococcus species are often used together as starter cultures in the manufacture of fermented dairy products; their complementary activity endows fermented dairy products with special flavor qualities37.
Generally, in our study, a positive association between some volatile compounds and the main bacterial genera, among PPs, it has been found. In some cases, negative correlations were identified, not due to the absence of a bacterial strain, which leads to an increase of some volatiles compounds, but by the fact that the decrease in the quantity of one bacterial group is associated with the presence of others, positively associated to some volatiles compounds.
The volatiles that showed high and significant (P-value <0.05) Pearson’s correlations with bacterial microbiome composition are reported in Figure 5. Among the total volatiles, most ethyl-esters, along with ethanol, are correlated with Leuconostoc, as observed by Pogačić et al.35, suggesting their potential combination as adjunct flavoring in cheese, associated with the development of fruity, sweet, and floral flavor notes, in particular in long ripened cheeses28. The 3-Hydroxybutan-2-one (or acetoin, derived by citrate metabolism) is strongly and positively associated with Lactococcus genera, confirming the results of Gallegos et al.38 and of Pogačić et al.35 where 3-Hydroxybutan-2-one and 3-methylbutan-1-ol were identified as relevant VOCs respectively for L. rhamnosus, L. paracasei, L. sakei and L. lactis subsp lactis, L. cremoris subsp. cremoris, L. paracasei subsp. paracasei. In our study 2-Heptanol, derived by the reduction of the corresponding ketone, showed a strong positive correlation with Lactobacillus as found by Dan et al.39, where this secondary alcohol was detected in samples fermented by different strains of Lactobacillus. Hexanal, a result of lipid oxidation, and its corresponding alcohol 1-hexanol, showed a strong positive correlation with Streptococcus, confirming the findings of Dan et al.40 where this compound was present in milk fermented by S. thermophilus. The only (fatty) acid that showed a significant positive correlation with Leuconostoc genus was propionic acid, derived by carbohydrate metabolism, as demonstrated by Porcellato et al.41 and Østlie et al.42, where Leuconostoc have been shown to dominate the cheese microbiota in the later stages of ripening with added propionic acid bacteria.
Fatty acid composition
Among the three principal milk macro-components, protein, sugar and fat, the latter is the most susceptible to changes in composition. The most part of fatty acids, in fact, is not produced de novo in the mammary gland, but derives from the bloodstream, after the ruminal activity on feeding matter. When feeding changes, like along the summer grazing period, a change in fatty acid composition takes place. Several studies were carried out on this aspect, in particular with the aim to compare winter versus summer productions or mountain extensive production versus intensive indoor production43-48.
The fatty acid composition (Figure 6) for all the PPs reflects that of a typical mountain cheese produced during grazing period, the most important change being the higher amount of polyunsaturated fatty acids (Figure 6E-F) and branched chain fatty acids (Figure 6C) deriving from pastures rich in dicotyledons. The level of total CLA varied from 1.06% to 1.85% (Supplementary Files Table S4), in agreement with data reported for other PDO cheeses produced during grazing period15,17,48, when no maize silage is given to cattle. Mean value for CLA (% FA) obtained in our study is more than twice (overall Least Mean Square ± Standard Error Means, 1.48 ±0.05, Figure 6G) respect to plain cheeses, where cattle are fed concentrates44,47. For the same reason, similar findings were noticed for BCFA (overall LMS ± SEM, 3.27 ± 0.05, (Figure 6C) and n-3 PUFA (overall LMS ± SEM, 0.80 ± 0.05, Figure 6E-F), confirming data of other authors49,50.
Among producers, PP3 and PP5 showed a higher abundance of SFA and PUFA, and a lower abundance of MUFA (Fig 6A, E-F, D; respectively). Since, during the production season, herds move to different grazing areas at diverse altitudes, characterized by different composition and phenological stage of herbs, a significant difference in the fatty acid composition among the PPs (Fig 6) and within the same PP during the three grazing periods (Supplementary Files Table S4) was detected. Producer PP3 that grazed the herd at higher altitudes produces cheeses with a higher content of CLA and a lower n6/n3 ratio (Figure 6G-H). The PP1, PP5 and PP6 showed huge variability among different GPs (Supplementary Files Table S4), indicating that their pastures differed during the summer season, especially from GPI and GPII, instead PP2 and PP3, showed less variability. PP5 and PP6 belonged to the same valley, the Bitto di Albaredo valleys (VA).
Fatty acid composition and microbiome association
Although feeding has a paramount importance in fatty acid composition of milk and cheese, it must be noticed that a certain association was revealed between fatty acid composition and the microbiota. Moreover, it must be underlined the influence of feeding on rumen microbiota20, that, in turn, influences cheese microbiota.
The FA composition of the HR cheese showed, for some variables, a significant Pearson’s correlations coefficient with the bacterial microbiota constituted by the five genera of Streptococcus, Lactococcus, Leuconostoc, Macrococcus, and Serratia. The abundance of Streptococcus showed a significant positive association (P-value < 0.05) with SFA, SCFA, PUFA n-6, PUFA n-3 (R = 0.30, R = 0.34, R = 0.27, R = 0.30; respectively) and a significant negative correlation with MUFA (R = − 0.45). Lactococcus, Leuconostoc, Macrococcus and Serratia abundance was negatively associated (P-value < 0.05) with SFA (R= − 0.33, R= − 0.46, R= − 0.30; R= − 0.33 respectively) and SCFA (R= − 0.31, R= − 0.28, R= − 0.31, R= − 0.28; respectively).
Leuconostoc and Macrococcus also demonstrated a significant positive association (P-value < 0.05) with MUFA (R= 0.31, R= 0.35; respectively). Finally, the Serratia abundance was positively associated (P-value < 0.05) with BCFA (R= 0.47) and CLA (R= 0.38).
Terpenoids composition
A big difference between mountain and plain vegetation is the content of terpenoids, an important class of compounds derived by isoprene units emitted by mountain herbs, mainly dicotyledons, in defense to predators. Dicotyledons are the main constituents of summer pasture and, among these, Apiaceae, Lamiaceae, and Asteraceae families contain high levels of terpenoids, while Poaceae along with Fabaceae, mainly constituents of the indoor feeding from plain, do not51.
Ten terpenoids resulted statistically different (P < 0.01) among the PPs (Table 2). A similar composition was noticed in other PDO mountain cheeses as ‘‘Fontina Valle d’Aosta” and ‘‘Bitto” 16,17,26. The richness (both qualitatively and quantitatively) of terpenoid moiety found in all the HR cheese samples, is typical of a mountain production during summer grazing period, where these volatile lipophilic compounds are transferred to the milk and finally to the cheese51. α-pinene, camphene and β-pinene were the most abundant compounds in HR cheeses. Their presence is confirmed in many plants from alpine pasture17,51.
In our study, a difference in terpenoid levels was detected among PPs (Table 2), suggesting that these compound can discriminate productions between different pastures15. Panseri et al.26 in a study comparing cheeses from different producers of PDO Bitto, a neighborly produced cheese, didn’t found such difference. The producers PP3, PP5 and PP6, settled at the highest grazing areas (2,000, 1,950, 1,930 m a.s.l., respectively), produce cheeses characterized by the higher levels of total terpene contents. Among PPs, PP5 showed the highest abundance of terpenoids. The terpenoids that best differentiated the PPs were allo-ocimene, α-terpinolene, α-pinene, and δ-3-carene. The latter, especially, was directly related to the altitude of grazing within the same producer PP5 (Supplementary Files Table S5).