Massive Survey on Bacterial–Bacteriophages Biodiversity and Quality of Natural Whey Starter Cultures in Trentingrana Cheese Production

This study focused on the microbial and bacteriophages identification and characterization in cheese-production facilities that use natural whey starter (NWS) cultures for Trentingrana production. Bacterial and phage screening was carried out on cooked not acidified whey and NWS samples isolated from six dairy factories, for 4 consecutive days in four different months. By means of a combined approach, using plate counts, bacterial isolation, and metataxonomic analysis Lactobacillus helveticus was found occurring as the dominant species in NWS cultures and Levilactobacillus brevis as codominant in the cheese factories where the temperature of NWS production was mainly lower than 40°C, suggesting that the variability in the parameters of the NWS culture preparation could differently modulate the bacterial species in NWS cultures. Using turbidity test approach on 303 bacterial isolates from the NWS cultures, 120 distinct phages were identified. L. helveticus phage contamination of NWS cultures was revealed in most of the analyzed samples, but despite the great recovery of bacteriophage contamination cases, the microbial quality of NWS cultures was high. Our results support the presence of natural bacteriophage resistance mechanisms in L. helveticus. The use of NWS cultures probably creates an ideal environment for the proliferation of different L. helveticus strains balanced with their phages without a clear dominance. It is evident, from this study, that the presence of a high biodiversity of NWS bacterial strains is relevant to avoid phages dominance in NWS cultures and consequently to keep a good acidification ability.

Introduction

Traditional and artisanal cheese productions are often based on fermentation processes carried on by defined starter cultures used to achieve defined and typical flavors and/or textures (Gatti et al., 2014; Gobbetti et al., 2018b). Natural whey starters (NWSs) are traditionally used in Italian long ripened hard cheese production as for Parmigiano Reggiano, Grana Padano, and Trentingrana (Bertani et al., 2020). Refreshed daily in the dairy factory from whey collected at the end of the cheese-making process (Bertani et al., 2020), these starters are mainly characterized by thermophilic lactic acid bacteria (LABs) such as Lactobacillus helveticus, Lactobacillus delbrueckii ssp. bulgaricus and L. delbrueckii ssp. lactis (Rossetti et al., 2008).

Together with raw milk quality, bacteriophages (or phages) may represent a dealing factor able to negatively affect fitness and performance of dairy starter cultures (Carminati et al., 2011). Phages are viruses able to infect bacteria; they are present in all the environments, including niches related to human activities as in the case of dairy factories (Carminati et al., 2011). Despite intense efforts as adapted factory design, sanitations, adequate ventilation, and culture rotation, a complete phage eradication in the dairy industry remains a utopian goal, and phage infection of starter LABs is still the most common cause of slow and/or incomplete fermentation (Guglielmotti et al., 2012). Consequently, high pH values and residual lactose may promote growth of pathogenic or spoilage bacteria, negatively affecting the quality or the yield of the final product (Josephsen and Neve, 1998). The most permanent source of phage in dairy is the raw milk, and considering that Trentingrana dairy factories use raw milk in cheese production, they can rapidly grow up to high concentration in NWS cultures (Guglielmotti et al., 2012).

Previous works showed the coexistence in NWS for Grana Padano of phage ecologically related to bacterial strains belonging to species such as L. helveticus and L. delbrueckii (Zago et al., 2006, 2008). As in other ecosystems, bacteriophages play an ecological role within NWS culture, acting as a biological pressure agent following the “kill-the-winner” hypothesis. NWS culture bacteriophages lead the natural selection of phage-resistant bacterial strains, where the fastest-growing bacteria are inhibited in dominating the NWS community, thus preserving the overall technological performances of the NWS cultures (Zago et al., 2008; Gobbetti et al., 2018a). Therefore, the comprehension of bacteria and phage biodiversity dynamic in NWS cultures is of prominent importance, and a good management of NWS processing for the Grana-like cheese production should consider the overall players in this ecological niche.

It has been demonstrated that small changes in the technological parameters, such as curd cooking temperature, titratable acidity, and pH, could affect the bacterial consortium present in NWS cultures (Rossetti et al., 2008). Despite the findings provided by these studies, the microbial and bacteriophage biodiversity of NWS culture is not fully known. To the best of our knowledge, this is the first study reporting the dynamics of bacteria and phage biodiversity in NWS by means of a massive sampling system in cheese factories following different NWS technologies of production. Therefore, in this study, we focused on NWS, collected just before addition to the vat milk, and cooked not acidified whey (cNAW), collected at the end of curdle cooking.

The samples were collected from six selected Trentingrana dairies, over 1 year of production, for 144 samples: 72 cNAW and 72 NWS samples.

This study was performed to understand the dynamics of bacteria and phages of NWS cultures over time, so in this work, we (i) characterized the dynamic of bacterial and phage community of the NWS collected to increase the existing knowledge and (ii) investigated if and how the cheese factory technology of NWS production affects the bacterial and phage community in the NWS cultures. The final goal is to maintain the high quality of Trentingrana production, avoiding economic loss associated with failure in milk fermentation and low cheese yield.

Materials and Methods

Collection of cNAW and NWS Samples, and Determination of Temperatures and pH

cNAW and NWS culture samples used in Trentingrana manufacturing were collected for four consecutive days in four different months (February, May, August, and November 2018) from six dairy factories (A–F) located in the province of Trento (north eastern of Italy) and operating into the Trentingrana Protected Designation of Origin cheese area of production (total of 144 samples).

The Trentingrana production had stopped at cheese factory C for July and August months, so NWS and cNAW samples were collected in September only for cheese factory C.

NWS titratable acidity was measured and expressed as Soxhlet–Henkel degrees (°SH); the pH was determined using a pH electrode (Crison Instruments, Barcelona, Spain). Samples were shipped to the laboratory under liquid nitrogen and stored at −80°C before the analysis.

The temperature of the room and of the NWS tank has been determined using Testo 175 T2 (Testo Ltd., Alton, United Kingdom). Temperature was recorded each minute for the overnight fermentation of the NWS.

Microbiological Counts and Isolation

All the samples were diluted in sterile peptone water and plated onto whey agar medium (WAM) according to Gatti et al. (2003) for cultivating thermophilic lactobacilli in anaerobic conditions for 72 h at 45°C. All culture media and anaerobic system were purchased from Oxoid (ThermoFisher Scientific, Milan, Italy).

Eight colonies grown onto WAM plates were randomly isolated from NWS samples (plates with a number of colonies in the range of 10–300). Each isolate was purified by subsequent culturing in whey broth (the broth version of the WAM used for plate counting). Pure cultures were stored at −80°C in glycerol (40% vol/vol) stocks. Cell morphology was determined by microscopic observation, Gram characterization was performed applying the KOH method (Gregersen, 1978), and catalase activity was tested by adding 5% H₂O₂ drops on the colonies.

DNA Extraction and Genotypic Identification of the NWS Bacteria

All bacterial isolates were grown 48 h in the whey broth culture at 45°C before DNA extraction. The bacterial DNA was isolated using Quick-gDNA^TM MicroPrep (Zymo Research, Italy) following the manufacturer’s instructions.

Randomly amplified polymorphic DNA–polymerase chain reaction (RAPD-PCR) was carried out in a total volume of 25 μL using primer M13 (Rossetti and Giraffa, 2005). Cluster analysis of DNA patterns was carried out using GelCompar II-BioNumerics software (package version 6.0; Applied Maths, bioMérieux, Belgium), exploiting the unweighted pair group method arithmetic averages. Similarity of PCR fingerprint profiles was calculated based on Pearson product–moment correlation coefficient. The threshold breakpoint value was fixed to 80%; isolates with similarity coefficient higher than 80% were classified into the same cluster, according to Gatti et al. (2003).

One isolate representative of each biotype was genotypically identified by 16S rRNA gene analysis. A fragment of the 16S rRNA gene was amplified using the primers 27F (5′-GAGAGTTTGATCCTGGCTCAG-3′) and 1495R (5′-CTACGGCTACCTTGTTACGA-3′), designed by Grifoni et al. (1995). The PCR products were purified using the Exo-SAP-IT^TM kit (USB Co., Cleveland, OH, United States) and sequenced in an ABI PRISM 3100 sequencer (Applied Biosystems, Foster City, CA, United States), using the BigDye Terminator v1.1 cycle sequencing kit (Applied Biosystems). The obtained sequences were compared using the BLAST algorithm available on the National Center for Biotechnology Information (NCBI, United States). All amplifications were run in a T100^TM ThermalCycler (Bio-Rad Laboratories, Hercules, CA, United States).

DNA Extraction, MiSeq Library Preparation, and Illumina Sequencing

Illumina analysis was performed on 24 cNAW and 24 NWS samples (one cNAW and NWS sample for each month and cheese factory). Three milliliters of sample was centrifuged at 3,200 × g for 15 min at 4°C. The genomic DNA was extracted from the pellet using the Power Food^TM Microbial DNA Isolation Kit (Qiagen, Hilden, Germany) according to the manufacturer’s instructions. All DNA samples were purified by PowerClean DNA Clean-up Kit (Qiagen). The DNA quality and concentrations were determined by NanoDrop^TM 8000 Microvolume UV-Vis spectrophotometer (ThermoFisher Scientific).

Amplicon library preparation, quality and quantification of pooled libraries, and pair-end sequencing using the Illumina MiSeq system (Illumina, United States) were carried out at the Sequencing Platform in Fondazione Edmund Mach (FEM, San Michele all’Adige, Italy). Briefly, for each sample, a 464-nucleotide sequence of the V3–V4 region (Baker et al., 2003; Claesson et al., 2010), of the 16S rRNA gene (Escherichia coli positions 341 to 805) was amplified. Unique barcodes were attached before the forward primers to facilitate the pooling and subsequent differentiation of samples. To prevent preferential sequencing of the smaller amplicons, the amplicons were cleaned using the Agencourt AMPure kit (Beckman Coulter, Brea, CA, United States) according to the manufacturer’s instructions; subsequently, DNA concentrations of the amplicons were determined using the Quant-iT PicoGreen dsDNA kit (Invitrogen, ThermoFisher Scientific) following the manufacturer’s instructions. In order to ensure the absence of primer dimers and to assay the purity, the generated amplicon libraries quality was evaluated by a Bioanalyzer 2100 (Agilent, Palo Alto, CA, United States) using the High Sensitivity DNA Kit (Agilent). Following the quantitation, cleaned amplicons were mixed and combined in equimolar ratios.

Illumina Data Analysis and Sequences Identification by QIIME2

Raw paired-end FASTQ files were demultiplexed using idemp¹ and imported into Quantitative Insights Into Microbial Ecology (Qiime2, version 2020.8). Sequences were quality filtered, trimmed, denoised, and merged using DADA2 (Callahan et al., 2016). Chimeric sequences had been identified and removed via the consensus method in DADA2. Representative bacterial sequences had been aligned with MAFFT and used for phylogenetic reconstruction in FastTree using plugins alignment and phylogeny (K. Price et al., 2009; Katoh and Standley, 2013). α- and β-diversity metrics had been calculated using the core-diversity plugin within QIIME2 and visualized by emperor (Vazquez-Baeza et al., 2013). Bacterial taxonomic and compositional analyses were carried on by using plugins feature-classifier². A pretrained naive Bayes classifier based on the Greengenes 13_8 99% operational taxonomic unit (OTU) database, which had been previously trimmed to the V4 region of 16S rDNA, bound by the 341F/805R primer pair, was applied to paired-end sequence reads to generate taxonomy tables. The data generated by MiSeq Illumina sequencing were deposited in the NCBI Sequence Read Archive and are available under Ac. Number PRJNA695135³ from sample SAMN17602423 to SAMN17602468.

Enrichment, Isolation, and Purification of Phages

All the 72 samples of NWS cultures were screened for the presence of bacteriophages. The related cNAW samples were used to isolate the virulent phages, which used the NWS isolated bacteria as host. Phage enrichment was performed as previously reported by Zago et al. (2005). Briefly 1 g of CaCO₃ was included in tubes of 40-mL cNAW sample, acidified with 0.5 mL of sodium maleate buffer (maleic acid 0.05 M titrated with 0.2 N NaOH until pH 5.15 is reached), and then inoculated with 1 mL of NWS from the same sampling day and dairy. After anaerobic incubation at 45°C for 72 h, the tubes were centrifuged at 4,000 revolutions/min for 15 min at 4°C. Supernatants obtained from the enrichment step were filtered with 45-μm membrane filter pore size (Merck–Millipore, Darmstadt, Germany) and checked for the presence of phages by means of turbidity test.

Tubes of WAM-Ca broth (WAM supplemented with 10 mM CaCl₂) were added with 5% of the previously filtered supernatants and inoculated with 2% mL of host culture in the logarithmic growth phase. The tubes were incubated at 45°C and observed each 6 h until cell lysis occurred. The lysis was detected by the higher clearness of the broth, compared to controls (tubes without addition of the supernatants), which were more turbid due to bacteria growth. After 24 h, where no lysis was detected, the bacterial cultures were used for 2% inoculation of new tubes of WAM-Ca broth. After four rounds of inoculations without lysis, the bacterial isolate was considered as not suitable for the phage isolation. Where lysis was observed, the tubes were then filtered through 0.45-μm membrane filter pore size (Merck–Millipore), set to pH 7, and stored at 4°C until use, for no more than 7 days. Phage stocks were stored at −80°C with 40% vol/vol of glycerol.

Titer Determination of Phage

Phage titer was expressed as the number of plaque-forming units per mL of sample (pfu/mL). The suspension of the phage was diluted sevenfold by using phosphate-buffered saline added with CaCl₂ 10 mM. One hundred thirty microliters of the host bacteria in the logarithmic growth phase was added to 5 mL of semisolid WAM-Ca (WAM-Ca, 4% wt/vol agar). The mixture was poured on the surface of WAM-Ca plates and then inoculated with 10 μL of phage dilutions in seven plate sectors. After anaerobic incubation at 45°C, clear plaques were checked. The phage titer was calculated using the formula titer (pfu/mL) = plaque number × dilution factor × 100.

Detection of the Lysogenic State

Each isolated bacterial biotype, characterized as a phage host, was tested for the detection of prophage, adapting the method of Cochran and Paul (1998). Overnight WAM cell cultures were split into two aliquots: the first aliquot was added with mitomycin C (MytC⁺, 1 μg/mL); the second aliquot was not added with MytC⁺ and used as control. After 24 h of incubation at room temperature, in the dark, induced phages were detected as follows: 100 μL of each cell suspension was used to inoculate 7 mL of soft WAM (WAM, 0.7% wt/vol agar) plated on a layer of 1.8% (wt/vol) agar and incubated for 24 h at the optimal conditions for each strain. The growth reduction onto MytC⁺ plates, after comparison with control plates, indicated the presence of phage.

If not specified, all chemical compounds were purchased from Sigma–Aldrich (St. Louis, MO, United States).

Statistical Analysis

A normality test (Shapiro-Wilk W) was performed, as well as a nonparametric test (Kruskal--Wallis) analyzing the day of collection as the independent variable and the microbial plate counts as the dependent variable. All the tests on plate counts were performed using the STATISTICA data analysis software system, version 9.1 (StatSoft, Inc, 2010⁴).

Differences in diversity indices (OTUs number and Shannon diversity index) of different samples were tested by Kruskal–Wallis test by a plug-in implemented in QIIME2. The overall structural changes of bacterial community were visualized by principal coordinate analysis (PCoA) based on both weighted and unweighted UniFrac distance matrices. The statistical significance was assessed via the nonparametric PERMANOVA (permutational multivariate analysis of variance) by means of plug-in implemented in QIIME2. For the differential abundance test, the taxonomy information for each OTU sequence was provided using ANCOM method (Mandal et al., 2015) by means of plug-in implemented in QIIME2.

Results

Thermophilic Bacteria Counts in Trentingrana NWS and cNAW Samples and Technological Process of NWS Production

The WAM counts of thermophilic bacteria in NWS and cNAW samples during the monitoring from February to December are shown in Table 1. NWS thermophiles varied without significant differences in the same and among the different dairy factories, ranging between 7.3 and 8.8 log colony-forming units (cfu)/mL. Only in November, at dairy factory D, NWS cultures’ averaged value was lower than 7.3 (6.7 log cfu/mL). The dairy factory C was the only one showing significant differences in the NWS thermophilic counts that were lower in February (7.5 log cfu/mL) and higher in August (8.6 log cfu/mL). More variability was observed in the cNAW thermophilic counts that ranged between 4.8 and 7.1 log cfu/mL. With the exceptions of dairy factories B and F, cNAW thermophilic counts were showing significant differences during the year without a well-defined trend.

TABLE 1

Table 1. Thermophilic lactic acid bacteria (TLAB) counts onto WAM of Trentingrana NWS and cNAW.

The NWS pH mean values ranged between 3.3 and 3.6 and the titratable acidity between 28°SH and 32°SH / 50 mL without any significant difference for month or cheese factory (data not shown).

As the final acidification activity in NWS cultures is directly associated with a slow temperature decrease from 55°C to 40°C during the overnight process of NWS production (Di Cagno et al., 2006), the NWS production cycle was split into two stages: a first thermophilic stage, when the cNAW is cooled from cooking temperature (55°C) to 40°C, and a second mesophilic stage, after the first cNAW cooling, when the temperature was in the range between 20°C and 40°C (Table 2).

TABLE 2

Table 2. Length in minutes of the first thermophilic and second mesophilic stage of the overnight NWS culture production cycle.

At dairy plants A, B, and C, the first thermophilic was at least 2 h 25 min shorter than the second mesophilic stage, conversely at dairy plant D, the first was approximately 7 h longer than the second stage. At cheese factories E and F, there was no relevant difference between the first and second stage of NWS culture production cycle.

Identification of Bacteria From NWS Cultures

A total of 576 isolates were picked up from WAM agar plates. A RAPD-PCR fingerprinting had been performed (i) to preliminarily recognize the bacterial isolates purified from colonies onto the same plate originated from genetically alike single mother cells and (ii) to investigate the microbial diversity beyond the species level, establishing the number of different biotypes inside cheese factory and over time. One isolate for each biotype was selected for further identification by 16S rRNA gene sequencing.

Of the 576 isolates, 45 did not grow after plate isolation and were discarded, 383 had been identified as L. helveticus, 120 as Levilactobacillus brevis, 13 as Lacticaseibacillus paracasei, eight as Lactiplantibacillus plantarum, and the other seven as Limosilactobacillus fermentum (Figure 1).

FIGURE 1

Figure 1. Distribution of L. helveticus, L. brevis, L. paracasei, L. plantarum, and L. fermentum isolates in NWS samples collected from six dairy factories labeled from A to F four times in 2018 year: February, May, August, and November. At cheese factory C, the samples were collected at the end of September because in August the Trentingrana production stopped.

Cluster analysis of RAPD-PCR profiles showed that, at a similarity level of 85%, L. helveticus and L. brevis, accounting for the 95% of all the isolates, grouped, respectively, into 243 and 51 different biotypes. L. paracasei, L. plantarum, and L. fermentum grouped into two, three, and five different biotypes, respectively (Table 3).

TABLE 3

Table 3. Distribution of L. helveticus, L. brevis, L. paracasei, L. plantarum, and L. fermentum biotypes (n° BioT) in NWS samples collected from six dairy factories labeled from A to F four times a year from February to November 2018.

Lactobacillus helveticus was always dominant in NWS samples of all cheese factories, with the exception of cheese factory A. In particular, its relative presence ranged between 58 and 96% and was higher than 90% in NWS samples at cheese factories E and F. Conversely, L. brevis presence ranged between 0 and 42% with the exception of factory A where its relative presence was 51% and dominated NWS samples from May to November (Figure 1). L. paracasei and L. plantarum isolates were found only in NWS samples collected in February, and the higher presence was at cheese factories D and C, respectively. L. fermentum isolates were found only in NWS samples collected in May and more frequently at cheese factory D.

The distribution of the isolates into the NWS cultures was mainly related with the dairy plant rather than with the month of production. In fact, considering only the month of NWS production, there was no great difference in the distribution of L. helveticus and L. brevis isolates whose relative presence ranged between 62 and 89% and 11 and 34%, respectively.

Characteristics of the Sequencing Data, Diversity Analysis, and Differential Abundance Analysis of the Bacterial Community in cNAW and NWS Samples

With the exception of two NWS samples, the extracted DNA was always successfully amplified in the bacterial V3–V4 16S rRNA gene region. A total of 2,082,145 paired-end sequences were obtained.

To address the hypothesis that species richness and biodiversity vary with sample source (cNAW or NWS sample), month (February, May, August, or November), and cheese factory (A, B, C, D, E, and F), the intragroup diversity estimation (α diversity) was calculated, using both the number of observed OTUs and Shannon diversity indexes. Both α-diversity indexes were not significantly different between cNAW and NWS samples (Kruskal–Wallis, P > 0.1), and no significant difference was found for month or dairy of sampling (data not shown).

In order to assess the amount of variation in bacteria composition among the samples, we calculated the phylogenetic β diversity based on both weighted and unweighted UniFrac distances. The PCoA plot onto weighted UniFrac distance matrix (Figure 2A) shows a little separation of the bacterial populations between cNAW (rings) and NWS samples (squares), visible on axis 2 explaining 18.3% of total variation. Conversely, the PCoA plot onto unweighted UniFrac distance matrix (Figure 2B) showed, with few exceptions, all cNAW and NWS samples overlapped in the same cloud on the left side of the graph. As the unweighted UniFrac distance accounts for the presence/absence of OTUs, so does the weighted UniFrac distance for abundance, although this result suggests that the individual microbial abundance more than the composition in species drives the distance among cNAW and NSW samples.

FIGURE 2

Figure 2. Principal coordinate analysis (PCoA) based on weighted (A) and unweighted (B) UniFrac distances of bacterial OTUs. (A) The first, second, and third axes explained 81.66% of the total variance. (B) The first, second, and third axes explained 69.97% of the total variance. Different symbols represent the different sampling sources: rings for cNAW and squares for NWS samples. Different colors represent the different sampling months: February = red, May = blue, August = orange and November = green.

A PERMANOVA was performed to explore the effects and significance of the variables (Table 4): sample source (cNAW and NWS), cheese factory (A, B, C, D, E, and F), and month of sampling (February, May, August, and November). The test revealed that sample source affected microbial communities (P < 0.05); the month of sampling was slightly significant in cNAW and NSW samples, and in particular, the samples collected in November showed a significantly different bacterial community. Cheese factory was the main factor affecting microbial population diversity in all the samples, and the cheese factory D showed a significantly higher effect on cNAW and NSW samples (Table 4).

TABLE 4

Table 4. PERMANOVA analysis (999 permutations) results for bacterial communities based on weighted UniFrac distances, respectively.

Considering only the OTUs whose relative abundance in each sample was greater than 0.01% (Table 5), Firmicutes was the dominant phylum, able to describe always more than 99% of the bacterial microbiota in all the samples. L. helveticus was always the dominant species over Lactobacillus spp. and Streptococcus spp., constituting more than 58% of the bacterial population in all the samples (Table 5). OTUs identified as Chryseobacterium, Staphylococcus equorum, Lactococcus, Acinetobacter guillouiae, and Acinetobacter johnsonii phylotypes were present in traces (<0.01%) and not in all the samples.

TABLE 5

Table 5. Relative abundances of L. helveticus, Lactobacillus, and Streptococcus spp. of bacterial sequences from cNAW and NWS samples using Illumina MiSeq.

To identify taxonomic groups driving differences among the bacterial community, a differential abundance test by ANCOM method had been performed. No difference was found considering the month of sampling. Some taxa had been found differentially abundant among the cheese factories; in particular, both Lactobacillus and Streptococcus spp. were informative for characterizing the microbial community of cheese factories A and B.

Dissemination and Titer of Phages in NWS Samples and Lysogenic State Analysis

All the 72 NWS samples collected from the six local cheese factories had been phage screened. The bacterial biotypes isolated from the different NWS cultures (Table 3) represented a reliable set of isolates suitable as hosts for phage detection in NWS samples; consequently, for each NWS sample screening, a set of putative bacterial hosts was selected, grouping all the isolates representative for each different biotype found in each NWS sample.

One hundred twenty phages were readily detected in 41 of the 72 NWS samples (Table 6) and were all lytic for their L. helveticus host isolate, with the exception of two phages lytic for L. fermentum. To detect the presence of prophage in their genomes, the 120 host bacterial isolates were treated with sublethal concentrations of MytC⁺. In none of the isolates, the MytC⁺ induced a sharp decrease in plate counts, confirming the lytic nature of all the phages.

TABLE 6

Table 6. Distribution of phage isolates (n° Is.) from NWS samples collected from six dairy factories labeled from A to F, four times a year from February to November 2018.

Phage titers were determined in all the 41 positive NWS samples. For L. helveticus phages, titers ranged between 2 × 10⁶ and 9 × 10⁸ pfu/mL with the exception of all the 12 phages from dairy factory C and six phages from dairy factory F, whose titers were in a range between 2 × 10² and 4 × 10⁵ pfu/mL.

Discussion

The NWS cultures used in Grana-like cheese production originate from artisanal back-slopping practices and face several selective pressures during the production process such as heat treatments and strong lactic acid fermentation. The bacterial consortia are characterized by a lower species complexity but high degree of strains diversity (Erkus et al., 2013), and bacteriophages have a regulatory role in the bacterial population dynamics through density-dependent predation (Rodriguez-Valera et al., 2009).

This study was carried out to investigate the bacterial and phage dynamics of cNAW and NWS production analyzing samples collected in 72 different days of Trentingrana production in six dairy plants, over a year period, by means of a combined approach, using plate counts, bacterial isolation, metataxonomic analysis, and phage isolation.

The NWS culture thermophilic counts, acidity, and pH were not significantly different for month (with exception of cheese factory C) or cheese factory, and their values were according to previously studied NWS cultures for Trentingrana (Franciosi et al., 2012) and Parmigiano Reggiano (Coppola et al., 2000). The homogeneity of the microbial counts throughout the overall sampling year could be the consequence of the consolidated process of NWS production from the cNAW, which guarantees for each dairy plant a microbial homogeneity of the NWS cultures. At cheese factory C, the thermophilic counts of NWS cultures were significantly higher in September and November. This cheese factory was the only one to stop the Trentingrana production in July and August. In September, the production of NWS cultures started again with the cheese. Therefore, this difference in thermophilic counts at dairy plant C could be a consequence of the total change in the NWS cultures that in September were no more linked to the NWS cultures of the previous months.

The sequencing of the 16S ribosomal gene of the isolates dominating the thermophilic microbial population in NWS cultures revealed the presence of five species, of which L. helveticus and L. brevis were the two dominant ones (95% of the thermophilic analyzed isolates). L. paracasei, L. plantarum, and L. fermentum were present in few NWS samples collected only in February and May. L. helveticus is well known as dominant in NWS cultures for Grana-like cheeses (Rossetti et al., 2008; Morandi et al., 2019; Bertani et al., 2020). L. fermentum had already been isolated in low abundance in NWS cultures of Grana-like cheeses (Rossetti et al., 2008; Morandi et al., 2019; Bertani et al., 2020) and the presence of two biotypes of L. paracasei and three of L. plantarum could be explained as an occasionally NWS culture contamination that occurred only in the February month of sampling. By contrast, the dominance of L. brevis together with L. helveticus in the thermophilic community of the NWS cultures in three cheese factories (L. brevis constituted > 30% of the isolates in A, B, and C) and during the overall year was surprising. The presence of 51 biotypes of L. brevis could not be explained as a contamination. L. brevis is one of the most common nonstarter LAB species found in the dairy industry (Settanni and Moschetti, 2010) and was already found in traces among the mesophilic bacteria of NWS cultures used for Parmigiano Reggiano (Coppola et al., 2000). To our knowledge, this is the first time this species had been isolated as one of the codominant thermophilic bacteria in NWS for Grana-like cheeses. L. brevis was in higher abundance in NWS from cheese factories A, B, and C where the overnight NWS production was dominated by a mesophilic stage with temperature of less than 40°C that is more selective for L. brevis than L. helveticus. As L. brevis growth could be strongly inhibited by the presence of homofermentative bacteria such as L. helveticus (Laleye et al., 1989), we speculate that a longer thermophilic stage greater than 40°C could be desirable because it is more selective for the growth of L. helveticus over L. brevis during the NWS culture production. L. brevis is an obligate heterofermentative bacteria producing CO₂ from lactose and known to be a possible cause of early gas production in cheeses (Michael and Mullan, 2000); therefore, it is better to limit its abundance in NWS cultures.

Cluster analysis of RAPD-PCR bacterial isolates profiles showed a higher genotypic diversity within L. helveticus than L. brevis isolates; in fact, the 383 isolates of L. helveticus grouped into 243 biotypes, whereas the 120 L. brevis into 51 biotypes. The better adaptability of the L. helveticus species to the stressing conditions caused by the strongly acidic environment of the NWS culture may explain its higher strain diversity and consequently the lower L. brevis diversity at strain level. Many L. helveticus strains were probably able to adapt to acidity and high temperature of NWS cultures, and conversely, fewer L. brevis biotypes were resistant to these same stress conditions. The qualitative distribution of the biotypes into the different cultures was mainly dairy plant specific rather than correlated with the month of production. This is confirming our previous speculation—that the dominant LAB species of the NWS cultures are modulated by the variability in the parameters of the culture preparation that are different in each dairy plant. NWS cultures characterized by a low species-level complexity but a higher strain-level diversity were already described in other studies performed on dairy starter cultures strains (Frantzen et al., 2018; Schmid et al., 2018; Somerville et al., 2019). The differences in NWS culture production associated with the dairy plant characteristics (farmers, dairy factory operators, and environment) may affect the selection of NWS culture biotypes. Previous studies have also highlighted that NWS species and biotypes show different ability to adapt to dairy ecosystem (Giraffa and Neviani, 1999; Coppola et al., 2000; Moser et al., 2018).

Regarding the relative composition of the bacterial community, the α indices did not show significant differences in the biodiversity of the microbial communities of cNAW and NWS samples. The β-diversity analysis showed that the same bacterial species were characterizing both NWS and cNAW samples but with different relative abundances, and the dairy factory had a higher significant effect than the month of sampling on the abundance of the different species inside the microbial community. Dairy factory was already considered the most important source of variation of the Trentingrana cheese quality index (Bittante et al., 2011), and our results are in agreement with this previous study. The pairwise comparison showed that the species and relative abundances of cNAW and NWS cultures were very similar in all the dairy plants with the exception of the cheese factory D, whose NWS samples showed the highest relative abundance of L. helveticus and the lowest of Lactobacillus spp. and Streptococcus thermophilus. The high abundance of L. helveticus could be due to the NWS cycle of production: dairy plant D was the only cheese factory where the NWS cultures were produced with a thermophilic stage significantly longer than the mesophilic one (approximately 7 h longer), confirming the importance of NWS technology of production in modulating the microbial community of the NWS cultures. The analysis of both metataxonomic and quantitative data obtained by plate counts and isolation suggested that L. helveticus decreased in relative abundance but increased in amount from cNAW to NWS. This was in agreement with Bertani et al. (2020), showing that L. helveticus growth was mainly favored, in the acidic conditions occurring during the overnight fermentation of NWS. In Trentingrana dairy plants A, B, and C, the NWS cultures showed a high abundance of L. brevis isolates adapted to survive in thermophilic environment. Their dominance was also confirmed by metataxonomic data. Lactobacillus spp. relative abundance was very high (29–39%) when compared to previous metataxonomic studies on NWS for Grana-like cheeses (Bertani et al., 2020) where Lactobacillus spp. was never higher than 1%. L. delbrueckii was neither isolated nor sequenced in our study, and the abundance of S. thermophilus was never higher than 2.6% in NWS cultures. Both results were in agreement with a previous study by Morandi et al. (2019) on NWS cultures for Trentingrana production collected from 2014 to 2016, when a drastic reduction of L. delbrueckii and S. thermophilus isolates had been observed.

Ubiquitously distributed among all ecosystems, bacteriophages are commonly present in dairy plants (Pujato et al., 2019). It has been widely shown that all the most common starter bacteria species used in dairy industry, such as S. thermophilus, Lactococcus, and Lactobacillus spp. (Zago et al., 2006, 2015; Oliveira et al., 2018), are suitable hosts for phages, with consequent possible impairment in fermentations performance.

Almost all the 120 investigated phages were infecting the L. helveticus isolates, and 54 were able to form lysis plaque. This low ability in plaques forming was already observed in L. helveticus phages (Carminati et al., 2011; Zago et al., 2015). All the 120 phages were isolated over 1 year from NWS samples suggesting a simultaneous presence of phages and different sensitive bacteria within the cultures, as shown in previous studies (Zago et al., 2005, 2015). Although the widespread presence of lytic phages, all the NWS cultures were successful in their acidification activity, suggesting that the phages virulence could be counteracted by the presence of many bacteria, with different phage sensitivity, representing a natural way to control phage overwhelm. The individual resistance of the NWS isolates to phage predation could be one of the factors influencing NWS culture diversity as previously confirmed by Spus et al. (2015). In this study, the phage presence was not equally distributed among dairy plants; in fact, approximately 58% of the phages have been isolated from NWS cultures produced in only two dairy plants (E and F). It is remarkable that the same dairy plants E and F showed the lower biodiversity in the number of L. helveticus biotypes, confirming that the lower the bacterial biodiversity, the higher could be the phage presence. Few bacterial biotypes have been already observed in dairy starter cultures, where the phage predation operated mainly at strain rather than at species level (Erkus et al., 2013). The high biodiversity of dairy bacteria biotypes may prevent massive phage multiplication in NWS cultures and explain the recovery of a lower number of phages from the dairy plants with the higher biodiversity in L. helveticus biotypes. We speculate that, if phage activity changed in NWS culture, a population of few bacteria biotypes could be more susceptible to a phage predation with the consequent loss of NWS activity and generation of a defective cheese (Madera et al., 2004). Conversely, NWS cultures rich in different biotypes have a better chance to include phage-resistant bacteria able to counteract changes in phages community.

These considerations highlight the importance in assessing the strain biodiversity in NWS cultures: the higher the biodiversity, the higher the probability of a natural selection of bacterial strains resistant to phage predation and thus higher preservation of the overall activity of the NWS cultures, as also stated by Carminati et al. (2011).

In conclusion, this study added further knowledge on the microbial composition and diversity of NWS cultures for Grana-like cheeses in relation with the phage isolation. We speculate the importance of the technological process of NWS production and in particular of the thermophilic and mesophilic stages as drivers of the bacterial species dominant in NWS cultures. In addition, we found a possible correlation among bacterial starter biodiversity and the number of recovered lytic phages.

The abundance of phages in NWS cultures underlines again the importance of phage control strategies in the dairy industry. Ongoing studies are in progress to evaluate the roles of bacterial biodiversity, phages, and NWS production technology on the species composition of Trentingrana NWS cultures.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/, PRJNA695135.

Author Contributions

EF and AMe devised the study. EF and AMa drafted the manuscript. EF, AMa, NC, AG, IC, and MR performed sampling plane and the experiments. MZ contributed to the study with her know-how on phages. EF, AMe, and KT provided resources and intellectual input that supported the study. All authors contributed to the article and approved the submitted version.

Funding

This work was funded by the Trento province Region Rural Development Program (PSR).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The authors would like to thank all the operators of the dairy facilities for providing precious information about the NWS production.

Footnotes

1. ^ https://github.com/yhwu/idemp/blob/master/idemp.cpp
2. ^ https://github.com/qiime2/q2-feature-classifier
3. ^ https://www.ncbi.nlm.nih.gov/bioproject/695135
4. ^ www.statsoft.com

References

Baker, G. C., Smith, J. J., and Cowan, D. A. (2003). Review and re-analysis of domain-specific 16S primers. J. Microbiol. Methods 55, 541–555. doi: 10.1016/j.mimet.2003.08.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Bertani, G., Levante, A., Lazzi, C., Bottari, B., Gatti, M., and Neviani, E. (2020). Dynamics of a natural bacterial community under technological and environmental pressures: the case of natural whey starter for Parmigiano Reggiano cheese. Food Res. Int. 129:108860. doi: 10.1016/j.foodres.2019.108860

PubMed Abstract | CrossRef Full Text | Google Scholar

Bittante, G., Cechinato, A., Cologna, N., Penasa, M., Tiezzi, F., and De Marchi, M. (2011). Factors affecting the incidence of first-quality wheels of Trentingrana cheese. J. Dairy Sci. 94, 3700–3707. doi: 10.3168/jds.2010-3746

PubMed Abstract | CrossRef Full Text | Google Scholar

Callahan, B. J., McMurdie, P. J., Rosen, M. J., Han, A. W., Johnson, A. J., and Holmes, S. P. (2016). DADA2: high-resolution sample inference from Illumina amplicon data. Nat. Methods 13, 581–583. doi: 10.1038/nmeth.3869

PubMed Abstract | CrossRef Full Text | Google Scholar

Carminati, D., Zago, M., and Giraffa, G. (2011). “Ecological aspects of phage contamination in natural whey and milk starters,” in Bacteriophages in Dairy Processing, eds A. L. Quiberoni and J. A. Reinheimer (New York, NY: Nova Science Publishers, Inc), 79–97.

Google Scholar

Claesson, M. J., Wang, Q., O’Sullivan, O., Greene-Diniz, R., Cole, J. R., Ross, R. P., et al. (2010). Comparison of two next-generation sequencing technologies for resolving highly complex microbiota composition using tandem variable 16S rRNA gene regions. Nucleic Acids Res. 38:e200. doi: 10.1093/nar/gkq873

PubMed Abstract | CrossRef Full Text | Google Scholar

Cochran, P. K., and Paul, J. H. (1998). Seasonal abundance of lysogenic bacteria in a subtropical estuary. Appl. Environ. Microbiol. 64, 2308–2312. doi: 10.1128/aem.64.6.2308-2312.1998

PubMed Abstract | CrossRef Full Text | Google Scholar

Coppola, R., Nanni, M., Iorizzo, M., Sorrentino, A., Sorrentino, E., Chiavari, C., et al. (2000). Microbiological characteristics of Parmigiano Reggiano cheese during the cheesemaking and the first months of the ripening. Lait 80, 479–490. doi: 10.1051/lait:2000139

CrossRef Full Text | Google Scholar

Di Cagno, R., De Angelis, M., Limitone, A., Fox, P. F., and Gobbetti, M. (2006). Response of Lactobacillus helveticus PR4 to heat stress during propagation in cheese whey with a gradient of decreasing temperatures. Appl. Environ. Microbiol. 72, 4503–4514. doi: 10.1128/aem.01829-05

PubMed Abstract | CrossRef Full Text | Google Scholar

Erkus, O., de Jager, V. C., Spus, M., van Alen-Boerrigter, I. J., van Rijswijck, I. M., Hazelwood, L., et al. (2013). Multifactorial diversity sustains microbial community stability. ISME J 7, 2126–2136. doi: 10.1038/ismej.2013.108

PubMed Abstract | CrossRef Full Text | Google Scholar

Franciosi, E., Gardini, F., Monfredini, L., Tabanelli, G., Fabris, A., Endrizzi, I., et al. (2012). Does milk treatment before cheesemaking affect microbial and chemical traits of ripened cheese? Grana Trentino as a case study. J. Dairy Sci. 95, 5485–5494. doi: 10.3168/jds.2011-4693

PubMed Abstract | CrossRef Full Text | Google Scholar

Frantzen, C. A., Kleppen, H. P., and Holo, H. (2018). Lactococcus lactis diversity in undefined mixed dairy starter cultures as revealed by comparative genome analyses and targeted amplicon sequencing of epsD. Appl. Environ. Microbiol. 84:e02199-17. doi: 10.1128/AEM.02199-17

PubMed Abstract | CrossRef Full Text | Google Scholar

Gatti, M., Bottari, B., Lazzi, C., Neviani, E., and Mucchetti, G. (2014). Invited review: microbial evolution in raw-milk, long-ripened cheeses produced using undefined natural whey starters. J. Dairy Sci. 97, 573–591. doi: 10.3168/jds.2013-7187

PubMed Abstract | CrossRef Full Text | Google Scholar

Gatti, M., Lazzi, C., Rossetti, L., Mucchetti, G., and Neviani, E. (2003). Biodiversity in Lactobacillus helveticus strains present in natural whey starter use for Parmigiano Reggiano cheese. J. Appl. Microbiol. 95, 463–470. doi: 10.1046/j.1365-2672.2003.01997.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Giraffa, G., and Neviani, E. (1999). Different Lactobacillus helveticus strain populations dominate during Grana Padano cheesemaking. Food Microbiol. 16, 205–210. doi: 10.1006/fmic.1998.0230

CrossRef Full Text | Google Scholar

Gobbetti, M., Di Cagno, R., Calasso, M., Neviani, E., Fox, P., and De Angelis, M. (2018a). Drivers that establish and assembly the lactic acid bacteria biota in cheeses. Trends Food Sci. Tech. 78, 244–254. doi: 10.1016/j.tifs.2018.06.010

CrossRef Full Text | Google Scholar

Gobbetti, M., Neviani, E., and Fox, P. (2018b). The Cheeses of Italy: Science and Technology. Cham: Springer International Publishing. doi: 10.1007/978-3-319-89854-4

CrossRef Full Text | Google Scholar

Gregersen, T. (1978). Rapid method for distinction of Gram-positive bacteria. Appl. Microbiol. Biothechnol. 5, 123–127. doi: 10.1007/bf00498806

CrossRef Full Text | Google Scholar

Grifoni, A., Bazzicalupo, M., Di Serio, C., Fancelli, S., and Fani, R. (1995). Identification of Azospirillum strains by restriction fragment length polymorphism of the 16S rDNA and the histidine operon. FEMS Microbiol. Lett. 127, 85–91. doi: 10.1111/j.1574-6968.1995.tb07454.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Guglielmotti, D. M., Mercanti, D. J., Reinheimer, J. A., and Quiberoni, A. (2012). Review: efficiency of physical and chemical treatments on the inactivation of dairy bacteriophages. Front. Microbiol. 2:282. doi: 10.3389/fmicb.2011.00282

PubMed Abstract | CrossRef Full Text | Google Scholar

Josephsen, J., and Neve, H. (1998). “Bacteriophages and lactic acid bacteria,” in Lactic Acid Bacteria, eds S. Salminen and A. von Wright (New York: Marcel Dekker Inc), 385–436.

Google Scholar

Katoh, K., and Standley, D. M. (2013). MAFFT multiple sequence alignment software version 7: improvements in performance and usability. Mol. Biol. Evol. 30, 772–780. doi: 10.1093/molbev/mst010

PubMed Abstract | CrossRef Full Text | Google Scholar

Laleye, L. C., Simard, R. E., Lee, B. H., and Holley, R. A. (1989). Control of heterofermentative bacteria during cheddar cheese maturation by homofermentative Lactobacillus Starters. J Dairy Sci. 72, 3134–3142. doi: 10.3168/jds.s0022-0302(89)79471-6

CrossRef Full Text | Google Scholar

Madera, C., Monjardin, C., and Suarez, J. E. (2004). Milk contamination and resistance to processing conditions determine the fate of Lactococcus lactis bacteriophages in dairies. Appl. Environ. Microbiol. 70, 7365–7371. doi: 10.1128/aem.70.12.7365-7371.2004

PubMed Abstract | CrossRef Full Text | Google Scholar

Mandal, S., Van Treuren, W., White, R. A., Eggesbø, M., Knight, R., and Peddada, S. D. (2015). Analysis of composition of microbiomes: a novel method for studying microbial composition. Microb. Ecol. Health Dis. 26:27663. doi: 10.3402/mehd.v26.27663

PubMed Abstract | CrossRef Full Text | Google Scholar

Michael, W., and Mullan, A. (2000). Causes and control of early gas production in cheddar cheese. Int. J. Dairy Tech. 53, 63–68. doi: 10.1111/j.1471-0307.2000.tb02660.x

CrossRef Full Text | Google Scholar

Morandi, S., Battelli, G., Silvetti, T., Goss, A., Cologna, N., and Brasca, M. (2019). How the biodiversity loss in natural whey culture is affecting ripened cheese quality? The case of Trentingrana cheese. LWT Food Sci. Tech. 115:108480. doi: 10.1016/j.lwt.2019.108480

CrossRef Full Text | Google Scholar

Moser, A., Schafroth, K., Meile, L., Egger, L., Badertscher, R., and Irmler, S. (2018). Population dynamics of Lactobacillus helveticus in Swiss Gruyère-type cheese manufactured with natural whey cultures. Front. Microbiol. 9:637. doi: 10.3389/fmicb.2018.00637

PubMed Abstract | CrossRef Full Text | Google Scholar

Oliveira, J., Mahony, J., Hanemaaijer, L., Kouwen, T. R. H. M., and van Sinderen, D. (2018). Biodiversity of bacteriophages infecting Lactococcus lactis starter cultures. J Dairy Sci 101, 96–105. doi: 10.3168/jds.2017-13403

PubMed Abstract | CrossRef Full Text | Google Scholar

Price, M. N., Dehal, P. S., and Arkin, A. P. (2009). FastTree: computing large minimum evolution trees with profiles instead of a distance matrix. Mol. Biol. Evol. 26, 1641–1650. doi: 10.1093/molbev/msp077

PubMed Abstract | CrossRef Full Text | Google Scholar

Pujato, S., Quiberoni, A., and Mercanti, D. (2019). Bacteriophages on dairy foods. J. Appl. Microbiol. 126, 14–30. doi: 10.1111/jam.14062

PubMed Abstract | CrossRef Full Text | Google Scholar

Rodriguez-Valera, F., Martin-Cuadrado, A.-B., Rodriguez-Brito, B., Pašiæ, L., Thingstad, T. F., Rohwer, F., et al. (2009). Explaining microbial population genomics through phage predation. Nat. Rev. Microbiol. 7, 828–836. doi: 10.1038/nrmicro2235

PubMed Abstract | CrossRef Full Text | Google Scholar

Rossetti, L., Fornasari, M. E., Gatti, M., Lazzi, C., Neviani, E., and Giraffa, G. (2008). Grana Padano cheese whey starters: microbial composition and strain distribution. Int. J. Food Microbiol. 127, 168–171. doi: 10.1016/j.ijfoodmicro.2008.06.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Rossetti, L., and Giraffa, G. (2005). Rapid identification of dairy lactic acid bacteria by M13-generated, RAPD-PCR fingerprint databases. J. Microbiol. Methods 63, 135–144. doi: 10.1016/j.mimet.2005.03.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmid, M., Muri, J., Melidis, D., Varadarajan, A. R., Somerville, V., Wicki, A., et al. (2018). Comparative genomics of completely sequenced Lactobacillus helveticus genomes provides insights into strain-specific genes and resolves metagenomics data down to the strain level. Front. Microbiol. 9:63. doi: 10.3389/fmicb.2018.00063

PubMed Abstract | CrossRef Full Text | Google Scholar

Settanni, L., and Moschetti, G. (2010). Non-starter lactic acid bacteria used to improve cheese quality and provide health bene?ts. Food Microbiol. 27, 691–697. doi: 10.1016/j.fm.2010.05.023

PubMed Abstract | CrossRef Full Text | Google Scholar

Somerville, V., Lutz, S., Schmid, M., Frei, D., Moser, A., Irmler, S., et al. (2019). Long-read based de novo assembly of low-complexity metagenome samples results in finished genomes and reveals insights into strain diversity and an active phage system. BMC Microbiol. 19:143. doi: 10.1186/s12866-019-1500-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Spus, M., Li, M., Alexeeva, S., Wolkers-Rooijackers, J. C., Zwietering, M. H., Abee, T., et al. (2015). Strain diversity and phage resistance in complex dairy starter cultures. J. Dairy Sci. 98, 5173–5182. doi: 10.3168/jds.2015-9535

PubMed Abstract | CrossRef Full Text | Google Scholar

Vazquez-Baeza, Y., Pirrung, M., Gonzalez, A., and Knight, R. (2013). EMPeror: a tool for visualizing high-throughput microbial community data. Gigascience 2:16.

Google Scholar

Zago, M., Bonvini, B., Rossetti, L., Meucci, A., Giraffa, G., and Carminati, D. (2015). Biodiversity of Lactobacillus helveticus bacteriophages isolated from cheese whey starters. J. Dairy Res. 82, 242–247. doi: 10.1017/s0022029915000151

PubMed Abstract | CrossRef Full Text | Google Scholar

Zago, M., Comaschi, L., Neviani, E., and Carminati, D. (2005). Investigation on the presence of bacteriophages in natural whey starters used for the production of Italian long-ripened cheeses. Milchwissenschaft 60, 171–174.

Google Scholar

Zago, M., De Lorentiis, A., Carminati, D., Comaschi, L., and Giraffa, G. (2006). Detection and identification of Lactobacillus delbrueckii subsp. lactis bacteriophages by PCR. J. Dairy Res. 73, 146–153. doi: 10.1017/s0022029905001524

PubMed Abstract | CrossRef Full Text | Google Scholar

Zago, M., Rossetti, L., Reinheimer, J., Carminati, D., and Giraffa, G. (2008). Detection and identification of Lactobacillus helveticus bacteriophages by PCR. J. Dairy Res. 75, 196–201. doi: 10.1017/S0022029908003208

PubMed Abstract | CrossRef Full Text | Google Scholar