We used five fungal pathogens, including Fusarium oxysporum, F. solani, F. chlamydosporum, F. tricinctum and Alternaria alternata, that were previously identified as root-rot-causing pathogens in wolfberry plants and they were preserved in glycerol at −80°C (Uwaremwe et al., 2021). The five fungi were cultivated at 25°C for 7 days on potato dextrose agar (PDA) that included 20% potato infusion, 2% dextrose, and 1.5% agar obtained from Qingdao Hope Bio-Technology Co., Ltd., in China. Afterward, a small block of mycelium agar was cut and placed into the center of a fresh PDA plate. One of the BCAs used in the current study was Bacillus amyloliquefaciens strain FBZ42, a commercial strain from the company ABiTEP GmbH that was donated by the Bacillus Genetic Stock Center (BGSC). This strain has shown both plant growth promotion and disease suppression potential for different plants (Kröber et al., 2014). The second bacterial BCA, HSB1, was an uncharacterized endophytic bacterium we isolated from wolfberry root tissues while conducting this study. Briefly, wolfberry roots were cleaned using 75% ethanol for 30 s and immediately transferred to 3% sodium hypochlorite for 5 min, and finally washed three times with sterile water. They were finally cut into small, thin blocks 0.5 cm × 0.5 cm, and placed on NA (Nutrient Agar) culture medium obtained from Qingdao Hope Bio-Technology Co., Ltd., in China, followed by incubation at 37°C. The FZB42 and HSB1 cells were grown and maintained in Luria Bertani liquid medium at 4°C for further experiment.

In the plate confrontation assay, a small block of agar covered with fungal mycelia was excised and placed onto the center of a fresh PDA plate and incubated for three additional days to ensure fungal colonization of the new plate. Afterward, FZB42 and HSB1 were added to Luria Bertani liquid medium (30 mL) and incubated for overnight at 37°C with shaking at 150 rpm. Then, 5 μl of either FZB42 or HSB1 cells were point inoculated on the PDA plate 2.5 cm away from the fungus. Their antifungal activity was evaluated by comparing fungal mycelium growth in the presence of each bacterial BCA, and plates that were inoculated with fungi alone (control plates), after 7 days of incubation at 25°C as reported in Cao et al. (2018).

Strain HSB1 was identified using 16S rRNA gene sequences as previously reported (Mignard and Flandrois, 2006; Janda and Abbott, 2007), after it had already displayed antifungal activity. The resulting 16S rRNA gene sequence was compared in a BLAST search to those in the National Library of Medicine (Bethesda, United States) database (Altschul et al., 1997). Phylogenetic analysis of HSB1 was performed using MEGA 5.5 (Kumar et al., 2008) and the relationships between HSB1 and other Bacillus sequences were analyzed using the neighbor-joining method (Saitou and Nei, 1987; Uwaremwe et al., 2021). Bootstrap values for the neighbor-joining tree were calculated for 1,000 replicates (Felsenstein, 1985; Uwaremwe et al., 2021).

In addition to assessing their antifungal abilities, FZB42 and HSB1 were tested twice for their potential to promote plant growth of Chinese wolfberry seedlings under laboratory conditions. The first assay, conducted from August to October 2019, compared the plant growth promotion potentials of HSB1 and FZB42. To prepare bacterial suspensions, each BCA was first added to Luria Bertani liquid medium (30 mL) and incubated for 12 h at 35°C with shaking at 150 rpm until the logarithmic growth phase was reached. Afterward, all cells were harvested by centrifugation at 5,000 rpm for 5 min, the supernatant was discarded, and the pellets were washed and resuspended in sterile distilled water to obtain initial bacterial population densities of 6 × 10⁹ (FZB42) and 4 × 10⁹ (HSB1) colony forming units (CFU) mL^–1. Finally, sterilized distilled water was used to make the final suspension at 5 × 10⁸ CFU mL^–1 for both HSB1 and FZB42 following the protocol of Cheng et al. (2019). Lycium barbarum “Ningxia N1” seeds obtained from the Institute of Plant Protection, Ningxia Academy of Agricultural and Forestry Sciences, China, were planted in three large plastic pots (40 cm × 60 cm) containing autoclaved pindstrup substrate (pH 5.5–6) obtained from market. Each pot was amended with 2 L of bacterial suspension, with FZB42 as treatment 1, HSB1 as treatment 2, and water control as treatment 3. When each plant had 3–4 leaves, seedlings were individually transferred to small plastic pots (3 seedlings/pot). In each pot was a mixture of autoclaved sand obtained from the Tengger desert in China (latitude 37°30′ to 40° north × longitude 102°20′ to 106° east), soil obtained from Gaolan county in Gansu province (latitude 36°05′–36°51′ north × longitude 103°32′–104°14′ east) and pindstrup substrate in a 1:1:1 ratio (vol/vol/vol), that was irrigated with 50 ml of each treatment suspension (FZB42, HSB1 or water). Root length (cm) and stem height (cm) were recorded for inoculated and non-inoculated seedlings at 60, 68, and 70 days after planting. The second assay was conducted from May through the end of June 2020, and involved HSB1 alone. L. barbarum “Ningxia N1” seeds were planted in a large plastic pot containing autoclaved pindstrup substrate as in the first assay, and the pot was amended with 2 L of water. When each had 3–4 leaves, seedlings were individually transferred to small plastic pots (1 seedling/pot) filled with a mixture of autoclaved sand, soil and pindstrup substrate (1:1:1 ratio, vol/vol/vol), and immediately irrigated with water. When seedlings reached 35 days of growth, HSB1 was inoculated at a concentration of 2 × 10⁷ CFU mL^–1. Shoot weight (g) and root weight (g) were recorded for inoculated and non-inoculated seedlings every 7 days after inoculation. The two experiments were conducted according to a randomized complete design composed of three replicates for each treatment and water control.

In addition to the in vitro biocontrol experiments, and in planta biocontrol experiment was conducted whereby HSB1 and FZB42 were used as BCAs against a strain of F. oxysporum that was previously isolated from Chinese wolfberry and identified based on ITS and TEF (Genbank reference: MN959986 and MT811807) (Uwaremwe et al., 2021). L. Barbarum seedlings were obtained from seed germination. Briefly, seeds were surface sterilized with 75% ethanol for 30 s, sodium hypochlorite for 5 min, and finally washed with sterile distilled water five times. They were planted in two different large plastic pots containing autoclaved pindstrup substrate amended with tap water. All pots were kept in laboratory conditions at a temperature of 25°C, with humidity between 75 and 90%, and an alternating cycle of 16 h light / 8 h dark having a total light intensity of 800 μMol/m2/s. After 20 days, germinated seedlings were individually transplanted into different plug trays containing a mixture of soil, sand and pindstrup substrate in equal ratio (1:1:1 ratio, vol/vol/vol) autoclaved two times to ensure complete disinfection. At the 4–5 leaf developmental stage, individual seedlings were transplanted into their own large plastic pots containing autoclaved soil, sand and pindstrup substrate.

For inoculum preparation, HSB1 and FZB42 were prepared as previously described in the plant growth promotion experiment. To obtain the fungal inoculum, F. oxysporum was cultured on petri dishes containing PDA and incubated at 25°C for 10–15 days. A conidial suspension was prepared by pouring 30 mL of sterile distilled water into each of the petri dishes and dislodging spores with a sterile toothbrush. The initial concentration of conidia in the suspension was determined using a hemocytometer, and the final inoculum concentration was adjusted to 5 × 10⁷ conidia ml^–1. Two methods of inoculation were used for this experiment. The first method involved inoculating seedlings with F. oxysporum 5 days before either FZB42 or HSB1 was introduced. The second inoculation procedure method involved inoculating FZB42 or HSB1 first, followed by a supplementary inoculation after 7 days to ensure their colonization as reported in Gadhave et al. (2018). After another 10 days, the BCA-treated seedlings were inoculated with F. oxysporum. In total, eight treatments were used: (1) CK, untreated seedlings (Control), (2) F. oxysporum alone, (3) FZB42 alone, (4) HSB1 alone, (5) F. oxysporum + FZB42, (6) F. oxysporum + HSB1, (7) FZB42 + F. oxysporum and (8) HSB1 + F. oxysporum. Seedlings were inoculated by pouring 40 ml of the prepared inoculum onto the soil surface. The control treatments consisted of an equivalent volume of sterile distilled water. This experiment was conducted according to a randomized complete design composed of three replicates for each treatment and water control.

Disease incidence (DI) was calculated according to the formula developed by Sharma and Kolte (1994) as the percentage of infected seedlings out of the total of all treated seedlings for each treatment, according to the following formula:

Disease severity (DS) of the foliage was evaluated using a rating scale from 0 to 4 as reported in Han et al. (2019). Based on different stages of root rot, 0 = seedlings with no symptoms, 1 = leaf yellowing, 2 = necrosis, 3 = wilting and 4 = leaf loss. The DS values were obtained from the averages of these scores. Biocontrol efficiency (BE) was calculated according to a formula by Wang et al. (2016) and Cheng et al. (2019) as follows:

Twenty-five days after F. oxysporum inoculation, rhizosphere soils for the three replicate pots were collected from all treatments (24 samples in total) and sieved (2 mm). Briefly, the roots were lightly shaken to remove loosely attached soil. The soil that was still tightly adhering to the roots was harvested as rhizosphere soil and frozen at −80°C for DNA extraction following the protocol of Wu et al. (2016). Total soil DNA was extracted using a Qiagen DNeasy PowerSoil Kit following all steps provided in the kit handbook¹. Genomic DNA concentration and purity were measured using a Qubit fluorometer (Thermo Fisher Scientific, United States).

Bacterial community composition was assessed by sequencing the V1-V9 region of the 16S rRNA gene using PCR primers 27F (5′- AGRGTTTGATYNTGGCTCAG-3′) and 1492R (5′- TASGGHTACCTTGTTASGACTT-3′) as reported in Větrovsky and Baldrian (2013). PCR conditions were initiated at 95°C for 5 min, followed by 25 cycles of denaturation at 95°C for 30 s, annealing at 50°C for 30 s, and extension at 72°C for 1 min, followed by a final elongation at 72°C for 7 min, and then hold at 4°C. The PCR products were pooled and visualized on 1% agarose gels, purified using a MinElute PCR Purification Kit according to the manufacturer’s instructions, and quantified using QuantiFluorTM-ST (Promega, United States). High-throughput sequencing was carried out on the PacBio Sequel II platform (BioMarker Technologies Co., Ltd., China). The original subreads were first corrected to generate circular consensus sequences (CCS) (SMRT Link, version 8.0), and then Lima software (v1.7.0) was used to identify the CCS of different samples through their barcoded sequences, and UCHIME1 (version 8.1) was used to remove the chimera bodies for high quality CCS sequences. Using USEARCH 4 (version 10.22) with a cut-off of 97% similarity, the Operational Taxonomic Units (OTUs) were clustered and the taxonomic classifications were performed using RDP Classifier (Version 2.2, based on Bergey’s taxonomy) with the classification threshold set at 0.5.

Putative bacterial metagenomic functions were inferred using a phylogenetic investigation of communities by reconstruction with unobserved states (PICRUSt) on the 16S rRNA gene abundance data as reported by Langille et al. (2013). Using functions within the PICRUSt pipeline, the OTU-table was normalized and used for metagenome inferences involving the KEGG (Kyoto Encyclopedia of Genes and Genomes) orthologs (KOs). The predicted functions were then collapsed into hierarchical KEGG pathways using the “categorize by function” step in the PICRUSt pipeline as performed by Wilkinson et al. (2017).

Plant growth promotion data were analyzed using the analysis of variance (ANOVA) procedure of SAS 8.1 software (SAS Institute Inc., Cary, NC, United States). Differences between treatments were assessed at each time point by Fisher’s protected least significance difference (LSD) test at 0.05 levels. Kruskal–Wallis one-way analysis of variance by ranks was used for comparing DI and DS. The rarefaction curve, corresponding to observed OTUs at different sequencing depths, was examined using QIIME software to determine whether the depth was reasonable. Chao1 and abundance-based coverage estimator (ACE) indices were used to calculate the evenness of each sample (Bokulich et al., 2013), while the Shannon and Simpson indices were used to measure diversity (Hong et al., 2015). Beta diversity among samples was determined by principal component analysis (PCA) using R software². Significant differences in bacterial community composition between paired samples were determined using the Metastats analysis and Mothur program which counted taxa in five classified levels (Lu et al., 2016). The BE values were not subjected to statistical analysis.

16S rRNA was used to identify the genus and species of HSB1. This sequence has been accessioned to GenBank (MT626060). The BLAST result showed 99.93% homology with a strain of Bacillus amyloliquefaciens. The phylogenetic analysis showed HSB1 clustering with other B. amyloliquefaciens, B. methylotrophicus, and B. velezensis strains, which are known to share identity and a most recent common ancestor (Figure 1a). B. mojavensis strains were used as out-group (Figure 1b).

The potential of the two B. amyloliquefaciens strains, HSB1 and FBZ42, to inhibit the five root rot fungal pathogens (F. oxysporum, F. solani, F. chlamydosporum, F. tricinctum and A. alternata) was assessed using dual culture technique. The results showed that both HSB1 and FZB42 inhibited mycelial growth of all five fungal pathogens compared to the control (petri plates without bacterial served as control (Figure 2). All five fungal pathogens were inhibited (up to 100%). Due to its faster growth rate compared to the Fusarium strains tested, A. alternata was the first to show inhibition from exposure to our BCAs (Figures 2A–E).

This set of experiments showed that the two bacterial strains, FZB42 and HSB1, had different effects on the development of wolfberry seedlings. Both strains promoted plant growth in some way compared to the CK. In the first assay comparing HSB1 to FZB42, all seedlings were growing at same rate during first 60 days (Supplementary Figures 1A–C). In the days after, both FZB42 and HSB1 increased stem length compared to the water control. At 68 and 78 day time points, the greatest stem length was observed in HSB1-treated seedlings (Figure 3A). Measurements showed that both FZB42 and HSB1 increased root length in 60-day-old seedlings compared to the control. However, HSB1 had significantly increased the root length compared to both the FZB42 and control treatments, with FZB42-treated roots having shorter root length than even the water control treatment (Figure 3B). In the second assay, HSB1-treated seedlings were healthier and taller compared to control treatment (Supplementary Figures 1D,E). Results showed that HSB1 increased wolfberry seedlings shoot weight at the 7, 14, and 21 days post-inoculation time points compared to the control. There was a significant difference between HSB1 and control shoot weights (Figures 4A,B). HSB1also increased seedling root weight compared to control treatment. There was a significant difference between HSB1 and control at 14 and 21 day time after inoculation (Figures 4C,D).

As a baseline, DI and DS were assessed at 15, 21, and 25 day post-inoculation with F. oxysporum alone. The highest DI (42%) and DS (1.7) were recorded for day 15 (Table 1 and Figure 5). In the foliage, the disease usually was recognized by the yellowing, necrosis and wilting of leaves, followed by complete leaf loss (Table 1 and Supplementary Figure 2). At the beginning, infected seedlings were characterized by severe yellowing of leaves and the disease symptoms continually increased throughout the experiment (Supplementary Figures 2A,B and Table 1). Assessment conducted at 21 and 25 day after F. oxysporum inoculation revealed that the greatest respective DI (65 and 84%) and DS (1.28 and 1.45) values were also recorded in F. oxysporum alone (Figures 5B,C,E,F). Unexpectedly, no DI and DS were recorded in treatments where F. oxysporum was inoculated before HSB1 (FO + HSB1) or where FZB42 was inoculated before F. oxysporum (FZB42 + FO) treatments throughout the experiment (BE = 100%); all seedlings had zero disease symptoms in the foliage (Supplementary Figures 2G,H). Additionally, control plants treated with water (CK), FZB42 alone, as well as HSB1 alone did not also show any disease symptoms (Supplementary Figures 2E,F,I). Brown lesions were observed on root surfaces of seedlings treated with FO alone or HSB1 + FO (Supplementary Figures 3A,C), whereas no lesions were observed on root surfaces of seedlings having the other treatments (Supplementary Figures 3B,D–H).

A total of 211,372 16S rRNA V1-V9 gene sequences were analyzed across 24 rhizosphere soil samples, with an average of 8,807 ± 959 sequences per soil sample (Supplementary Table 1). Based on a threshold of 97% shared nucleotide identity, these sequences were grouped by OTU, yielding 8,366, 8,259, 7,901, 7,484, 7,145, and 6,623 OTUs at the phylum, class, order, family, genus, and species levels, respectively (Table 2). The sequencing depth was analyzed to identify new taxa. A rarefaction curve analysis at 3% dissimilarity for the bacterial community revealed that the sharp of the curve was increasing and depth did not reach saturation, indicating that a greater sequencing depth was needed (Figure 6A). However, the data were sufficient for showing differences among the treatments and suggested that BCA application increased bacterial diversity. The richness indices (ACE and Chao1), and diversity indices (Shannon and Simpson) were further estimated and presented in Table 3. Chao1 and ACE indices were higher in HSB1 + FO and FZB42 + FO treatments compared to FO alone. The lowest values were recorded in HSB1 alone and FO + HSB1 treatments. In addition, Simpson index was higher in HSB1 alone and FO + HSB1 treatments compared to the FO alone treatment and CK, whereas the Shannon index was higher for FO + FZB42, FZB42 + FO, and HSB1 + FO treatments. PCA comparison of changes in soil bacterial community across different treatment groups showed that the first two principle components could explain 42.82 and 17.92% of the total variation. The bacterial communities differed according to treatment, indicating the effects of F. oxysporum and the two Bacillus strains on bacterial composition in the wolfberry rhizosphere (Figure 6B).

The average number of microbial groups at phylum, class, order, family, genus and species levels were 24, 34, 73, 105, 192, and 239, respectively (Supplementary Table 2). All samples showed similar phylum and genus composition, but differed in terms of the relative abundances of other taxonomic groups (Figure 7). Of the 10 most abundant phyla across all samples, Proteobacteria and Bdellovibrionota were the least abundant. Although phylum Proteobacteria predominated, each treatment affected the relative abundance of this phylum in the overall rhizosphere bacterial composition (Figure 7A). The 10 most abundant genera included Massilia followed by Arenimonas, Pelomonas, Gemmatimonas, Vicinamibacter, Comamonas, Pseudoxanthomonas, Pseudomonas, Pedosphaera and Piscinibacter (Figure 7B). Consistent with our phylum level observations, relative abundances of Massilia differed by treatment (Figure 6B and Supplementary Table 3). Finally, LEfSe analysis was used to detect taxa with significantly different abundances between the FO alone, FO + HSB1, FO + FZB42, HSB1 alone, and FZB42 + FO treatments. The most differentially abundant bacterial taxa in the rhizosphere samples belonged to the Proteobacteria phylum (Supplementary Figure 4).

KEGG ortholog predictions were performed on the 16S rRNA data using PICRUSt. We conducted comparisons between different treatments, including CK vs. HSB1, F. oxysporum vs. HSB1, and HSB1 vs. F. oxysporum + FZB42 (Figure 7). Seven pathways related to lipid transport and metabolism, transcription, energy production and conversion, amino acid transport and metabolism, inorganic ion transport and metabolism, secondary metabolite biosynthesis, transport and catabolism, and “function unknown” were overrepresented in the HSB1 alone sample compared to the CK (Figure 8A). Comparison between F. oxysporum and HSB1 showed that pathways related to membrane transport, cellular community-prokaryotes, and lipid metabolism were overrepresented in the HSB1 alone sample (Figure 8B). In comparing HSB1 to FO + FZB42, pathways related to metabolism of cofactors and vitamins, translation, glycan biosynthesis and metabolism, replication and repair, global and overview maps were overrepresented in the FO + FZB42 sample (Figure 8C).