Analysis of fungal dynamic changes in the natural fermentation broth of ‘Hongyang’ kiwifruit

Materials

Fresh ‘Hongyang’ kiwifruits (see Fig. 1) were collected from Aziying, Panlong District, Kunming, Yunnan Province, China. Sterilized medical rubber gloves were used to pick fresh fruits to avoid the contamination of the samples. The average single fruit of ‘Hongyang’ kiwifruit weighs 75 g, with red pulp scattered around the central column of the cross-section.

Yeast in naturally fermented ‘Hongyang’ kiwifruit juice

Five hundred grams of ripe, fresh, and unsterilized fruits were chopped with the sterilized knife on the ultra-clean workbench. They were evenly placed in three sterile one-liter Erlenmeyer flasks, added 167 mL ddH₂O, and fermented at 28 °C. Two milliliters of natural fermentation broth samples on the first day (D1T), the third day (D3T), the fifth day (D5T), and the fifteenth day (D15T), were taken from each Erlenmeyer flask respectively and were placed in a centrifuge tube. Samples were centrifuged in a centrifuge (6000 rpm, 10 min) and washed with phosphate buffer saline (PBS) buffer (1/15 mol/L KH₂PO₄, and 1/15 mol/L Na₂HPO₄. 2H₂O, pH 7.4) three times. Then cell pellets obtained were stored in liquid nitrogen for high-throughput sequencing.

Yeast in enrichment culture

After taking and chopping 500g of ‘Hongyang’ kiwifruit as in the previous paragraph, the fruits chopped were evenly divided into three one-liter Erlenmeyer flasks containing 167 mL YPD liquid medium (yeast extract 10 g/L, peptone 20 g/L, glucose 20 g/L, agar 20 g/L, pH 7.0, 121 °C, sterilization for 15 minutes, stored at 4 °C for later use). Fungi were enriched and cultivated at 28 °C until the liquid medium became turbid, then sampled, and named CK1, CK2, and CK3. Three replicate samples were centrifuged in a centrifuge (6000 rpm, 10 min) and washed with PBS buffer three times, and stored in liquid nitrogen for high-throughput sequencing.

High-throughput sequencing and bioinformatics analysis

Results and analysis

Species composition

A total of 1,308,736 sequences were obtained, and the average length of the sequences was 219 bp. According to the 95% similarity, the samples obtained 582 OTUs. In D1T, D3T, D5T, D15T, and CK samples, a total of 131 genera and 178 species were found (See Annex 1). Among the samples from natural fermentation, D1T had 112 genera and 146 species, D3T had 36 genera and 42 species, D5T had 38 genera and 47 species, D15T had 31 genera and 37 species. In the enriched culture samples, CK had 34 genera and 42 species (See Annex 1). As the fermentation of kiwifruit continued, the number of yeast species continued to decrease, and alpha diversity metrics also showed a decreased trend in diversity (see Table 1). In addition, the number of non-Saccharomyces species obtained through YPD medium enrichment culture (CK sample) was lower than those of D1T, D3T, and D5T (p-Value = 0.000942). See Annex 1 for the relationships between species composition of different samples. The Krona analysis results of all samples are shown in Fig. 2. The yeasts accounted for 80% of all the species. Among the yeasts, the Pichia genus accounted for the most significant proportion (51%), and it was the dominant genus in the samples (abundance value = 0.409795, p-Value = 0.000201). The top ten species in terms of richness for each sample were plotted in a bar chart (see Fig. 3). The dominant fungi in D1T were A. pullulans (abundance value = 0.117850, p-Value = 0.002684). However, as the fermentation continued, the content of A. pullulans and other fungi continued to decrease. On the 5th day of fermentation, A. pullulans disappeared. The dominant fungi were P. terricola (abundance value = 0.429040, p-Value = 0.000021) and P. kluyveri (abundance value = 0.392138, p-Value = 0.000000). Subsequently, the advantages of P. terricola and P. kluyveri continued to decrease. On the 15th day, the proportions of P. terricola (abundance value = 0.176059, p-Value = 0.000021) and P. kluyveri (abundance value = 0.176925, p-Value = 0.000000) dropped respectively. On the contrary, as the fermentation continued, P. membranifaciens (abundance value = 0.573930, p-Value = 0.001871) gradually became the dominant species. The amount of P. kluyveri showed an increasing trend, followed by a decrease. In the CK, the dominant species were P. kluyveri (abundance value = 0.176925, p-Value = 0.000000). The highest abundance of P. kluyveri was found in the CK, enriched in culture with YPD.

A heat map was drawn for the abundance data of the top 10 species in average abundance to reflect the correlation of fungi among samples and display the distribution trend of species abundance in each sample. The results are shown in Fig. 3. The relationship between each sample and each fungus could be seen from the heat map. The species diversity was the highest on the first day. As the natural fermentation progressed, the diversity gradually decreased. The results of the heat map analysis (see Fig. 4) were consistent with the species composition analysis.

Alpha diversity analysis

The Alpha diversity index analysis results among the samples are shown in Table 1 and Fig. 5. The Goods coverage values of D1T, D3T, D5T, D15T, and CK in the sample were all about 1 (see Table 1), indicating that the sequencing depth of the sequencing results of all samples this time had covered all species, and the sequencing results represented the yeast in the sample.

The results of the Alpha diversity index analysis revealed that in the natural fermentation process of kiwifruit, the Shannon index of D1T was the highest (2.424493, p-Value = 0.000000), and Simpson was relatively low (0.601790, p-Value = 0.000003), reflecting the high species diversity of yeast in the sample D1T. The Shannon index value of D15T was the lowest (1.948073, p-Value = 0.000000), and the Simpson index value was the highest (0.581287, p-Value = 0.000003), reflecting the lowest species diversity in D15T. In the four samples of naturally fermented juice D1T, D3T, D5T, and D15T, the Shannon index first decreased and then increased, and the Simpson index first increased and then decreased, indicating that as the natural fermentation continued, the yeast species diversity first increased and then decreased. In the enriched culture broth, the Shannon index and Simpson index were smaller than those in the natural fermentation process (Shannon index = 0.148268, Simpson index = 0.029788, and p-Value = 0.000000), indicating that the yeast diversity in the kiwifruit enriched culture broth was lower than that in the natural fermentation process.

In the natural fermentation process, the Chao1 index and observed species index of D1T were the highest (Chao1 index = 147.559333, observed species index = 146.600000, and p-Value = 0.000001). Among D3T, D5T, and D15T, although the Chao1 index and observed species index of D3T was the lowest (Chao1 index = 32.572333, observed species index = 31.366667, and p-Value = 0.000000), indices of D5T was the highest (Chao1 index = 37.960167, observed species index = 36.333333, and p-Value = 0.000000), and the overall trend was increasing and then decreasing, but compared to D1T stability means that in the natural fermentation process, the number of species would drop sharply at the beginning of the fermentation and then stabilize. Compared with natural fermentation broth, the Chao1 index and observed species index in the enriched culture broth were lower than those in the natural fermentation broth (Chao1 index = 32.357000, observed species index = 30.366667, and p-Value = 0.000000), indicating that the number of yeast species in the enriched culture broth of kiwifruit was lower than the number of yeast species in the natural fermentation broth.

Beta diversity analysis

The Beta diversity analysis results are shown in Fig. 6. The species composition of the CK samples was far away from the species composition of the natural fermentation samples, indicating that there were not many species in the CK samples, and the colony structure of CK samples was quite different from that of the natural fermentation ones.

The similarities among samples were displayed in the form of a hierarchical tree. As shown in Fig. 7, the species composition distance between D3T and D5T samples was the closest (p-Value = 0.000147), indicating that the species composition between the two samples was the most similar, and P. kluyveri occupied a relatively large proportion (D3T abundance value = 0.363838, D5T abundance value = 0.392138, p-Value = 0.000000). On the other hand, the distance between the species composition of the D1T sample and other samples was the furthest (p-Value = 0.000080), indicating that the species composition of D1T was different from other samples. The difference between them was enormous. This result was consistent with the results of PCoA and NMDS analysis.

Functional potential prediction

KEGG database divided natural fermentation samples’ metabolic pathways into five categories: biosynthesis, degradation/utilization/assimilation, generation of metabolite and energy, glycan pathways, and metabolic clusters. By comparing the abundance of five kinds of samples, the results are shown in Fig. 8. The high proportion of biosynthesis and precursor metabolite and energy generation indicated that the natural fermentation process was active and produced many by-products. Furthermore, the fatty acid and lipid degradation pathways in the generation of precursor metabolite and energy and degradation/utilization/assimilation were directly related to the production of aromatic ester compounds during fermentation.

To explore the reasons for the changes of metabolic pathways, MetaCyc metabolic pathway was also predicted in this study, took the enriched cultured sample CK as the control group, and analyzed the differences of metabolic pathways in each sample of natural fermentation. Then the metabolic differential pathways following the changes of metabolic pathways in each sample were compared, to clarify which pathways affected the changes of KEGG pathways. The results are shown in Fig. 9.

There is a direct relationship between metabolic pathways and colony structure. It was found that the pathway changes were regulated by Pichia, unidentified, Alternaria, and unclassified Ascomycota; Pichia affected most of the pathways (see Figs. S1–S6).