Impact Factor 4.235 | CiteScore 6.4
More on impact ›


Front. Microbiol., 09 July 2020 |

Complementary Use of Cultivation and High-Throughput Amplicon Sequencing Reveals High Biodiversity Within Raw Milk Microbiota

Franziska Breitenwieser1†, Etienne V. Doll2, Thomas Clavel3,4, Siegfried Scherer2,4 and Mareike Wenning2,5*
  • 1Milchprüfring Baden-Württemberg e.V., Kirchheim unter Teck, Germany
  • 2Chair for Microbial Ecology, Weihenstephan School of Life Sciences, Technische Universität München, Freising, Germany
  • 3Functional Microbiome Research Group, Institute of Medical Microbiology, RWTH University Hospital, Aachen, Germany
  • 4ZIEL – Institute for Food and Health, Technische Universität München, Freising, Germany
  • 5Bavarian Health and Food Safety Authority, Oberschleißheim, Germany

Raw milk microbiota are complex communities with a significant impact on the hygienic, sensory and technological quality of milk products. However, there is a lack of knowledge on factors determining their composition. In the present study, four bulk tank milk samples of two farms at two different time points were analyzed in detail for their microbiota using cultivation and 16S rRNA amplicon sequencing. Diversity in samples from the first time point was assessed via cultivation of 500 aerobic mesophilic bacterial isolates in each sample. A high biodiversity of 70 and 110 species per sample was determined, of which 25–28% corresponded to yet unknown taxa. The isolates were dominated by Gram-positive members of the genera Staphylococcus, Corynebacterium, Streptococcus, or Janibacter, whilst Chryseobacterium and Acinetobacter were most abundant among the Gram-negative taxa. At the second time point, samples of the same farms were analyzed via both cultivation (1,500 individual colonies each) and high-throughput 16S rRNA gene amplicon sequencing. The latter revealed a threefold higher biodiversity at the genus level, as anaerobic or fastidious species were also detected. However, cultivation identified genera not captured by sequencing, indicating that both approaches are complementary. Using amplicon sequencing, the relative abundance of a few genera was distorted, which seems to be an artifact of sample preparation. Therefore, attention needs to be paid to the library preparation procedure with special emphasis on cell lysis and PCR.


The indigenous microbiota of raw milk has re-gained attention during the last years for economic and safety reasons. Consumers increasingly prefer unpasteurized products such as raw milk cheese or raw milk bought in retail sale or directly at the farm. This increases the risk of zoonotic diseases due to the presence of pathogenic bacteria like Listeria monocytogenes, Salmonella spp.,Campylobacter spp. or Shiga toxin-producing Escherichia coli (Giacometti et al., 2013; Gould et al., 2014; Christidis et al., 2016; Artursson et al., 2018; Willis et al., 2018). In the United States, unpasteurized milk and cheese, which are consumed by a minority of the population only, caused 96% of outbreak-related illnesses because of contaminated dairy products between 2009 and 2014 (Costard et al., 2017). Health risks are largely diminished by heat treatment of milk, but spoilage as a result of heat-resistant microorganisms or enzymes can still occur. Spore-formers are responsible for the major part of spoilage in pasteurized milk (Huck et al., 2008; Schmidt et al., 2012; Doll et al., 2017). Even UHT milk, although microbiologically stable, can be indirectly affected by the raw milk microbiota. Heat-resistant peptidases primarily produced by Gram-negative bacteria such as Pseudomonas or Serratia (Marchand et al., 2009; von Neubeck et al., 2015; Machado et al., 2017) impact the sensory and textural properties of UHT milk during elongated shelf-life and even very low residual activities may lead to bitterness and age gelation (Stoeckel et al., 2016). Thus, different product types are affected by different parts of the microbiota that are not eliminated during processing. However, besides all negative implications, the raw milk microbiota is valuable and particularly appreciated, for instance for obtaining individual and flavor-rich cheeses (Callon et al., 2005; Montel et al., 2014).

To improve consumer safety, maintain tasteful product quality, achieve a longer shelf-life, and optimize the use of beneficial bacteria in manufacturing of raw milk dairy products, the complex inter-play of farming practices, animal physiology and composition of raw milk microbiota needs to be understood. This will enable the development of housing and milking strategies effective in reducing undesired bacteria whilst maintaining a stable and beneficial microbiota. Different factors influencing the milk microbial communities have been examined over the last decades (Desmasures and Gueguen, 1997; Christiansson et al., 1999; Vacheyrou et al., 2011; Mallet et al., 2012; Verdier-Metz et al., 2012; Miller et al., 2015), but these studies relied all on laborious culture-dependent methodologies and, therefore, were limited in size, sensitivity, and focus. The emergence of high-throughput sequencing techniques has paved the way for much more comprehensive investigations, which led to a multitude of analyses published during the last few years (Quigley et al., 2013; Oikonomou et al., 2014; Ganda et al., 2016; Kable et al., 2016; Bonsaglia et al., 2017; Doyle et al., 2017a, b; Rodrigues et al., 2017; Fretin et al., 2018; Li et al., 2018; Lima et al., 2018). However, as demonstrated for other microbial communities such as those in the gut of mammals, sequencing approaches do have their own limitations and cultivation is still a valuable tool (Browne et al., 2016; Clavel et al., 2016; Clooney et al., 2016). To be able to link the outcomes of new investigations with former results, it is necessary to use culture- and sequence-based approaches concomitantly, thereby revealing outcome similarities and divergences between methodologies. In addition, the power and limits of culture-based approaches to detect microbial biodiversity in raw milk was never studied in detail.

The main objective of this study was a deep analysis of raw milk biodiversity by using a combination of cultivation, focusing on aerobic taxa including food-borne pathogens and potential spoilage organisms, and high-throughput 16S rRNA gene amplicon analysis.

Materials and Methods

Collection of Milk Samples

Bulk tank milk of two farms (A and B) in Southern Germany with a herd size of about 45 animals was collected, each on two occasions. The first time, both farms were sampled simultaneously in January (samples 1 and 2). Samples were automatically taken from the bulk tank containing four to six successive milkings using the sampling device of the milk collecting vehicle, transported to the laboratory under refrigerated conditions (maximum 6°C) within 24 h and processed within 1 h after arrival. The second time, each farm was sampled independently, farm A (sample 3) 26 months and farm B (sample 4) 18 months after the first sampling. Milk was collected directly from the bulk tank containing one and two milkings, respectively. The coupling pipe was disinfected with 70% ethanol and rinsed with milk from the tank before samples were collected. Milk was transported under refrigerated conditions and analysis started 3–4 h after sample collection. Part of the milk was used immediately for the culture-dependent analysis, whereas 150 mL were frozen and stored at −80°C for molecular work.

Bacterial Counts

Total aerobic counts (TAC) of the mesophilic raw milk microbiota were determined using the spread-plate method. Decimal dilutions were plated on tryptic soy agar (TSA, Oxoid) with 10 plates per dilution and incubated at 30°C for 5 days.

Culture-Dependent Analysis of Raw Milk Microbiota

In the first approach (sample 1 and 2), 500 colonies were randomly isolated from agar plates used for TAC determination. All colonies were selected from the identical dilution step with plates showing well-separated colonies. Starting with the first petri dish, all colonies present on this dish were picked and sub-cultured on TSA. The procedure was repeated until a total of 500 colonies was obtained. Isolates were identified using FTIR spectroscopy and gene sequencing for representative isolates. In the second approach (samples 3 and 4), the same isolation procedure was followed with the exception that larger petri dishes (Ø = 15 cm) were used to apply larger sample volumes and obtain higher numbers of isolates: 1,500 colonies were selected per sample and all were identified by gene sequencing.

Identification of Isolates by FTIR Spectroscopy

Sample preparation and recording of spectra was performed as described previously (Wenning et al., 2014; von Neubeck et al., 2015). Identification of spectra was conducted using in-house reference databases containing approximately 8,000 spectra of 1,000 species using the parameters described by Kümmerle et al. (1998) and Oberreuter et al. (2002). To reduce redundancy, FTIR spectra of all isolates were compared for each sample by hierarchical cluster analysis (HCA) as described previously (von Neubeck et al., 2015). Representative isolates were selected and subsequently identified by gene sequencing.

Sanger Sequencing

For identification of most bacterial isolates, the 16S rRNA gene was used. The rpoB gene (Mellmann et al., 2006) was used for Staphylococcus spp. and the 26S rRNA gene for yeasts (Kurtzman and Robnett, 2003). Cell lysis, PCR, and sequencing were performed as described previously (von Neubeck et al., 2015). rpoB and 26S RNA gene sequences were identified using the NCBI database, 16S rRNA gene sequences using EzBiocloud (Yoon et al., 2017). Thresholds for species and genus assignment of 16S rRNA gene sequences were 98.65 (Kim et al., 2014) and 95% sequence similarity, respectively; 95% similarity was used for species assignment of rpoB sequences.

High-Throughput 16S rRNA Gene Amplicon Analysis of Raw Milk Microbiota

Milk samples 3 (farm A) and 4 (farm B) subjected to culture-based analysis of 1,500 colonies were also used for amplicon sequencing.

Concentration of Microbial Cells in Raw Milk

As microbial loads in fresh raw milk are comparatively low and did not exceed log 5 cfu/mL in samples 3 and 4, bacterial cells were concentrated prior to DNA extraction. A volume of 150 mL was centrifuged at 8,000 × g for 20 min at 4°C to pellet bacterial cells. The supernatant consisting of a fat layer and skim milk was carefully removed, except for 10 mL skim milk, which were left to re-suspend the pellet. The suspension was transferred to a 50 mL tube and additionally centrifuged at 5,000 × g for 10 min at 4°C. Again, the fat layer was carefully removed.

Removal of Sedimented Casein

Sedimented casein was removed according to Murphy et al. (2002). The cell suspension was divided into 10 aliquots of 1 mL and 300 μL 0.5 M EDTA (pH 8.0) as well as 200 μL TE buffer (pH 7.6) were added to each aliquot. Within approximately 1 min, the casein micelles disintegrated due to chelating of calcium ions indicated by clarification of the solution, which was centrifuged at 16,000 × g for 1 min at room temperature. The supernatant was discarded and the pellet re-suspended in 100 μL Ringer’s solution. All ten subsamples were then pooled and centrifuged again. The supernatant was discarded and the pellet re-suspended in 350 μl quarter strength Ringer’s solution.

DNA Extraction and Enrichment of Bacterial DNA

The PathoProof DNA Extraction Kit (Thermo Fisher Scientific), including enzymatic lysis of microbial cells, was applied according to the manufacturer’s instructions to extract DNA, producing a final volume of 100 μL DNA extract per sample. As bulk tank milk contains high amounts of eukaryotic DNA originating from somatic cells of the cow, enrichment of bacterial DNA using the Looxster Enrichment Kit (Analytik Jena, Germany) was performed according to the manufacturer’s instructions resulting in 30 μL of DNA for each sample.

Determination of Bacterial DNA Concentrations

The concentration of total DNA was determined using the Qubit® 2.0 fluorometer (Life Technologies Co.) using the dsDNA HS Assay Kit. In a further step, the concentration of bacterial DNA in both extracts was determined using quantitative real-time PCR. PCR mixtures contained 4 μL Phusion® Buffer HF (Thermo Fisher Scientific Inc.), 1 μL dNTPs (20 nM), 2 μL each of primer 515F (5′-GTGCCAGCMGCGCGGTAA) and 806R (5′-GGACTACHVGGGTWTCTAAT) (10 pmol/μL), 0.1 μL Phusion High-Fidelity DNA Polymerase (Thermo Fisher Scientific Inc.), 1 μL SYBR-Green (diluted 1:50,000), and 2 μL of DNA extract (diluted 1:10) in a final volume of 20 μL. After initial denaturation at 98°C for 30 s 40 cycles of 5 s denaturation at 98°C, 10 s primer annealing at 52.5°C, and 10 s elongation at 72°C were run. A DNA standard composed of somatic DNA extracted from bulk tank milk and different fractions (10, 1, and 0.1%) of a bacterial DNA extract of a Pseudomonas and Enterococcus raw milk isolate was used for quantification of bacterial DNA.

Library Preparation for Amplicon Sequencing

DNA quantification resulted in a total DNA concentration of 15 ng/μL for sample 3 of farm A and 10 ng/μL for sample 4 of farm B. Sample 3 contained 4.8% (0.7 ng/μL) bacterial DNA, sample 4 2.5% (0.25 ng/μL). DNA extracts were then standardized to contain 0.2 ng/μL bacterial DNA. Library preparation followed a two-step protocol (Berry et al., 2011) based on amplification of the V3–V4 region of the 16S rRNA gene using primers 341F (5′-CCTACGGGNGGCWGCAG) and 785R (5′-GGATTAGATACCCBDGTAGTC) (Klindworth et al., 2013). Due to the low concentrations of template microbial DNA, the PCR protocol was adjusted and the number of replicates as well as the number of cycles augmented. For the first PCR step, eight parallel PCRs (four barcodes in duplicate) were run for each sample, each reaction using 3.5 μL DNA extract containing 0.7 ng bacterial DNA. After initial denaturation at 98°C for 30 s, 30 cycles of 5 s denaturation at 98°C, 10 s primer annealing at 55°C, and 10 s elongation at 72°C were run. The second PCR step added barcodes (dual combinatorial indexing) and Illumina adaptors to the amplified fragments and was done using the same protocol but 2 μl of PCR product of step 1 as template and was run only for 10 additional cycles. A negative control (PCR blank) was included using PCR-grade water as template.

Purification of Amplified DNA and Sequencing

Duplicate PCRs for any given samples (i.e., having the same barcodes) were transferred into a 1.5 mL tube and purified using Agencourt AMPure XP Beads (Beckman Coulter Inc.). DNA concentrations were determined via QubitTM 2.0 fluorometer and adjusted to 2 nM. All libraries were pooled and sequencing was performed in paired-end mode (2 × 275 cycles) using a MiSeq platform (Illumina).

Analysis of Sequencing Data

Raw reads were processed using the IMNGS pipeline (Lagkouvardos et al., 2016) based on UPARSE (Edgar, 2013). After demultiplexing, forward and reverse reads were merged and trimmed by five nucleotides on each end. Chimera filtering was done using UCHIME (Edgar et al., 2011). Operational Taxonomic Units (OTU) clustering (97% identity) was performed by USEARCH 8.0 (Edgar, 2010) and OTUs occurring at a relative abundance <0.05% in all samples were discarded. Taxonomical identification was assigned using the RDP classifier (Wang et al., 2007) and diversity analyses were done in Rhea (Lagkouvardos et al., 2017). Relative abundances of OTUs were normalized to account for differences in sequence depth and α-diversity was assessed on the basis of species richness as well as Shannon and Simpson diversity indices. As automated identification may be imprecise and often terminates at higher ranks, 487 OTU sequences not assigned to a genus were manually identified at the genus level (95% similarity) using EzBiocloud (Yoon et al., 2017). OTUs for manual identification were selected according to the following criteria: all OTUs that were identified only to the phylum level or above, all OTUs belonging to families detected by culturing, and the 100 most abundant OTUs in each sample. ß-diversity was computed for both the culture-dependent and -independent results based on generalized UniFrac distances (Chen et al., 2012) and visualized in a phylogram based on the Ward’s minimum variance method (Murtagh and Legendre, 2014).


Bulk tank milk of two farms in Southern Germany was analyzed in detail for microbial biodiversity. In a first culture-dependent approach, 500 colonies per sample were isolated and identified. In a second approach, milk of the same farms was analyzed based to 1,500 colonies per sample combined with high-throughput sequencing.

Culture-Dependent Analysis of 500 Isolates

The sample from farm A (sample 1) contained six, the one of farm B (sample 2) four successive milkings and both had a TAC of 4.7 log cfu/mL. We identified 495 (sample 1) and 472 (sample 2) isolates to the species level revealing a high biodiversity. Sample 2 contained 110 species assigned to 57 genera, whereas sample 1 was a little less diverse with 71 species belonging to 36 genera (Figure 1 and Table 1). Most of the diversity was represented by only a minor fraction of isolates, as the majority of species occurred with a relative abundance <1%. In sample 1, 49 species (69%) were represented by 16% of isolates and in sample 2, 83 species (75%) accounted for 29% of isolates (Figure 1 and Supplementary Tables S1, S2). Both microbiota were dominated by Gram-positive taxa (Figure 1 and Supplementary Tables S1, S2). Isolates from Gram-negative species accounted for less than 10% in sample 1 and little more than 20% in sample two. Many species in sample 1 belonged to Staphylococcus and Corynebacterium (8 and 10, respectively, Supplementary Table S1), many of which were among the most abundant species. In sample 2, 14 species belonged to Corynebacterium, but the genus with the second most species was Chryseobacterium (Supplementary Table S2). The fact that 40% of all species were detected with only a single isolate (relative abundance 0.2%) indicates that the biodiversity present was insufficiently covered and a substantial fraction of species remained undiscovered.


Figure 1. Composition of the microbiota in two bulk tank milk samples based on each 500 isolates. Only species with a relative abundance >1% of the total isolates are listed. Species occurring in only one sample are indicated in green or blue, species detected in both samples are given in red. S1, sample 1; S2, sample 2; Brachyb., Brachybacterium; Bro., Brochotrix; C., Corynebacterium; Chr., Chryseobacterium; D., Dietzia; I., Isoptericola; J., Janibacter; Ln., Leuconostoc; Lc., Lactococcus; M., Microbacterium; O., Ornithinimicrobium; P., Propionibacterium; St., Staphylococcus; Str., Streptococcus; sp. nov., hitherto undescribed species.


Table 1. Microbial biodiversity of samples 1–4 based on the cultivation of isolates on TSA.

Culture-Dependent Analysis of 1,500 Isolates

In order to capture a larger fraction of the high biodiversity present in bulk tank milk of these two farms, a second in-depth aerobic cultivation approach was undertaken (1,500 isolates per sample). Milk of farm A (sample 3, one milking) was collected 26 months later than sample 1, that of farm B (sample 4, two milkings) 18 months later than sample 2. Sample 3 had a TAC of 4.3 log cfu/mL and sample 4 contained 4.4 log cfu/mL (Table 1). It turned out that sample 4 comprised a high fraction of yeasts (27% of isolates) mainly belonging to the genera Candida and Yarrowia. After eliminating these from the data set, 1,102 bacterial isolates remained.

As expected, the increased sample size covered a higher biodiversity (Figure 2 and Supplementary Tables S3, S4). Samples 1 and 3 from farm A showed a nearly identical slope in rarefaction analysis (Figure 2), but 58% more species were detected in sample 3 because of the higher depth of analysis. The slope of sample 4 from farm B (1,100 colonies) was higher than for samples from farm A indicating a higher biodiversity, which is also expressed by species richness as well as Shannon and Simpson index (Table 1). However, the slope of sample 4 (only bacterial isolates, no yeast) was considerably lower than that of sample 2 (500 colonies). Although the number of isolates identified was more than doubled, only 14% more species were detected. While sample 3 was dominated by staphylococci as already observed for sample 1 (Supplementary Table S3), a shift in the dominating genera in sample 4 was recognized (Supplementary Table S4). Beside the highly abundant yeasts (data not shown), Microbacterium was predominant (25.6% of isolates), whereas it represented only 5.2% in sample 2. Further genera with high relative abundance in both samples were Corynebacterium, Streptococcus, and Aerococcus. Pseudoxanthomonas and Stenotrophomonas were detected as the most prevalent Gram-negative genera.


Figure 2. Rarefaction analysis of bulk tank milk microbiota of farm A (blue) and B (green) analyzed by cultivation. Solid lines, samples 1 and 2 (500 isolates); dotted lines, samples 3 (1,500 isolates), and 4 (1,100 isolates).

As the number of species with a relative abundance >1% remained constant (farm A) or even dropped slightly (farm B, Table 1) the increase in biodiversity was exclusively found among the rare species with a relative abundance <1%. Accordingly, the Simpson index did not increase with the larger sampling depth and rarefaction curves still display considerable slope (Figure 2), indicating that much more individuals need to be identified in order to better cover the detectable biodiversity.

Molecular Analysis Based on 16S rRNA Gene Amplicon Sequencing

To study the biodiversity of microbial communities in milk samples 3 and 4 more comprehensively, a cultivation-independent approach was also followed. In order to test for reproducibility of the PCR-based library preparation, four replicates of each sample (same DNA extraction) were analyzed, each replicate consisting of eight single PCR reactions pooled together. For sample 4 (farm B), all four replicates worked and yielded a similar number of reads; for sample 3 (farm A) one replicate failed in sequencing. After data pre-processing, 30,000 high-quality reads were obtained for sample 3 and nearly 42,000 reads for sample 4. The reproducibility was very high and only minor variations were observed at the compositional level between different replicates of the same sample (Figure 3), demonstrating the robustness of the analysis.


Figure 3. Genus level composition of samples 3 and 4 determined by amplicon sequencing. The DNA extract of each sample was analyzed by 4 PCR replicates (one of sample 3 failed in sequencing), which are shown separately. Genera are sorted in alphabetical order (bottum-up), reads with no genus assignment are shown at the very top in gray. For clarity, only genera with a relative abundance >1.5% in at least one replicate are listed; missing fractions to 100% are genera with a relative abundance <1.5%.

The two farm samples displayed a similar biodiversity, both concerning the number of OTUs and genera identified (Table 2). As many as 144 genera and 75% of the 624 OTUs detected occurred in both samples. A considerable number of reads (33% in sample 3 and 21% in sample 4) belonging to 287 OTUS were not assigned to a genus using the RDP classifier. Manual identification of OTUs using EzBiocloud was able to classify nearly a quarter of those (69 OTUs). Nonetheless, 23 and 14% of reads were identified only to the level of phylum, order or family (Figure 3).


Table 2. Comparison of microbial biodiversity of samples 3 and 4 as detected by 16S rRNA amplicon sequencing and cultivation.

Although the microbiota resembled each other in terms of overall taxa occurrence, there were marked discrepancies in relative abundances, already at the phylum level. The most abundant genus in sample 4 was Caryophanon (15%) belonging to the Firmicutes that was only rarely detected in sample 3 (0.3%). Sample 3, in contrast, contained higher fractions of Chryseobacterium (Bacteroidetes), Janibacter (Actinobacteria), and Staphylococcus (Firmicutes) whereas sample 4 harbored more Kocuria (Actinobacteria) as well as Aerococcus and Streptococcus (Firmicutes). Altogether, sample 3 was characterized by a higher fraction of Bacteroidetes, whereas sample 4 contained more Firmicutes. Still, some abundant genera such as Bacillus, Brachybacterium, Clostridium, Corynebacterium, Facklamia, Jeotgalicoccus, and Pseudomonas occurred in similar fractions in both samples (Figure 3).

Comparison of Culture- and Sequencing-Based Biodiversity

The results obtained by the two different approaches exhibit many discrepancies (Figures 4, 5 and Table 2). Amplicon sequencing detected a higher biodiversity, as 4–5 times more OTUs than species, 2.8–3.4 times more genera and 2–3 times more families were identified. This is also reflected by higher Shannon and Simpson indices (Table 2). However, there were only few highly abundant genera (>1%) and this favored the detection of a large biodiversity. Additionally, they sum up to a relatively small fraction of reads (53 and 66%, respectively), whereas in the cultivation-dependent approach they represented 90% of isolates. This allows for the discovery of many more minor genera and agrees with the trend already observed for the cultivation-dependent analysis of samples 2 and 4. In addition, sample 3 exhibited a higher biodiversity than sample 4 in the amplicon-based approach, although fewer reads were obtained (Table 2). This clearly contradicts the results obtained by the cultivation-dependent analysis. Rarefaction curves of amplicon analyses flattened for both samples indicating that the largest fraction of biodiversity has been detected within the approximately 10,000 reads of one PCR replicate (Figure 4). However, combining the three or four replicate PCRs for each sample still increased the biodiversity by 10% for sample 3 (3 replicates) and 17% for sample 4 (4 replicates).


Figure 4. Rarefaction analysis of microbiota of bulk tank milk samples 3 (blue) and 4 (green) analyzed by a cultivation-based isolation (3_I, 4_I) or amplicon sequencing (3_A, 4_A). Individual PCR replicates are shown for the amplicon approach.


Figure 5. Comparison of microbiota in samples 3 (farm A) and 4 (farm B) as detected by the cultivation (3-I, 4-I) and molecular (3-A, 4-A) approaches. (A) Relative abundance of families. Only families >1% in at least one analysis are listed. i.s., incertae sedis; n.d., not detected. (B) Relative abundance of phyla. C., Candidatus. (C) Phylogram based on the Ward’s minimum variance method showing the hierarchical clustering of samples. Genus-level data were used to calculate distances; only the 57 genera detected at least once by both the cultivation-dependent and –independent analyses were included.

Many genera not detected using the culture-dependent approach belonged to families that grow anaerobically or are not culturable using TSA under aerobic conditions such as Lachnospiraceae, Rhodocyclaceae, Ruminococcaceae, Saccharibacteria incertae sedis, Sphaerobacteraceae, or Clostridiales incertae sedis (Supplementary Table S5). These missing taxa lead to the much lower slope of rarefaction curves of the culture-based analysis (Figure 4). The biggest discrepancy was observed in sample 4, where Caryophanon belonging to the Planococcaceae was detected as the most abundant genus by amplicon sequencing accounting for almost 15% of reads, but was not isolated a single time in the culture-dependent analysis due to its fastidious growth requirements. However, although the sequencing approach revealed a largely higher biodiversity, there were still 24 genera not detected, but found by culturing. Two thirds of these belonged to the Actinobacteria and 21% to the Firmicutes (Supplementary Table S6). Eight genera belonged to six families that were not identified by amplicon sequencing. Jonesiaceae and Morganellaceae were detected exclusively by culturing in sample 3, Mycobacteriaceae and Tsukamurellaceae only in sample 4 and in both samples Paenibacillaceae and Propionibacteriaceae were overlooked in the molecular analysis. In summary, none of the approaches alone is able to cover and detect the full biodiversity present in the samples.

It was furthermore observed that the relative abundance of some families is massively distorted depending on the type of analysis (Figure 5A). Microbacteriaceae and Xanthomonadaceae were considerably underrepresented in the amplicon approach in both samples as well as Staphylococcaceae or Gordoniaceae in sample 3, whereas Flavobacteriaceae or Intrasporangiaceae were overrepresented. The misestimation of Microbacteriaceae and Gordoniaceae (Actinobacteria) as well as Staphylococcaceae and Planococcaceae (Firmicutes) lead to substantial differences observed already at the phylum level (Figure 5B). To correlate the similarities between the two samples and both techniques and to visualize methodology-driven discrepancies, a phylogram was calculated (Figure 5C). It was based on those 57 genera that have been determined at least once using each approach in order to emphasize the distribution of relative abundances rather than the differences in biodiversity. It thus included 50 and 59% of reads of sample 3 and 4, respectively, the remaining OTUs were not considered. Data obtained by amplicon sequencing showed the highest similarities and the culture-based analyses were more distantly related. Differences introduced by methodology, thus, exceeded sample-specific differences. But besides all discrepancies observed, many of those families identified in a sample by both techniques were captured with relatively similar fractions (Figure 5A and Supplementary Table S5) revealing substantial overlap in results between the two techniques.


The four milk samples analyzed in this study contained a microbiota typical for fresh bovine raw milk dominated by Gram-positive bacteria (Delbes et al., 2007; Fricker et al., 2011; Curone et al., 2018; Lima et al., 2018), although the bulk tank milk of samples 1 and 2 had been stored at the farm for 3 or 2 days before collection (combining six or four milkings, respectively), which is thought to favor the growth of Gram-negative psychrotolerant genera such as Pseudomonas. The major phyla detected were Firmicutes, Actinobacteria, Proteobacteria, and Bacteroidetes as already found for bulk tank milk or milk from single healthy animals in other studies (Ganda et al., 2016; Rodrigues et al., 2017; Fretin et al., 2018; Lima et al., 2018). Kable et al. (2016) analyzed the microbiome of almost 900 samples of tanker trucks covering 200 farms in California and found 29 taxa belonging to the core microbiome present in all samples. 18 were identified down to genus level, the remaining eleven only to the family level. Of the 18 taxa identified at genus level four (Butyrivibrio, Yaniella, Mycoplasma, 5-7N15) were absent in samples 3 and 4 in this study, indicating that regional characteristics may influence the milk microbiota and detected core microbiomes will be significantly affected by various factors, some of which are probably still unknown. Numerous studies highlighted the influence of farm practices, herd management, milking hygiene (Verdier-Metz et al., 2009; Vacheyrou et al., 2011; Mallet et al., 2012; Miller et al., 2015; Fretin et al., 2018; Metzger et al., 2018), season (Doyle et al., 2017b; Li et al., 2018), lactation stage (Doyle et al., 2017a; Curone et al., 2018) or breed (Curone et al., 2018) on the bovine milk microbiome, but many more studies are needed to elucidate the complex relationships between different factors associated with milk production and milk microbiota. Because of this multifactorial influence, outcomes of many studies regarding species composition differ largely. In the present study, we found Staphylococcus, Corynebacterium, Caryophanon, or Streptococcus as leading genera in the amplicon approach, whereas Li et al. (2018) found the four genera Bacillus, Acinetobacter, Pseudomonas, and Lactococcus making up more than half of all reads.

With the advance of next-generation sequencing techniques, huge data sets can be generated and very comprehensive studies have become feasible. To be able to relate these new data to outcomes of former studies based on the cultivation and identification of isolates, comparative analyses are necessary. To ensure comparability, an adequate depth of culture-based analysis is required, which is laborious. In the first part of the present study, a large biodiversity could be detected in bulk tank milk samples using a culture-dependent analysis. Although only the aerobic mesophilic bacteria were analyzed, more than 100 species were identified out of 500 isolates in sample 2. For samples 3 and 4, dedicated to analysis by both culture-dependent and –independent techniques, the number of isolates was increased to 1,500, restricting the analysis to the aerobic bacteria growing on TSA in order to strengthen the validity for taxa related to food quality and safety and to retain a high sampling depth for the chosen culture condition. In contrast, 16S rRNA gene amplicon sequencing detects prokaryotic taxa regardless of growth requirements and, in addition, was based on a larger sample volume. For DNA extraction 150 mL milk were centrifuged to concentrate bacterial cells. In the culture dependent approach, however, the 1,500 colonies picked from agar plates were contained in <0.1 mL of sample. This difference is reflected by the three times higher number of genera identified among sequence reads vs. isolates. As the sample preparation did not include a step using propidium monoazide (PMA) to block DNA from dead cells, part of the detected biodiversity may not belong to the viable fraction of the microbiota. However, Kable et al. (2019) found equivalent cell numbers as well as no detectable differences in biodiversity between raw milk samples with and without PMA treatment. We therefore conclude that raw milk contains only a minor fraction of dead cells which does not significantly contribute to the biodiversity detected.

Nevertheless, species richness in sequencing data may be overestimated, as 54 of all OTUs (8.6%) were classified only to the phylum level or above. These unknown taxa detected may indeed exist, or may represent sequencing artifacts, which can never be completely removed during processing (Edgar, 2013; Allali et al., 2017). Increasing threshold levels for relative abundance filtering would contribute to diminishing this problem, but would also exclude real taxa. In the present study, a relative abundance threshold of 0.1% would have detected 29 genera less than the 0.05% filter threshold applied. This equals 16% of all genera detected. Six of those were also detected by culturing and, thus, are evidentially no artifacts. Diakite et al. (2019) reported similar findings for the analysis of fecal microbiota. They detected 27 bacterial species using cultivation that were represented by less than five sequence reads at an average sequencing depth of approximately 100,000 reads per sample (equaling < 0.005% of reads). Relative abundance filtering should therefore be combined with an additional parameter such as identification to a certain taxonomic level to improve efficiency of filtering artifacts.

To the contrary, 24 genera were not detected by sequencing, but found using cultivation. The majority of these culture-specific taxa were of very low relative abundance and often represented by only one isolate. However, some occurred more often such as Agromyces, Citricoccus, Jonesia, Luteococcus, Tessaracoccus, and Mycobacterium, representing 0.35–0.75% of all isolates in the respective sample. The high sensitivity of culture-based approaches concerning the detection of species has also been shown by other studies. Diakite et al. (2019) applied a large-scale culture approach for fecal samples using 58 different conditions and detected 494 species while sequencing was able to uncover merely 362 species (73%). Analyzing the taxonomic assignment of OTUs in the present study revealed that Actinobacteria were overrepresented among those 24 genera not detected (66%) while they made up only 44% of genera found by culturing and 22% in amplicon sequencing. Three of the 24 genera (Brevibacillus, Paenibacillus, and Rummeliibacillus) belonged to spore forming bacteria and DNA from endospores will remain undetected in molecular analyses.

It is very likely that missing or underestimating genera is at least partially methodology driven, choice of primers and DNA extraction from Gram-positive bacteria probably being the most critical steps (Klindworth et al., 2013; Parada et al., 2016; Panek et al., 2018; Vaidya et al., 2018). We found that four of the 24 taxa that were not detected (among them Luteococcus and Tessaracoccus) as well as Microbacterium, the most underestimated genus in the amplicon-based approach, exhibited a mismatch in the reverse primer pairing. All five genera belonged to the families Microbacteriaceae or Propionibacteriaceae (Supplementary Table S6) and showed the same substitution at the third position of the 785R primer (C instead of A). However, all other 20 undetected genera as well as Gordonia and Staphylococcus, that were greatly underestimated, showed a perfect primer match.

Numerous studies confirm that the method of DNA extraction leads to shifts in the microbiome most probably due to differential lysis of various taxa (Biesbroek et al., 2012; Hart et al., 2015; Panek et al., 2018; Vaidya et al., 2018). Actinobacteria are particularly challenging in cell lysis due to their strong cell walls. In accordance to the present study, Fretin et al. (2018) also missed many Actinobacteria and particularly Microbacterium species when comparing culture-dependent and -independent analyses of raw milk samples. Similar findings were observed by Treven et al. (2019) for human milk microbiota. Besides, Laursen et al. (2017) reported a negative correlation of detected relative abundance with genomic G+C content, which was largely diminished by elongating the denaturation step in the PCR from 30 to 120 seconds. In the present study, the denaturation step lasted 30 s and could also be a reason why Microbacterium and other genera belonging to the Actinobacteria were greatly underestimated.

Another reason why some taxa were missed may be their resistance to lytic enzymes that renders them difficult to lyse. This is well known for many species e.g., belonging to the Firmicutes such as staphylococci and often originates from modifications of the peptidoglycan structures (Davis and Weiser, 2011). The PathoProof Mastitis Kit used for DNA extraction in the present study does not include a bead-beating step, which strengthens the dependency on enzymatic lysis. In the end, there are multiple reasons for discrepancies between cultivation-dependent and –independent analyses, which may, however, be at least partially overcome when eliminating possible sources of deficiencies in the sample preparation steps by appropriately adaptating the protocol for DNA extraction and library preparation, e.g., by adding a bead-beating step.

Species identified by culture covered between 20 and 24% of OTUs detected in this study and represented only a minor fraction of biodiversity. During the last years, several studies demonstrated that the cultured fraction of microbiota can be substantially improved by applying numerous culture media and conditions in parallel and increasing sampling effort, referred to as culturomics (Lagier et al., 2012) or culture-enriched molecular profiling (Lau et al., 2016) and may even largely exceed the number of OTUs detected by culture-independent sequencing (Lau et al., 2016; Diakite et al., 2019). Besides the complementary nature of the two approaches in estimating biodiversity, capturing unknown species by culture is necessary to allow for species characterizations, to reveal possible interactions in the microbiome or to determine implications for food processing (in case of food microbiomes). The core microbiome of bovine milk analyzed by Kable et al. (2016) using amplicon sequencing included almost 40% of unclassified taxa. However, the relevance of many of these remains unclear. Our study is the first using an intense and in-depth culture-based analysis of milk microbiota. Although only one common cultivation condition was applied, numerous potential novel species and genera have been detected. Across all four samples between 20 and 30% of species are hitherto unknown and they summed up to 15.5 and 23% of all isolates identified in samples 1 and 2, respectively. Culturing is therefore still a precious and essential effort in microbial taxonomy and ecology.


Raw milk microbiota contain a high biodiversity that can be relatively well covered by classical aerobic cultivation as long as hundreds of isolates are analyzed. Cultivation is valuable to unravel and characterize unknown taxa and make a deposit in culture collections to guarantee broad availability. Amplicon sequencing, in contrast, provides a higher resolution of analysis, but may also include the analysis of artifacts due to methodological bias. Distortion of relative abundance estimation of single genera was observed in the present study, which increases the risk of misinterpretations. Hence, protocols for library preparation need to be optimized, with special attention to cell lysis of Gram-positive bacteria and PCR conditions, and bioinformatic processing must be carefully thought through. As culture-based analysis of the raw milk microbiome including anaerobic and fastidious species is laborious, further comprehensive studies will require both technical improvements to enhance throughput and the use of molecular techniques.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The amplicon sequences of this study are available in the Sequence Read Archive under Bioproject PRJNA622551.

Author Contributions

FB and ED carried out experiments in the lab. FB performed identification of Sanger and amplicon sequences. TC supported evaluation of amplicon data and helped with manuscript writing. MW and SS designed the study. SS corrected the manuscript. MW was responsible for data evaluation and wrote the manuscript. All authors contributed to the article and approved the submitted version.

Conflict of Interest

FB was employed by Milchprüfring Baden-Württemberg e.V. during the course of the study.

The remaining 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.


We thank Angela Sachsenhauser and Caroline Ziegler from the ZIEL Core Facility Microbiome for technical assistance with high-throughput amplicon sequencing.

Supplementary Material

The Supplementary Material for this article can be found online at:


Allali, I., Arnold, J. W., Roach, J., Cadenas, M. B., Butz, N., Hassan, H. M., et al. (2017). A comparison of sequencing platforms and bioinformatics pipelines for compositional analysis of the gut microbiome. BMC Microbiol. 17:194. doi: 10.1186/s12866-017-1101-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Artursson, K., Schelin, J., Thisted Lambertz, S., Hansson, I., and Olsson Engvall, E. (2018). Foodborne pathogens in unpasteurized milk in Sweden. Int. J. Food Microbiol. 284, 120–127. doi: 10.1016/j.ijfoodmicro.2018.05.015

PubMed Abstract | CrossRef Full Text | Google Scholar

Berry, D., Ben Mahfoudh, K., Wagner, M., and Loy, A. (2011). Barcoded primers used in multiplex amplicon pyrosequencing bias amplification. Appl. Environ. Microbiol. 77, 7846–7849. doi: 10.1128/AEM.05220-11

PubMed Abstract | CrossRef Full Text | Google Scholar

Biesbroek, G., Sanders, E. A., Roeselers, G., Wang, X., Caspers, M. P., Trzcinski, K., et al. (2012). Deep sequencing analyses of low density microbial communities: working at the boundary of accurate microbiota detection. PLoS One 7:e32942. doi: 10.1371/journal.pone.0032942

PubMed Abstract | CrossRef Full Text | Google Scholar

Bonsaglia, E. C. R., Gomes, M. S., Canisso, I. F., Zhou, Z., Lima, S. F., Rall, V. L. M., et al. (2017). Milk microbiome and bacterial load following dry cow therapy without antibiotics in dairy cows with healthy mammary gland. Sci. Rep. 7:8067. doi: 10.1038/s41598-017-08790-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Browne, H. P., Forster, S. C., Anonye, B. O., Kumar, N., Neville, B. A., Stares, M. D., et al. (2016). Culturing of ‘unculturable’ human microbiota reveals novel taxa and extensive sporulation. Nature 533, 543–546. doi: 10.1038/nature17645

PubMed Abstract | CrossRef Full Text | Google Scholar

Callon, C., Berdague, J. L., Dufour, E., and Montel, M. C. (2005). The effect of raw milk microbial flora on the sensory characteristics of Salers-type cheeses. J. Dairy Sci. 88, 3840–3850. doi: 10.3168/jds.S0022-0302(05)73069-1

CrossRef Full Text | Google Scholar

Chen, J., Bittinger, K., Charlson, E. S., Hoffmann, C., Lewis, J., Wu, G. D., et al. (2012). Associating microbiome composition with environmental covariates using generalized UniFrac distances. Bioinformatics 28, 2106–2113. doi: 10.1093/bioinformatics/bts342

PubMed Abstract | CrossRef Full Text | Google Scholar

Christiansson, A., Bertilsson, J., and Svensson, B. (1999). Bacillus cereus spores in raw milk: factors affecting the contamination of milk during the grazing period. J. Dairy Sci. 82, 305–314. doi: 10.3168/jds.S0022-0302(99)75237-9

CrossRef Full Text | Google Scholar

Christidis, T., Pintar, K. D., Butler, A. J., Nesbitt, A., Thomas, M. K., Marshall, B., et al. (2016). Campylobacter spp. prevalence and levels in raw milk: a systematic review and meta-analysis. J. Food Prot. 79, 1775–1783. doi: 10.4315/0362-028X.JFP-15-480

PubMed Abstract | CrossRef Full Text | Google Scholar

Clavel, T., Lagkouvardos, I., and Hiergeist, A. (2016). Microbiome sequencing: challenges and opportunities for molecular medicine. Expert Rev. Mol. Diagn. 16, 795–805. doi: 10.1080/14737159.2016.1184574

PubMed Abstract | CrossRef Full Text | Google Scholar

Clooney, A. G., Fouhy, F., Sleator, R. D., O’Driscoll, A., Stanton, C., Cotter, P. D., et al. (2016). Comparing apples and oranges?: next generation sequencing and its impact on microbiome analysis. PLoS One 11:e0148028. doi: 10.1371/journal.pone.0148028

PubMed Abstract | CrossRef Full Text | Google Scholar

Costard, S., Espejo, L., Groenendaal, H., and Zagmutt, F. J. (2017). Outbreak-Related disease burden associated with consumption of unpasteurized cow’s milk and cheese, United States, 2009-2014. Emerg. Infect. Dis. 23, 957–964. doi: 10.3201/eid2306.151603

PubMed Abstract | CrossRef Full Text | Google Scholar

Curone, G., Filipe, J., Cremonesi, P., Trevisi, E., Amadori, M., Pollera, C., et al. (2018). What we have lost: mastitis resistance in Holstein Friesians and in a local cattle breed. Res. Vet. Sci. 116, 88–98. doi: 10.1016/j.rvsc.2017.11.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Davis, K. M., and Weiser, J. N. (2011). Modifications to the peptidoglycan backbone help bacteria to establish infection. Infect. Immun. 79, 562–570. doi: 10.1128/IAI.00651-10

PubMed Abstract | CrossRef Full Text | Google Scholar

Delbes, C., Ali-Mandjee, L., and Montel, M. C. (2007). Monitoring bacterial communities in raw milk and cheese by culture-dependent and -independent 16S rRNA gene-based analyses. Appl. Environ. Microbiol. 73, 1882–1891. doi: 10.1128/AEM.01716-06

PubMed Abstract | CrossRef Full Text | Google Scholar

Desmasures, N., and Gueguen, M. (1997). Monitoring the microbiology of high quality milk by monthly sampling over 2 years. J. Dairy Res. 64, 271–280. doi: 10.1017/S002202999600213

CrossRef Full Text | Google Scholar

Diakite, A., Dubourg, G., Dione, N., Afouda, P., Bellali, S., Ngom, I. I., et al. (2019). Extensive culturomics of 8 healthy samples enhances metagenomics efficiency. PLoS One 14:e0223543. doi: 10.1371/journal.pone.0223543

PubMed Abstract | CrossRef Full Text | Google Scholar

Doll, E. V., Scherer, S., and Wenning, M. (2017). Spoilage of microfiltered and pasteurized extended shelf life milk is mainly induced by psychrotolerant spore-forming bacteria that often originate from recontamination. Front. Microbiol. 8:135. doi: 10.3389/fmicb.2017.00135

PubMed Abstract | CrossRef Full Text | Google Scholar

Doyle, C. J., Gleeson, D., O’Toole, P. W., and Cotter, P. D. (2017a). High-throughput metataxonomic characterization of the raw milk microbiota identifies changes reflecting lactation stage and storage conditions. Int. J. Food Microbiol. 255, 1–6. doi: 10.1016/j.ijfoodmicro.2017.05.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Doyle, C. J., Gleeson, D., O’Toole, P. W., and Cotter, P. D. (2017b). Impacts of seasonal housing and teat preparation on raw milk microbiota: a high-throughput sequencing study. Appl. Environ. Microbiol. 83:e02694-16. doi: 10.1128/aem.02694-16

PubMed Abstract | CrossRef Full Text | Google Scholar

Edgar, R. C. (2010). Search and clustering orders of magnitude faster than BLAST. Bioinformatics 26, 2460–2461. doi: 10.1093/bioinformatics/btq461

PubMed Abstract | CrossRef Full Text | Google Scholar

Edgar, R. C. (2013). UPARSE: highly accurate OTU sequences from microbial amplicon reads. Nat. Methods 10, 996–998. doi: 10.1038/nmeth.2604

PubMed Abstract | CrossRef Full Text | Google Scholar

Edgar, R. C., Haas, B. J., Clemente, J. C., Quince, C., and Knight, R. (2011). UCHIME improves sensitivity and speed of chimera detection. Bioinformatics 27, 2194–2200. doi: 10.1093/bioinformatics/btr381

PubMed Abstract | CrossRef Full Text | Google Scholar

Fretin, M., Martin, B., Rifa, E., Isabelle, V. M., Pomies, D., Ferlay, A., et al. (2018). Bacterial community assembly from cow teat skin to ripened cheeses is influenced by grazing systems. Sci. Rep. 8:200. doi: 10.1038/s41598-017-18447-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Fricker, M., Skanseng, B., Rudi, K., Stessl, B., and Ehling-Schulz, M. (2011). Shift from farm to dairy tank milk microbiota revealed by a polyphasic approach is independent from geographical origin. Int. J. Food Microbiol. 145(Suppl 1), S24–S30. doi: 10.1016/j.ijfoodmicro.2010.08.025

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganda, E. K., Bisinotto, R. S., Lima, S. F., Kronauer, K., Decter, D. H., Oikonomou, G., et al. (2016). Longitudinal metagenomic profiling of bovine milk to assess the impact of intramammary treatment using a third-generation cephalosporin. Sci. Rep. 6:37565. doi: 10.1038/srep37565

PubMed Abstract | CrossRef Full Text | Google Scholar

Giacometti, F., Bonilauri, P., Serraino, A., Peli, A., Amatiste, S., Arrigoni, N., et al. (2013). Four-year monitoring of foodborne pathogens in raw milk sold by vending machines in Italy. J. Food Prot. 76, 1902–1907. doi: 10.4315/0362-028X.JFP-13-213

PubMed Abstract | CrossRef Full Text | Google Scholar

Gould, L. H., Mungai, E., and Behravesh, C. B. (2014). Outbreaks attributed to cheese: differences between outbreaks caused by unpasteurized and pasteurized dairy products, United States, 1998-2011. Foodborne Pathog. Dis. 11, 545–551. doi: 10.1089/fpd.2013.1650

PubMed Abstract | CrossRef Full Text | Google Scholar

Hart, M. L., Meyer, A., Johnson, P. J., and Ericsson, A. C. (2015). Comparative evaluation of DNA extraction methods from feces of multiple host species for downstream next-generation sequencing. PLoS One 10:e0143334. doi: 10.1371/journal.pone.0143334

PubMed Abstract | CrossRef Full Text | Google Scholar

Huck, J. R., Sonnen, M., and Boor, K. J. (2008). Tracking heat-resistant, cold-thriving fluid milk spoilage bacteria from farm to packaged product. J. Dairy Sci. 91, 1218–1228. doi: 10.3168/jds.2007-0697

PubMed Abstract | CrossRef Full Text | Google Scholar

Kable, M. E., Srisengfa, Y., Laird, M., Zaragoza, J., McLeod, J., Heidenreich, J., et al. (2016). The core and seasonal microbiota of raw bovine milk in tanker trucks and the impact of transfer to a milk processing facility. mBio 7:e00836-16. doi: 10.1128/mBio.00836-16

PubMed Abstract | CrossRef Full Text | Google Scholar

Kable, M. E., Srisengfa, Y., Xue, Z., Coates, L. C., and Marco, M. L. (2019). Viable and total bacterial populations udergo equipment- and time-dependent shifts during milk processing. Appl. Environ. Microbiol. 85:e00270-19. doi: 10.1128/AEM.00270-19

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, M., Oh, H. S., Park, S. C., and Chun, J. (2014). Towards a taxonomic coherence between average nucleotide identity and 16S rRNA gene sequence similarity for species demarcation of prokaryotes. Int. J. Syst. Evol. Microbiol. 64(Pt 2), 346–351. doi: 10.1099/ijs.0.059774-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Klindworth, A., Pruesse, E., Schweer, T., Peplies, J., Quast, C., Horn, M., et al. (2013). Evaluation of general 16S ribosomal RNA gene PCR primers for classical and next-generation sequencing-based diversity studies. Nucleic Acids Res. 41:e1. doi: 10.1093/nar/gks808

PubMed Abstract | CrossRef Full Text | Google Scholar

Kümmerle, M., Scherer, S., and Seiler, H. (1998). Rapid and reliable identification of food-borne yeasts by Fourier-transform infrared spectroscopy. Appl. Environ. Microbiol. 64, 2207–2214. doi: 10.1128/AEM.64.6.2207-2214.1998

CrossRef Full Text | Google Scholar

Kurtzman, C. P., and Robnett, C. J. (2003). Phylogenetic relationships among yeasts of the ‘Saccharomyces complex’ determined from multigene sequence analyses. FEMS Yeast Res. 3, 417–432. doi: 10.1016/S1567-1356(03)00012-6

CrossRef Full Text | Google Scholar

Lagier, J. C., Armougom, F., Million, M., Hugon, P., Pagnier, I., Robert, C., et al. (2012). Microbial culturomics: paradigm shift in the human gut microbiome study. Clin. Microbiol. Infect. 18, 1185–1193. doi: 10.1111/1469-0691.12023

PubMed Abstract | CrossRef Full Text | Google Scholar

Lagkouvardos, I., Fischer, S., Kumar, N., and Clavel, T. (2017). Rhea: a transparent and modular R pipeline for microbial profiling based on 16S rRNA gene amplicons. PeerJ 5:e2836. doi: 10.7717/peerj.2836

PubMed Abstract | CrossRef Full Text | Google Scholar

Lagkouvardos, I., Joseph, D., Kapfhammer, M., Giritli, S., Horn, M., Haller, D., et al. (2016). IMNGS: a comprehensive open resource of processed 16S rRNA microbial profiles for ecology and diversity studies. Sci. Rep. 6:33721. doi: 10.1038/srep33721

PubMed Abstract | CrossRef Full Text | Google Scholar

Lau, J. T., Whelan, F. J., Herath, I., Lee, C. H., Collins, S. M., Bercik, P., et al. (2016). Capturing the diversity of the human gut microbiota through culture-enriched molecular profiling. Genome Med. 8:72. doi: 10.1186/s13073-016-0327-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Laursen, M. F., Dalgaard, M. D., and Bahl, M. I. (2017). Genomic GC-content affects the accuracy of 16S rRNA gene sequencing based microbial profiling due to PCR Bias. Front. Microbiol. 8:1934. doi: 10.3389/fmicb.2017.01934

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, N., Wang, Y., You, C., Ren, J., Chen, W., Zheng, H., et al. (2018). Variation in raw milk microbiota throughout 12 months and the impact of weather conditions. Sci. Rep. 8:2371. doi: 10.1038/s41598-018-20862-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Lima, S. F., Bicalho, M. L. S., and Bicalho, R. C. (2018). Evaluation of milk sample fractions for characterization of milk microbiota from healthy and clinical mastitis cows. PLoS One 13:e0193671. doi: 10.1371/journal.pone.0193671

PubMed Abstract | CrossRef Full Text | Google Scholar

Machado, S. G., Bagliniere, F., Marchand, S., Van Coillie, E., Vanetti, M. C., De Block, J., et al. (2017). The biodiversity of the microbiota producing heat-resistant enzymes responsible for spoilage in processed bovine mlk and dairy products. Front. Microbiol. 8:302. doi: 10.3389/fmicb.2017.00302

PubMed Abstract | CrossRef Full Text | Google Scholar

Mallet, A., Guéguen, M., Kauffmann, F., Chesneau, C., Sesboué, A., and Desmasures, N. (2012). Quantitative and qualitative microbial analysis of raw milk reveals substantial diversity influenced by herd management practices. Int. Dairy J. 27, 13–21. doi: 10.1016/j.idairyj.2012.07.009

CrossRef Full Text | Google Scholar

Marchand, S., Heylen, K., Messens, W., Coudijzer, K., De Vos, P., Dewettinck, K., et al. (2009). Seasonal influence on heat-resistant proteolytic capacity of Pseudomonas lundensis and Pseudomonas fragi, predominant milk spoilers isolated from Belgian raw milk samples. Environ. Microbiol. 11, 467–482. doi: 10.1111/j.1462-2920.2008.01785.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Mellmann, A., Becker, K., von Eiff, C., Keckevoet, U., Schumann, P., and Harmsen, D. (2006). Sequencing and staphylococci identification. Emerg. Infect. Dis. 12, 333–336. doi: 10.3201/eid1202.050962

PubMed Abstract | CrossRef Full Text | Google Scholar

Metzger, S. A., Hernandez, L. L., Skarlupka, J. H., Suen, G., Walker, T. M., and Ruegg, P. L. (2018). Influence of sampling technique and bedding type on the milk microbiota: results of a pilot study. J. Dairy Sci. 101, 6346–6356. doi: 10.3168/jds.2017-14212

PubMed Abstract | CrossRef Full Text | Google Scholar

Miller, R. A., Kent, D. J., Boor, K. J., Martin, N. H., and Wiedmann, M. (2015). Different management practices are associated with mesophilic and thermophilic spore levels in bulk tank raw milk. J. Dairy Sci. 98, 4338–4351. doi: 10.3168/jds.2015-9406

PubMed Abstract | CrossRef Full Text | Google Scholar

Montel, M. C., Buchin, S., Mallet, A., Delbes-Paus, C., Vuitton, D. A., Desmasures, N., et al. (2014). Traditional cheeses: rich and diverse microbiota with associated benefits. Int. J. Food. Microbiol. 177, 136–154. doi: 10.1016/j.ijfoodmicro.2014.02.019

PubMed Abstract | CrossRef Full Text | Google Scholar

Murphy, M. A., Shariflou, M. R., and Moran, C. (2002). High quality genomic DNA extraction from large milk samples. J. Dairy Res. 69, 645–649. doi: 10.1017/s0022029902005848

PubMed Abstract | CrossRef Full Text | Google Scholar

Murtagh, F., and Legendre, P. (2014). Ward’s hierarchical agglomerative clustering method: which algorithms implement Ward’s criterion? J. Classif. 31, 274–295. doi: 10.1007/s00357-014-9161-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Oberreuter, H., Seiler, H., and Scherer, S. (2002). Identification of coryneform bacteria and related taxa by Fourier-transform infrared (FT-IR) spectroscopy. Int. J. Syst. Evol. Microbiol. 52(Pt 1), 91–100. doi: 10.1099/00207713-52-1-91

PubMed Abstract | CrossRef Full Text | Google Scholar

Oikonomou, G., Bicalho, M. L., Meira, E., Rossi, R. E., Foditsch, C., Machado, V. S., et al. (2014). Microbiota of cow’s milk; distinguishing healthy, sub-clinically and clinically diseased quarters. PLoS One 9:e85904. doi: 10.1371/journal.pone.0085904

PubMed Abstract | CrossRef Full Text | Google Scholar

Panek, M., Cipcic Paljetak, H., Baresic, A., Peric, M., Matijasic, M., Lojkic, I., et al. (2018). Methodology challenges in studying human gut microbiota - effects of collection, storage, DNA extraction and next generation sequencing technologies. Sci. Rep. 8:5143. doi: 10.1038/s41598-018-23296-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Parada, A. E., Needham, D. M., and Fuhrman, J. A. (2016). Every base matters: assessing small subunit rRNA primers for marine microbiomes with mock communities, time series and global field samples. Environ. Microbiol. 18, 1403–1414. doi: 10.1111/1462-2920.13023

PubMed Abstract | CrossRef Full Text | Google Scholar

Quigley, L., McCarthy, R., O’Sullivan, O., Beresford, T. P., Fitzgerald, G. F., Ross, R. P., et al. (2013). The microbial content of raw and pasteurized cow milk as determined by molecular approaches. J. Dairy Sci. 96, 4928–4937. doi: 10.3168/jds.2013-6688

PubMed Abstract | CrossRef Full Text | Google Scholar

Rodrigues, M. X., Lima, S. F., Canniatti-Brazaca, S. G., and Bicalho, R. C. (2017). The microbiome of bulk tank milk: characterization and associations with somatic cell count and bacterial count. J. Dairy Sci. 100, 2536–2552. doi: 10.3168/jds.2016-11540

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmidt, V. S. J., Kaufmann, V., Kulozik, U., Scherer, S., and Wenning, M. (2012). Microbial biodiversity, quality and shelf life of microfiltered and pasteurized extended shelf life (ESL) milk from Germany, Austria and Switzerland. Int. J. Food Microbiol. 154, 1–9. doi: 10.1016/j.ijfoodmicro.2011.12.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Stoeckel, M., Lidolt, M., Achberger, V., Glück, C., Krewinkel, M., Stressler, T., et al. (2016). Growth of Pseudomonas weihenstephanensis, Pseudomonas proteolytica and Pseudomonas sp. in raw milk: Impact of residual heat-stable enzyme activity on stability of UHT milk during shelf-life. Int. Dairy J. 59, 20–28. doi: 10.1016/j.idairyj.2016.02.045

CrossRef Full Text | Google Scholar

Treven, P., Mahnic, A., Rupnik, M., Golob, M., Pirs, T., Matijasic, B. B., et al. (2019). Evaluation of human milk microbiota by 16S rRNA gene Next-Generation Sequencing (NGS) and cultivation/MALDI-TOF mass spectrometry identification. Front. Microbiol. 10:2612. doi: 10.3389/fmicb.2019.02612

PubMed Abstract | CrossRef Full Text | Google Scholar

Vacheyrou, M., Normand, A. C., Guyot, P., Cassagne, C., Piarroux, R., and Bouton, Y. (2011). Cultivable microbial communities in raw cow milk and potential transfers from stables of sixteen French farms. Int. J. Food Microbiol. 146, 253–262. doi: 10.1016/j.ijfoodmicro.2011.02.033

PubMed Abstract | CrossRef Full Text | Google Scholar

Vaidya, J. D., van den Bogert, B., Edwards, J. E., Boekhorst, J., van Gastelen, S., Saccenti, E., et al. (2018). The effect of DNA extraction methods on observed microbial communities from fibrous and liquid rumen fractions of dairy cows. Front. Microbiol. 9:92. doi: 10.3389/fmicb.2018.00092

PubMed Abstract | CrossRef Full Text | Google Scholar

Verdier-Metz, I., Gagne, G., Bornes, S., Monsallier, F., Veisseire, P., Delbes-Paus, C., et al. (2012). Cow teat skin, a potential source of diverse microbial populations for cheese production. Appl. Environ. Microbiol. 78, 326–333. doi: 10.1128/AEM.06229-11

PubMed Abstract | CrossRef Full Text | Google Scholar

Verdier-Metz, I., Michel, V., Delbes, C., and Montel, M. C. (2009). Do milking practices influence the bacterial diversity of raw milk? Food Microbiol. 26, 305–310. doi: 10.1016/

PubMed Abstract | CrossRef Full Text | Google Scholar

von Neubeck, M., Baur, C., Krewinkel, M., Stoeckel, M., Kranz, B., Stressler, T., et al. (2015). Biodiversity of refrigerated raw milk microbiota and their enzymatic spoilage potential. Int. J. Food Microbiol. 211, 57–65. doi: 10.1016/j.ijfoodmicro.2015.07.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, Q., Garrity, G. M., Tiedje, J. M., and Cole, J. R. (2007). Naive Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl. Environ. Microbiol. 73, 5261–5267. doi: 10.1128/AEM.00062-07

PubMed Abstract | CrossRef Full Text | Google Scholar

Wenning, M., Breitenwieser, F., Konrad, R., Huber, I., Busch, U., and Scherer, S. (2014). Identification and differentiation of food-related bacteria: a comparison of FTIR spectroscopy and MALDI-TOF mass spectrometry. J. Microbiol. Methods 103, 44–52. doi: 10.1016/j.mimet.2014.05.011

PubMed Abstract | CrossRef Full Text | Google Scholar

Willis, C., Jorgensen, F., Aird, H., Elviss, N., Fox, A., Jenkins, C., et al. (2018). An assessment of the microbiological quality and safety of raw drinking milk on retail sale in England. J. Appl. Microbiol. 124, 535–546. doi: 10.1111/jam.13660

PubMed Abstract | CrossRef Full Text | Google Scholar

Yoon, S. H., Ha, S. M., Kwon, S., Lim, J., Kim, Y., Seo, H., et al. (2017). Introducing EzBioCloud: a taxonomically united database of 16S rRNA gene sequences and whole-genome assemblies. Int. J. Syst. Evol. Microbiol. 67, 1613–1617. doi: 10.1099/ijsem.0.001755

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: raw milk, biodiversity, microbiota, 16S rRNA gene amplicon sequencing, cultivation

Citation: Breitenwieser F, Doll EV, Clavel T, Scherer S and Wenning M (2020) Complementary Use of Cultivation and High-Throughput Amplicon Sequencing Reveals High Biodiversity Within Raw Milk Microbiota. Front. Microbiol. 11:1557. doi: 10.3389/fmicb.2020.01557

Received: 06 April 2020; Accepted: 16 June 2020;
Published: 09 July 2020.

Edited by:

Abd El-Latif Hesham, Assiut University, Egypt

Reviewed by:

Baltasar Mayo, Consejo Superior de Investigaciones Científicas (CSIC), Spain
Mirna Mrkonjic Fuka, University of Zagreb, Croatia

Copyright © 2020 Breitenwieser, Doll, Clavel, Scherer and Wenning. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Mareike Wenning,

Present address: Franziska Breitenwieser, Tetra Holdings GmbH, Stuttgart, Germany