A total of 90 tea samples were collected at different processing stages from four companies in Anhua, Hunan Province, and one processing facility in Meitan, Guizhou Province, China. The sampled fermentation stages were (in chronological order) as follows: fresh tea (9 samples), pile processing (14 samples), drying (14 samples), early stage golden flower bloom (2–12 h) (22 samples), full golden flower bloom (1–20 days) (23 samples), and post-fermentation (8 samples; Figure 1). The first three stages cover the processing of fresh tea leaves to raw dark tea, and the last three stages cover the processing of raw dark tea to Fu Bricks (Zheng et al., 2015). The pile fermentation and golden flower stages were each sampled at different time points, within each stage, to more accurately assess changes to the microbial community, which were averaged to obtain an overall effect during each stage. The samples were placed in sterile polyethylene bags and transported on ice to the laboratory. All the samples were portioned immediately, with one stored portion at 4°C for culturing of microorganisms and the other stored at −80°C for DNA extraction.

For each sample, 2.5 g of DNA was extracted with the MoBio PowerSoil^® DNA Isolation Kit (MO BIO Laboratories, Inc.) by following the manufacturer's protocols. The extracted DNA was stored at −20°C for further PCR amplification and barcoded amplicon sequencing. Primer sets ITS1-F (5'-CTTGGTCATTTAGAGGAAGTAA-3') and ITS2-R (5'-GCTGCGTTC TTCATCGATGC-3') were used for the amplification of fungal DNA, and primer sets V4-F (5'-GTGCCAGCMGCCGCGGTAA-3') and V4-R (5'-GGACTACHVGGGTWTCTAAT-3') were used for the amplification of bacterial DNA. Both the forward and reverse primers had 6-bp barcodes unique to each sample, which were used to permit multiplexing of samples. The DNA concentrations were adjusted to 20 ng/μl for each sample in the PCR. Each 50 μl PCR solution consisted of 25 μl 2x Taq PCR MasterMix (TIANGEN Company, Beijing), 1.5 μl sense primer (10 μM), 1.5 μl anti-sense primer (10 μM), 3 μl DNA template, and 19 μl ddH₂O. PCR amplifications were performed in triplicate reactions on a Veriti thermal cycler (Applied Biosystems). PCR amplifications were performed under the following conditions: initial denaturing at 94°C for 4 min, followed by 35 cycles of 94°C for 40 s, 54°C for 45 s (for ITS), or 52°C for 45 s (for V4), 72°C for 50 s, and a final extension step at 72°C for 10 min. Negative controls with no template DNA were also included in PCR amplifications. The products of the triplicate PCRs of each soil sample were pooled and subjected to gel electrophoresis, from which purified amplicons were harvested using the EasyPure Quick Gel Extraction Kit (TransGen, Beijing). The purified amplicons and products encoding different barcodes were combined into a single library for automated cluster generation and paired-end sequencing on the MiSeq system.

Primary data analysis (image analysis, base calling) was performed on a MiSeq instrument. The sequences were analyzed using the Quantitative Insights Into Microbial Ecology (QIIME) toolkit with default parameters for each step (Caporaso et al., 2010). Fungal and bacterial sequences with lengths shorter than 220 nucleotides were removed. De novo and reference-based chimera checking were performed, and sequences that were characterized as chimeric were also removed. Bacterial reads were binned into OTUs at the ≥97% sequence similarity level using the open-reference OTU picking protocol in UPARSE-pipeline (Edgar, 2013), and the most abundant sequences from each OTU were taken as representative sequences for the respective OTU. Taxonomic configuration of bacterial OTUs was performed by searching the Basic Local Alignment Search Tool (BLAST) for each representative sequence against a subset of the Silva database (Quast et al., 2013). The fungal sequences were subjected to ≥97% sequence similarity single-linkage clustering using GenBank (Benson et al., 2005) and UNITE (Abarenkov et al., 2012) as reference taxonomic databases.

For each sample, 10 g was homogenized with sterile distilled water using a blender at 12,000 rpm using three 30 s pulses. The homogenate was serially diluted in sterile distilled water and spread onto each of four types of media: MEA (2% malt extract, 0.1% peptone, 2% glucose, 1.8% agar; containing 100 mg/L chloramphenicol and 100 mg/L streptomycin), LA medium (BD Co.), YDP (1% yeast extract, 1% peptone, 2% glucose, 1.8% agar; containing 100 mg/L chloramphenicol), and Gauze No.1 (AOBOX BIOTECHNOLOGY; 25 mg/L nalidixic acid and 50 mg/L cycloheximide) (Liu et al., 1996; Rutiaga-Quiñones et al., 2012). Plates were prepared for each sample and medium in triplicate. The plates were incubated at 25°C for 2 days, and colonies were transferred to new plates to purify and preserve.

For the identification of microbes from pure cultures, morphological examinations were conducted with microscopy. These identifications were then followed by the implementation of DNA-based methods. The genomic DNA was extracted with a simple and rapid “thermolysis” method (Zhang et al., 2010). Briefly, a sterile toothpick was used to transfer mycelia or yeast cells from the colony into 100 μl pure water. The mixture was vortexed thoroughly and centrifuged at 10,000 g for 1 min. The supernatant was removed, and then, 100 μl of lysis solution [EasySelectTM Pichia Expression Kit (Invitrogen, USA)] was added into the microcentrifuge tube. The mixture was then incubated at 85°C for 20–30 min. The extract contained genomic DNA and was stored at −20°C for future use. Multiple primers are used for confirmation of our morphological identifications; they are listed in Supplementary Table 1, and they included the ITS region of fungi (White et al., 1990), the D1/D2 region of yeasts (Boekhout et al., 1995), the 16S rRNA gene of bacteria with degenerate bases (DeLong et al., 2006), and the 16S rRNA gene of actinomycetes without degenerate bases (Lane, 1991). Special primers were used for better differentiating yeast and actinomycetes species based on literature. Each 25 μl PCR solution consisted of 12.5 μl 2x Taq PCR MasterMix (TIANGEN Company, Beijing), 1 μl sense primer (10 μM), 1 μl anti-sense primer (10 μM), 2 μl DNA template, and 8.5 μl ddH₂O. PCR amplifications were performed on a Veriti Thermal Cycler (Applied Biosystems) under the following conditions: initial denaturing at 94°C for 3 min, followed by 35 cycles at 94°C for 30 s, 55°C for 30 s, and 72°C for 1 min, with a final extension step at 72°C for 10 min (fungi); initial denaturing at 94°C for 3 min, followed by 35 cycles at 94°C for 40 s, 55°C for 30 s, and 72°C for 80 s, with a final extension step at 72°C for 10 min (bacteria); initial denaturing at 94°C for 3 min, followed by 35 cycles at 94°C for 30 s, 55°C for 30 s and 72°C for 1 min, with a final extension step at 72°C for 10 min (yeast); initial denaturing at 94°C for 3 min, followed by 35 cycles at 94°C for 30 s, 55°C for 30 s, and 72°C for 90 s, with a final extension step at 72°C for 10 min (actinomycetes). The PCR products were sent to Omega Genetics Co. Ltd., China, for downstream sequencing, and the obtained sequences were uploaded to GenBank database to search for homologous sequences by BLAST.

Sequencing data were normalized using total sum scaling (Weiss et al., 2017). Culture-dependent data were shown as counts (CFU/g). All data were checked for homogeneity of variance, and normality was confirmed by inspection of residuals. The overall effects of manufacturing stages on the composition of bacterial and fungal communities were examined using permutational multivariate ANOVA (PERMANOVA) based on a Bray-Curtis dissimilarity matrix (Oksanen et al., 2006). Non-metric multidimensional scaling (NMDS) plots, also based on Bray-Curtis, were used to visualize the separations between manufacturing stages. Linear mixed models were used to examine the effects of the manufacturing process on the Simpson diversities of bacteria and fungi, and the relative abundances of bacterial phyla, fungal classes, A. cristatus, yeast, and lactic acid bacteria. In the model, processing stages were used as the fixed factors, and different sampling facilities were used as random factors. Then, post-hoc Tukey's tests were used for the pairwise comparisons between all processing stages. The linear mixed model and post-hoc Tukey's test described earlier were also used to examine the effects of processing stages on the relative abundances of core microbes (those that dominated in 80% of all samples or had relative abundances >1%). Network analysis based on Pearson's correlations (P < 0.05) was used to explore only the relationships between core microorganisms because too few sample sizes would lead to false-positive correlations. The network was constructed using Gephi-0.9.2. Given the likelihood that not all microbes present in our samples would be culturable during the culture-dependent experiments and would result in missing values in the culture-dependent dataset, statistical analysis may not be possible for the culture-dependent dataset. All analyses were performed in R (version 3.0.1, R Development Core Team, 2013).

Fermentation stage one explained 13% of the variances in the bacterial community and 17% in the fungal community (Table 1). Fermentation stage two did not significantly influence the microbial communities. The overall effects of fermentation stages explained 16% of the variances in the bacterial community and 15% in the fungal community (Table 1 and Figure 2). Fresh tea leaves contained more distinct fungal and bacterial communities compared with those present during pile fermentation (Figures 2A,B). Manufacturing stages significantly influenced the Simpson diversity of bacterial and fungal communities (Figures 2C,D). With regard to bacterial diversity, pile fermentation had the lowest level, while early stage golden flower bloom had the highest (Figure 2C). For fungi, the drying stage had the highest amount of diversity, while the Fu Brick tea product (i.e., the final stage of processing) had the lowest (Figure 2D).

Across all stages of the manufacturing process, the dominant phyla comprising bacterial communities included the Proteobacteria, Firmicutes, Bacteroidetes, Actinobacteria, and Acidobacteria (Figure 3A and Supplementary Table 2). Fresh tea leaves had the highest relative abundances of Proteobacteria and Acidobacteria and the lowest relative abundance of Firmicutes (Supplementary Table 2). Throughout both fermentation stages, however, the relative abundances of Proteobacteria and Acidobacteria decreased, while the relative abundance of Firmicutes increased (Supplementary Table 2). The final Fu Brick product contained the highest relative abundance of Firmicutes and the lowest relative abundances of Proteobacteria and Acidobacteria (Supplementary Table 2).

At the genus level, the bacterial communities on fresh tea leaves were dominated by Klebsiella, followed by Methylobacterium. The pile fermentation and drying stages increased the relative abundances of Klebsiella, Sphingomonas, Paenibacillus, and Streptomyces, which then decreased at the golden flower bloom stages (Table 2). Throughout the fermentation stages, the relative abundance of Methylobacterium decreased significantly and reached its lowest level in the Fu Brick product (Table 2). Conversely, the abundances of Bacillus and Pelomonas species were lowest in fresh tea leaves, but they significantly increased during pile fermentation and reached their highest levels in the Fu Brick product (Table 2).

The dominant classes within fungal communities throughout the manufacturing process consisted of Eurotiomycetes, Sordariomycetes, Saccharomycetes, and Dothideomycetes (Figure 3B and Supplementary Table 2). Fresh tea leaves contained the lowest relative abundance of Eurotiomycetes and the highest relative abundances of Dothideomycetes and Sordariomycetes (Supplementary Table 3). As subsequent processing stages progressed, the relative abundance of Eurotiomycetes increased and the relative abundances of Dothideomycetes and Sordariomycetes decreased, with Fu Bricks containing the highest relative abundance of Eurotiomycetes and the lowest relative abundances of Dothideomycetes and Sordariomycetes (Supplementary Table 3). The relative abundance of Saccharomycetes first increased, and became highest in the pile fermentation stage, and then decreased, with Fu Bricks containing the lowest relative abundance of Saccharomycetes (Supplementary Table 3).

So far, Aspergillus, especially A. cristatus, was the dominant fungal genus, which increased significantly across all processing stages from fresh tea leaves (low levels) to comprise 98% of the Fu Brick product (Table 2). Other fungal genera predominated in fresh tea leaves, such as Penicillium, Purpureocillium, and Cladosporium. As fermentation progressed, their relative abundances fluctuated. By the drying stage, Penicillium and Purpureocillium species had the lowest relative abundances (Table 2). Both fermentation stages caused consistent decreases in the relative abundance of Cladosporium fungi with the Fu Brick product containing the lowest relative abundance of this genus (Table 2). Despite being considered a dominant fungal class, no genus within the Saccharomycetes satisfied the two standards of core microbes for this study.

In this study, we classified the microbes as A. cristatus, yeasts, or lactic acid bacteria based on their colony forming units (cfu) (Figure 4). Overall, we found that the culture-based data at each processing stage supported the sequencing data: fresh tea leaves contained the lowest abundance of A. cristatus cfu, which increased along the manufacturing process to attain the greatest abundance in the Fu Brick product (Figures 4A,E). Yeasts were abundant in fresh tea leaves, but during fermentation, their numbers decreased. This fate continued through the final Fu Brick stage where they had the lowest abundance (Figures 4B,F). The culture-dependent data showed that all isolated bacteria cfu (including the lactic acid bacteria) were highest at the drying stage and lowest in the Fu Brick product (Figure 4G), while the sequencing data suggested that the relative abundances of lactic acid bacteria did not significantly fluctuate during the manufacturing process (Figure 4C).

Network analysis was conducted for the microbes that presented over 80% of all samples (Figure 5). Klebsiella, Staphylococcus, and Acinetobacter were in the core of the network, Klebsiella negatively correlated with most other microbes, while Staphylococcus and Acinetobacter had both positive and negative relationships with other microbes. Bacillus, Penicillium, and Streptomyces positively correlated with A. cristatus, but Lactobacillus, Klebsiella, and Cyberlindnera negatively correlated with this fungus.