Bacillus altitudinis GLB197 was examined for its antagonistic potential against fungal pathogens. Results demonstrated that GLB197 significantly inhibited the mycelia growth of nine fungal pathogen with a clear inhibition zone (Figures 1A,B). The endophyte exhibited broad-spectrum antagonistic activities. To understand the cell motility and biofilm formation of the endophytic strain, the swarming motilities and biofilm formation of GLB197 was investigated. The GLB197 exhibited excellent swarming motility and colonized more than half the plate after 5 h of growth (Figures 1C,D). Meantime, GLB197 formed colonies composed of dense wrinkles and compact structures on a MSgg plate and formed robust, wrinkled pellicles on MSgg liquid medium (Figures 1E–G). During the period of 0–48 h, the level of cells accumulated in the biofilm increased about fivefold for GLB197. A previous study has shown that GLB197 could recolonize in grapevine leaves 5 weeks later (Zhang et al., 2017).

In summary, we have isolated B. altitudinis GLB197 from grape leaf, showing broad-spectrum antagonistic activities. Above all, the endophytic strain GLB197 is a beneficial microbe with potential for biotechnological application.

The genome of B. altitudinis GLB197 consists of one chromosome and one plasmid encoding a total number of 3,781 coding DNA sequences (CDS). Among the predicted CDSs, 2,905 of them could be assigned a putative function, whereas 876 were predicted to encode hypothetical proteins. The total length of the chromosome genome was 3,733,835 bp, with a mean G + C content of 41.56%. This chromosome genome contained 3,770 CDS that made up 87.82% of genome. Furthermore, a total of 80 tRNA-coding genes and 24 rRNA genes were predicted in the chromosome sequence. Meanwhile, the GLB197 harbors a plasmid of 7,061 bp. The G + C content of the plasmid is 35.14% and lower than in the chromosome. The genomic features of B. altitudinis GLB197 are shown in Table 1 and Supplementary Figure 1.

The phylogenetic tree of the 17 Bacillus genomes was constructed based on the concatenation of the 1,226 core genes that were present in single copy in all genomes with ML method and rooted by B. cereus ATCC 10987. It displayed that GLB197 and some B. pumilus and B. altitudinis strains formed a big cluster (Figure 2A). Meantime, the phylogenetic analysis revealed that GLB197 was a sister group to W3, which was isolated from raw gallnut honey in China. The sister group of GLB197 and W3 was a clade sister of strains GR-8 and MTCC B6033 (Figure 2A). However, the results showed that they are not distinguished by their geographical origin and the mixed trend of phylogenetic clustering of strains isolated from a similar environment. The results indicated that the phylogenetic clustering of genomes was apparently different from the habitat-specific grouping of these strains.

Furthermore, to confirm the findings from phylogenetic analysis, we also calculated the ANI values of different strains. The results confirmed that there are three major groups of selected B. pumilus group strains (Figure 2B). Strains with ANI values >95% are considered to be the same species (Richter and Rossello-Mora, 2009). As shown in Figure 2B, 17 Bacillus strains clustered three parts. The ANI values of three cluster strains are quite striking. Additional genomic features of GLB197 and 16 B. pumilus group genomes, such as the sequence similarity, was analyzed (Supplementary Figure 2). The analysis indicated a close genetic relatedness of the strains belonging to B. pumilus group and showed that most regions within their genomes were conserved.

A pan-genome for the strain GLB197 and 16 sequenced B. altitudinis strains was determined by the Pan-Genomes Analysis Pipeline (PGAP) software, comparing with the translated CDS set, followed by clustering of orthologous proteins and the representatives of each orthologous cluster and strain-specific CDS in the total pan-genome. The total pan-genome for the 17 compared B. altitudinis strains encompasses 5,506 CDS. Among the 5,506 protein-coding genes, 3,157 genes, which accounted for 57.34% of the genes in the pan-genome of B. altitudinis, were represented in all genomes. The accessory gene families (2,349 genes) were smaller than the core gene families (3,157 genes). The accessory gene families were further classified into 1,150 dispensable genes and 1,199 strain-specific gene families. The smallest numbers of specific genes were encoded by B. altitudinis strains NJ-M2 and NJ-V2 with 21 and 22, respectively. The highest numbers of unique genes (157) was found in B. altitudinis P-10 (Figure 3A).

The genome size for pan-genome and core genome among the selected B. altitudinis genomes was plotted against the number of genomes. The pan-genome curve displayed an asymptotic trend, indicating that 17 genomes were insufficient to describe the complete gene repertoire of the B. altitudinis. The generated pan-genome curves of B. altitudinis are well-represented by the Heaps law mathematical functions: y = 624.45x^0.48 + 3033.86, where y refers to the pan-genome size while x refers to the number of sequenced genomes. According to these equations, the pan-genome size of B. altitudinis appeared to reach infinity when the number of genomes increase to infinity (Figure 3B). Therefore, our data suggest that B. altitudinis have open pan-genomes, which indicates that the species have infinite genomes. Analysis of the core genome was also asymptotic, with 3,157 core genes after the addition of the 17th genome. For these B. altitudinis, we estimated about 72 new genes detected when each additional genome is added from the mathematical equation (Supplementary Figure 3). The infinite pan-genome of B. altitudinis suggests the bacteria will keep acquiring new genes as they evolve independently over evolutionary time.

It is noteworthy that the fermentation broth and cell-free supernatant of B. altitudinis GLB197 showed a distinctive inhibition to grapevine downy mildew on leaf disk and the sporangium release (data not published). Furthermore, the genome sequence of GLB197 was mined for the presence of gene-encoding lipopeptides. Only one gene cluster involved in nonribosomally synthesized lipopeptides was identified in the genome of GLB197 (Supplementary Table 1). The gene cluster probably plays an important role in antimicrobial activity. Figure 4A shows the six genes lchAA, lchAB, lchAC, lgrD, lgrB, and lchAD which are organized as a lichenysin gene cluster arranged within an 8.4-kb genomic region in B. altitudinis GLB197. Meantime, the putative product structure is displayed in Figure 4B. The modular organization of NRPSs involved in biosynthesis of lichenysin in B. altitudinis GLB197 may obey the linear rule (Figure 4C). Similarity searches against parts of publicly available B. altitudinis genomes revealed that the lichenysin gene cluster is also composed of the six genes (Figure 4D and Supplementary Figure 4). The result showed that the gene cluster is highly conserved among B. altitudinis strains. However, the G + C contents of the lichenysin cluster are higher than the genome of B. altitudinis GLB197 (Figure 5). Furthermore, we compared the gene cluster encoding lipopeptides of surfactin family between B. altitudinis and other Bacillus strains. We found that Bacillus strains such as B. subtilis and B. atrophaeus displayed that the gene cluster of surfactin family is composed of four genes (Figure 4D and Supplementary Figure 5). Importantly, the results raise the possibility that the two genes (lgrD and lgrB) were acquired through HGT from other strains. To test this hypothesis, individual phylogentic tree were reconstructed for all of the six cluster genes and their respective homologs identified across different bacterial genomes. The phylogenetic trees based on each of the individual lchAA, lchAB, lchAC, lchAD, lgrD, and lgrB protein sequences (Supplementary Figures 6–11) show that lchABCD of Bacillus are nested with those of Paenibacillus, implying that Paenibacillus and Bacillus have a common lchABCD gene ancestor. Whereas, the phylogeny show that lgrB and lgrD proteins of B. altitudinis fall into the Brevibacillus lineages, supporting that the ancestor of lgrB and lgrD may originate from Brevibacillus. Above all, it is most likely that the lchABCD cluster was transferred from Paenibacillus strains and then B. altitudinis strains acquired the lgrD and lgrB via HGT from Brevibacillus. However, the function and structure need to further confirm.

Endophyte uses multiple strategies to interact with their host. The “PIFAR” open-access web-based tool allowed the identification of several genetic factors involved in plant-bacteria interactions (Manuel Martinez-Garcia et al., 2016). A summary of the genes identified in B. altitudinis GLB197 is provided in Table 2. Specific information about genetic factors can be found in Supplementary Table 2. A set of known mechanisms, which endophytes often use to interact with plants, are found in the genome of B. altitudinis GLB197 such as chemotaxis (cheA, cheW, and cheY), EPs (gpsX and galU), plant adaptation and protection (i.e., four genes for siderophores, one genes for plant hormones, and 11 genes for detoxification), adhesion and biofilm (i.e., attC and six genes for biofilm) (Supplementary Table 2). For instance, the EPS is necessary for rhizoplane and endosphere colonization of rice plants in Gluconacetobacter diazotrophicus (Meneses et al., 2011). The genome analysis revealed genes which might be important for plant-bacteria interaction.

The B. altitudinis GLB197 used in this study, which has been described previously, was isolated from grapevine leaves (Zhang et al., 2017). GLB197 was routinely grown at 37°C on Luria-Bertani (LB) broth or on solid LB medium supplemented with 1.5% agar. For assays of biofilm formation, MSgg medium was used. The recipe for MSgg is as follows: 5 mM potassium phosphate (pH 7.0), 100 mM morpholine propane sulfonic acid (MOPS) (pH 7.0), 2 mM MgCl₂, 700 mM CaCl₂, 50 mM MnCl₂, 50 mM FeCl₃, 1 mM ZnCl₂, 2 mM thiamine, 0.5% glycerol, 0.5% glutamic acid, 50 mg/ml tryptophan, 50 mg/ml threonine, and 50 mg/ml phenylalanine (Branda et al., 2001).

We performed a plate confrontation assay to test the activity of GLB197 against common fungal pathogens including Colletotrichum sp., Alternaria solani, Botrytis cinerea, Valsa mali, and Trichothecium roseum. In the plate confrontation assay, the fungi were cultivated on potato-dextrose-agar plates (PDA) at 28°C for 3–5 days. A 5-mm-diameter block of mycelium culture was cut and transferred into the center of a fresh PDA plates. After 1 day of incubation, 2 μl of the overnight cultures grown in LB medium was spotted on the PDA plate 2.5 cm away from the center, where the mycelium agar block was placed. The antifungal activity was evaluated by measuring diameter of mycelium after 5–7 days of incubation at 28°C. All of the plant pathogens are stored in the lab.

In brief, 5 ml LB liquid cultures was prepared with shaking (200 rpm) at 37°C to OD₆₀₀ = 0.8, 1 ml of cells were collected by centrifugation at 6,000 × g for 5 min, washed with phosphate-buffered saline (PBS, 137 mM NaCl, 2.7 mM KCl, 10 mM Na₂HPO₄, and 2 mM KH₂PO₄), and resuspended in 100 μl PBS. Swarming motility of GLB197 was tested using standard protocols with minor modifications. LB plates containing 0.7% agar were dried in a laminar flow hood for 20 min and then 3 μL of the cell suspension was spotted on the center of each plate (Chen et al., 2012). The plates are incubated at room temperature for 5 h to allow cell growth in order to clearly visualize the swarming zone. Then, the swarming agar plates were dried for another 2 h in a laminar flow hood. The diameter of the swarming zone was measured. The assay was performed with three independent experiments, each with at least three technical replicates.

To monitor pellicle formation, the biofilm formation was analyzed in MSgg medium. GLB197 was grown in LB medium at 37°C overnight. Then, 4 μl of culture was added to 4 ml of MSgg medium in 12-well plates and incubated statically at 28°C for up to 48 h. To quantify biofilm formation, the culture beneath the biofilm was drawn off carefully and the biofilm in each well was washed with 3 ml sterile saline three times and fixed with 2 ml of 99% (v/v) methanol for 15 min, followed by air-drying. The dried biofilms were stained with 2 ml 1% crystal violet (CV) for 10 min. Excess CV was then removed. The CV bound to the cells was dissolved in 5 ml of 33% (v/v) glacial acetic acid and diluted 200 times for OD₅₇₀ detection (Weng et al., 2013). The assay was performed with three independent experiments.

DNA concentration and purity were checked by a spectrophotometer (Nanodrop 2000, Thermo Scientific, Waltham, MA, United States), and the DNA integrity was analyzed by agarose gel electrophoresis. Then, the genomic DNA was sent to Tianjin Biochip Corporation (Tianjin, China), a commercial NGS service provider, for whole genome sequencing. The complete genome of GLB197 was sequenced by the Illumina HiSeq and Pacific Biosciences (PacBio, Menlo Park, CA, United States) platforms. Briefly, genomic DNA sample were fragmented using an ultrasonication approach (Covaris, Woburn, MA, United States) according to the manufacturer’s instruction. The 10-kb template library was constructed using DNA Template Prep Kit 2.0 with the 10-kb insert library protocol. The library was sequenced using the PacBio RSII Sequencer. Meantime, 1 μg DNA sample was also fragmented using an ultrasonication approach, size selected and end repaired. Each generated fragment was ligated to Illumina-specific adapter sequence, quantified, indexed, and sequenced on the Hiseq platform of Illumina. Quality-of-sequence reads were first analyzed using FastQC tool. Then, adaptor sequence removal, trimming, error correction, and assembly were performed using the HGAP software (Chin et al., 2013). Gene predictions were performed with Prokka version 1.11 which predicts coding DNA sequence (CDS) using Prodigal (Seemann, 2014). Annotation of the protein-coding sequence was conducted using the Basic Local Alignment Search Tool (BLAST) against the COG, Kyoto Encyclopedia of Genes and Genomes, and Interpro databases. The final annotated chromosome was plotted using CIRCOS to show the gene locations, GC skew, and GC content (Krzywinski et al., 2009).

We constructed a phylogenetic tree for B. altitudinis GLB197 and other 16 Bacillus genomes, and B. cereus ATCC 10987 was included as an outgroup. All sequences of Bacillus were downloaded from the National Center for Biotechnology information. A phylogenetic tree for the set of 18 genomes was inferred using a core-genome alignment concatenation approach. Multiple alignments of amino acid sequences were conducted using MAFFT (version 7.310) (Katoh and Standley, 2013), and conserved blocks from multiple alignments of test protein were selected by using Gblocks (Castresana, 2000). The maximum likelihood tree was constructed using the RAxML (version 8.2.10) software and the PROTGAMMALGX model with 100 bootstrap replicates. To evaluate the phylogeny of the B. altitudinis strains, ANI was calculated using the JSpecies software with MUMmer (NUCmer) alignment (Richter and Rossello-Mora, 2009). The phylogenetic tree was generated using MEGA (Kumar et al., 2008). Then, the heat map was visualized using R package as a confirmation. Moreover, pairwise genome alignment of all 16 B. altitudinis strains was performed and represented by BRIG 0.95. Single gene alignments were aligned with molecular evolutionary genetics analysis (MEGA). The neighbor-joining trees were constructed using the same software, and 500 bootstraps were done (Kumar et al., 2008).

To identify the core and strain-specific genes, a pan-genome analysis of B. altitudinis strains was carried out by the PGAP software which implements functional ortholog clustering using the amino acid sequences based on Gene Family (GF) method (Zhao et al., 2012). PGAP pipeline-based protein similarity method was used to detect a set of core orthologs from the 17 B. altitudinis strains, and the core orthologs were clustered with at least 50% protein sequence identity to each other and 50% overlap with the longest sequence with an e-value 1e−5. The dataset of shared genes among the 17 strains was defined as their core genomes, and the total set of genes within test genomes was defined as the pan-genome. Strain-specific genes were extracted from the orthologous table by using Perl script.

The core and pan-genomes, as well as estimated respective sizes and trajectories were made using models and regression algorithms proposed by Tettelin et al. (2005, 2008). The curve fitting of the pan-genome was performed using a power-law regression model based on Heaps law [y = Apanx^Bpan + Cpan] as previously described (Tettelin et al., 2005, 2008), where y denotes pan-genome size, x the genome number, and A_pan, B_pan, and C_pan are fitting parameters. Here, B_pan is equivalent to the parameter γ used by Tettelin et al. (2005, 2008) in estimating the open or closed nature of a pan-genome (Rasko et al., 2008). When 0 < B_pan < 1, the size of the pan-genome increases unboundedly with sequential addition of new genomes and can be considered open. Conversely, when B_pan < 0 or > 1, the pan-genome trajectory approaches a plateau as further genomes are added and can be considered closed. The curve fitting of core genome was performed using an exponential regression model [y = A_coree^(B_corex) + C_core] (Tettelin et al., 2005, 2008; Rasko et al., 2008). New gene plots were derived from the pan-genome showing the number of new “strain-specific” genes contributing to the pan-genome per additional sequenced strain as a function of the number of strains. Both the core and pan-genome were visualized through PanGP v1.0.1 software (Zhao et al., 2014). PanGP was run using default parameters generating distribution plots of (i) total genes, (ii) conserved genes, and (iii) new genes found upon progressive sampling of “n” genomes.

Secondary metabolite clustering was predicted using the antiSMASH website (Weber et al., 2015). Moreover, bioinformatics identification of genetic factors (i.e., adhesion, antibiotics, biofilm, detoxification, microbe-associated molecular patterns (MAMPs), multidrug resistance (MDRs), plant cell wall-degrading enzymes (PCWDEs), bacterial lipopolysaccharides (LPSs), siderophores, proteases, etc.) involved in plant-bacteria interaction was performed by implementing the PIFAR open-access, web-based tool (Manuel Martinez-Garcia et al., 2016).

This whole genome shotgun project of B. altitudinis GLB197 has been deposited at DDBJ/ENA/GenBank under the accession CP018574.1.