Variations of Phyllosphere and Rhizosphere Microbial Communities of Pinus koraiensis Infected by Bursaphelenchus xylophilus

Pine wood nematode, Bursaphelenchus xylophilus, as one of the greatest threats to pine trees, is spreading all over the world. Plant microorganisms play an important role in the pathogenesis of nematodes. The phyllosphere and rhizosphere bacterial and fungal communities associated with healthy Pinus koraiensis (PKa) and P. koraiensis infected by B. xylophilus at the early (PKb) and last (PKc) stages were analyzed. Our results demonstrated that pine wood nematode (PWD) could increase the phyllosphere bacterial Pielou_e, Shannon, and Simpson index; phyllosphere fungal Chao 1 index, as well as rhizosphere bacterial Pielou_e, Shannon, and Simpson index; and rhizosphere fungal Pielou_e, Shannon, and Simpson index. What’s more, slight shifts of the microbial diversity were observed at the early stage of infection, and the microbial diversity increased significantly as the symptoms of infection worsened. With the infection of B. xylophilus in P. koraiensis, Bradyrhizobium (rhizosphere bacteria), Massilia (phyllosphere bacteria), and Phaeosphaeriaceae (phyllosphere fungi) were the major contributors to the differences in community compositions among different treatments. With the infection of PWD, most of the bacterial groups tended to be co-excluding rather than co-occurring. These changes would correlate with microbial ability to suppress plant pathogen, enhancing the understanding of disease development and providing guidelines to pave the way for its possible management.


Introduction
Terrestrial ecosystems are confronted with more and more abiotic and bio-disturbances with the acceleration of global climate change. And extreme climate change aggravates the emergence of plant diseases and insect pests in global forest systems, which has brought great ecological and economic challenges [29]. In the forest ecosystem, the loss of leaves caused by diseases and pests leads to tree dieback and large-scale forest decline, resulting in changes of forest community structure [56], which in turn affects microbial communities [9] and ecosystem function [60]. Plant adaptation to the environment is the result of the integration of the plant itself and its microbiome, which is essential for maintaining the function of terrestrial ecosystems [12,15,51,80]. Over the years, pine trees have been confronted with a severe and devastating disease, pine wilt disease (PWD) mainly caused by Bursaphelenchus xylophilus [33], namely the pine wood nematode, which is one of the most serious conifer diseases worldwide, threatening several species of pine trees [78], and resulting in profound economic losses and adverse ecological environmental threat worldwide [67,78,85,95].
Trees infected with PWD have symptoms such as xylem deformation, resin duct disruption, and cortex and cambium tissue damage, which affect water transportation and conduction [26], resulting in the decline in photosynthesis, and ultimately leading to discoloration, wilting, and consequent death of host trees [27]. Given that, the biological characteristics of PWD [50,101], the dispersing vector [13,42,83], and the mechanism of PWD pathogenicity [63,103] have become research hot spots. In addition, there is growing evidence suggesting that during pathogenesis of B. xylophilus, plant microorganisms play important roles in host fitness [31,34,74,81,84]. The growing amount of data demonstrated that plant-related bacteria have beneficial effects, promote plant growth, and improve plant stress and disease resistance [36,68], particularly bacterial genera Trichoderma, Serratia, Bacillus, and Esteya which have nematicidal activity through mechanisms of parasitism or their production of toxic compounds [53,89]. In this perspective, there is a pressing need to illuminate the variations of plant microbiome during the disease development [61], advancing our understanding of the relationship between plant compartments and the microbial communities after B. xylophilus infection, and paving the way for its possible management.
In recent years, the endophytic microbial community of several B. xylophilus host pine trees, such as Pinus flexilis [11], Pinus contorta [7], Pinus pinaster (Proença et al., 2017a), and Pinus sylvestris [38], has been well documented. It is well established that pine endophytic bacterial diversity and composition play an important role in regulating plant response to PWN [3,54,61]. However, our understanding of the significant implications of phyllosphere and rhizosphere microorganisms remains limited. To the best of our knowledge, phyllosphere and rhizosphere microorganisms are two important components of plant microflora [39]. Previous studies have shown that the decline in plant healthy status or changes in growth conditions caused by host pathogens could affect the microbial community in leaves and roots of the host [23,76,85]. The phyllosphere microbiome interacts with the host plant affecting its health and function, and act as mutualists promoting plant growth and tolerance of environmental stressors [72]. The plant healthy status also can affect the exudates of the plant root, which is the essential factor affecting the rhizosphere microbial community [37,96]. The rhizosphere microbial community has an excellent ability to synthesize diverse secondary metabolites in response to different abiotic and biotic stresses, which is fundamental to the healthy growth of plants [6,54]. In addition, microorganisms also form a complex and diverse co-occurring network through direct or indirect interactions, and under the action of interspecies relationships such as symbiosis and competition, they produce positive and negative directions between each other, which are closely related to plant resistance [47]. In addition, previous researches collected the samples were only at one time point at the last stage of the disease [52]. Therefore, it is necessary to further investigate the microbial community at different stages of PWD after PWD occurrence and elucidate the relationship between the plant pathogen and host microbial community. P. koraiensis is widely distributed in Northeast China (Liaoning, Jilin, and Heilongjiang provinces), Japan, North Korea, South Korea, Russia, and other countries, with a total area of about 300,000 km 2 [59]. P. koraiensis, as a famous and valuable economic tree species, plays a momentous ecological environmental value. However, in P. koraiensis, a main host of B. xylophilus that is generally distributed in China, the plant microbiome information after PWD infection under field conditions has been scarcely studied. Here, the phyllosphere and rhizosphere bacterial and fungal communities in healthy P. koraiensis (PKa) and P. koraiensis naturally infected by B. xylophilus at the early stage (PKb) and at the last stage (PKc) of the disease were analyzed by sequencing 16S rDNA and ITS (internal transcribed spacer) rDNA using the high-throughput Illumina NovaSeq PE250 to uncover the differences in host microbial community potentially caused by PWD. In this study, we hypothesized that [2] PWD could increase the diversity of phyllosphere and rhizosphere microbial communities that differed as the infection of B. xylophilus progressed; [3] the host bacterial and fungal community differed as the infection of B. xylophilus progressed; and [4] with the infection of PWD, most of the microbial taxa tended to be co-excluding; besides, some microbial species unsuited towards living in infected pines disappeared and some species would present. This present study clarified the shifts of the host microbial community of P. koraiensis caused by PWD for the first time, so as to better understand the relationships between pathogens and the host microbial community at different stages of PWD after the occurrence of PWD. What's more, microbial community composition and microbial network can be used as biological indicators at different stages after the occurrence of PWD.

Overview of the Research Area
The study area is located at Dengta City, Liaoyang City, Liaoning Province, China (41°17′44″ N, 123°35′47″ E). The climate in this area is characterized as north temperate continental climate with an annual average temperature of 8.8 °C, annual average precipitation of 600 to 800 mm, and an annual average frost-free period of 140 to 160 days. The soil type is classified as Eutrochrepts soil [71]. Five fixed sites with the same site conditions and soil basic properties were selected, with an area of 1 ha. Eight healthy P. koraiensis trees (PKa), eight diseased P. koraiensis trees infected by B. xylophilus at the early stage (PKb), and eight diseased P. koraiensis trees at the last stage (PKc) were selected for sampling in each site. The distance between diseased and adjacent healthy trees was less than 15 m. The selection of healthy and P. koraiensis trees infected by B. xylophilus was made in accordance with the method described by Millberg et al. [57], in which the healthy trees were completely green needles and from which no B. xylophilus was isolated. The early stage of infection tree refers to needles that have become slightly wilted and browning of needles. The last stage of diseased trees was completely dead looking with brown needles (Fig. S1).

Sample Collection
For the leaf sampling, 24 needle samples were collected from the tops of the branches from three directions (120° as the boundary) from 8 trees, and mixed as one replicate at each site. A total of 15 leaf samples (3 types × 5 replicates) were collected from the study area. With regard to soil sampling, 24 rhizosphere soils were collected from three directions from 8 trees, and pooled together as one replication in each site, resulting in a total of 15 rhizosphere soil samples (3 types × 5 replicates). All the samples were put in the icebox and transported to the laboratory for subsequent analysis.
To analyze the microbial community on the leaves according to Kembel and Mueller [43], and Ren et al. [65], 10 g of leaf samples from each replicate was cut into pieces and transferred to a sterile triangle flask. 1:20 (leaf weight/ volume TE buffer =1:20) phosphate-buffered saline solution (20 ml, PBS, 0.01 M, pH 7.4) was added to each triangle flask. After sealed with a sterilized film, the samples were shaken on a shaker at 200 r/min for 30 min at room temperature, and the microbial cells were separated from the leaf surface. Vacuum filtration was used in a sterile environment, and microbes from the oscillating liquid were collected on a 0.22-μm microporous membrane placed into 2-ml sterile centrifuge tubes, then stored at −80 °C prior to DNA extraction. Fresh soil removed plant residues and stone was passed through a 2-mm sieve and immediately put into 2-ml sterile centrifuge tubes and frozen at −80 °C for later DNA extraction and high-throughput sequencing.

DNA Extraction, Amplification, and NovaSeq PE250 Sequencing
Genomic DNA was extracted from microporous membranes and 0.5 g soil using the FastDNA SPIN Kit for soil (MP Biomedicals, Santa Ana, CA, USA), in accordance with the manufacturer's instructions. The concentrations of DNA were measured using a NanoDrop ND-1000 spectrophotometer (Thermo Fisher Scientific, USA). The V3-V4 regions of the bacterial 16S rDNA gene were amplified and sequenced using the primer pairs 338F (5′-ACT CCT ACG GGA GGC AGC AG-3′) and 806R (5′-GGA CTA CHVGGG TWT CTAAT-3′) with barcode sequence. And ITS1 regions of the fungal ITS rDNA gene were amplified using the primer pairs ITS1F (5′-CTT GGT CAT TTA GAG GAA GTAA-3′) and ITS2 (5′-GCT GCG TTC TTC ATC GAT GC-3′) with barcode sequence [20]. All the PCR were carried out with 25-μl mixture, including 2 μl of dNTPs (2.5 mM); 2 μl DNA template (40-50 ng); 0.25 μl (5 U/μl) of Q5 High-Fidelity DNA Polymerase; 8.75 μl of ddH 2 O; 5 μl of Q5 High-Fidelity GC buffer (5×) and Q5 reaction buffer (5×), respectively, 1 μl (10 uM) of forward primer; and 1 μl (10 uM) of reverse primer. The following PCR thermal cycling conditions consisted of an initial denaturation step for 5 min at 98 °C, followed by 25 cycles of denaturation for 15 s at 98 °C, annealing for 30 s at 55 °C, and elongation of 72 °C for 30 s, with the final elongation step for 5 min at 72 °C. The PCR amplicons were further purified and quantified by using Agencourt AMPure Beads (Beckman Coulter, Indianapolis, IN) and PicoGreen dsDNA Assay Kit (Invitrogen, Carlsbad, CA, USA). PCR products for sequencing were carried out using an Illumina's NovaSeq PE250 platform at Shanghai Personal Biotechnology Co., Ltd., Shanghai, China. The high-throughput sequencing raw data of phyllosphere and rhizosphere microbes were uploaded in the NCBI database with the SRA accession numbers of PRJNA689361 and PRJNA689392.

Bioinformation Analysis
After removing primers and barcode sequences with cutadapt, quality filter, denoise, joint, and removal of chimeras, the high-quality sequences were finally obtained. Sequences with ≥97% similarity were assigned to the same OTU. The Silva Database for bacteria (Release132, http:// www. arbsilva. de) [64] and Unite Database for fungi (Release 8.0, https:// unite. ut. ee/) [46] were used for each representative sequence.

Statistical Analysis
One-way analysis of variance (ANOVA) with an LSD test was used to identify differences in microbial community richness (Chao 1 index, Observed species), diversity (Shannon index, Simpson index), and evenness (Pielou_e index) among different treatments. Venn diagrams were constructed using subsampled data to show the shared and unique OTUs in RStudio with the "Venn" package. Linear discriminant analysis Effect Size (LEfSe) in Galaxy software was employed to identify microbial lineages (from the phylum to genus) responsible for the differentiation of the phyllosphere and rhizosphere microbial communities caused by different treatments. Principal coordinate analysis (PCoA) based on the Bray-Curtis distance matrix, as one of the classical multidimensional scaling, was used to visualize the distinction of the microbial community structure. Heatmap plots of phyllosphere and rhizosphere microbial communities with the relative abundance of top 50 at the genus level were performed based on the Bray-Curtis distance matrix using RStudio with the package "Vegan". The co-occurrence patterns of OTUs from different treatments were evaluated by network analysis using the "psych" package in RStudio based on the Spearman rank correlation, and Gephi software was applied to visualize the networks with a Fruchterman-Reingold layout.

Changes in Phyllosphere and Rhizosphere Microbial Community Alpha Diversity
A total of 1,367,262 and 1,174,114 high-quality phyllosphere and rhizosphere bacterial sequences were generated across all samples after sequence denoising and quality filtering with the average number of sequences per sample 91,150 and 78,274, severally, which were assigned into 11,294 and 18,175 OTUs. The number of shared phyllosphere bacterial OTUs among PPKa, PPKb, and PPKc was 1306, and the unique OTUs of PPKa, PPKb, and PPKc were 2400, 3028, and 2862, respectively (Fig. 1A). The number of shared rhizosphere bacterial OTUs among RPKa, RPKb, and RPKc was 2297, and the unique OTUs of RPKa, RPKb, and RPKc were 3416, 3145, and 5437, respectively (Fig. 1B). Moreover, the shared OTUs among RPKa, RPKb, RPKc, PPKa, PPKb, and PPKc were 13 (Fig. S2A).
The fungal communities were further explored by highthroughput amplicon sequencing. Across all samples, we obtained a total of 1,318,977 and 1,316,090 high-quality phyllosphere and rhizosphere fungal sequences after sequence denoising and quality filtering with the average number of sequences per sample 87,739 and 87,931, severally, which were respectively grouped into 1272 and 3190 OTUs. The number of shared phyllosphere fungal OTUs among PPKa, PPKb, and PPKc was 276, and the number of unique OTUs of PPKa, PPKb, and PPKc was 599, 846, and 997, respectively (Fig. 1C). The number of shared rhizosphere fungal OTUs among RPKa, RPKb, and RPKc was 269, and the number of unique OTUs of RPKa, RPKb, and RPKc was 272, 260, and 245, severally (Fig. 1D). In addition, the shared OTUs among RPKa, RPKb, RPKc, PPKa, PPKb, and PPKc were 23 (Fig. S2B).

Variations in Phyllosphere and Rhizosphere Microbial Community Beta Diversity
It is well established that the microbial compositions from rhizosphere and phyllosphere samples formed distinct clusters (Fig. 2), which indicated that the plant compartment is a major selective force for the formation of plant-related microbial composition. The unconstrained principal coordinate analysis (PCoA) of the Bray-Curtis distance from all phyllosphere and rhizosphere samples based on the OTU data detected 64.9% of the total variance among bacterial communities, with the first and second axes explaining 57.2% and 7.7% of the variance, respectively ( Fig. 2A). PCoA based on the OTU data detected 67.0% of the total variance of fungal communities, with the first and second axes explaining 55.7% and 11.3% of the variance, respectively (Fig. 2B). As expected, the infection of B. xylophilus had a profound effect on plant microbe. The results demonstrated that rhizosphere bacterial community (Fig. S3A), rhizosphere fungal community (Fig. S3C), and phyllosphere fungal community (Fig. S3D) from PKa, PKb, and PKc formed three distinct clusters, especially along the PCoA1.

Microbiological Information Network and co-Occurrence Analysis
In order to further disentangle complex microbe-microbe interactions, we created association networks of phyllosphere and rhizosphere bacterial and fungal communities from OTU data ( Fig. 7; Table S1). Total nodes of phyllosphere and rhizosphere bacterial community association network in PKc existed the highest, followed by PKb and PKa, and total nodes of phyllosphere and rhizosphere fungal community association network in PKb existed the highest, followed by PKc and PKa, both indicating that the OTUs of the ecological network increased after infection ( Fig. 7; Table S1). Graph density in the network of PKb, a key topological property to describe how well a node is connected with its neighbors, showed higher than PKa and PKc, suggestive of more intensive microbial coupling at the early stage of infection ( Fig. 7; Table S1). Except for phyllosphere fungi, positive links showed decreased with the infection of PWD, and at the last stage, positive links existed the lowest ( Fig. 7; Table S2), demonstrating that most of the microbial taxa tended to be co-excluding rather than co-occurring.

Microbial Community Diversity Response to Different Samples
Plant microorganisms play critical roles in ecosystem function, sustainable restoration and management, as well as health of many plant species ( [ 18,45,58,61,77]. As climate change and human activity disrupt natural environments and microbial processes, there is essential to further explore the variations of microbe-microbe interactions [30] and microbe-host interactions [17]. We investigated the microbial community of phyllosphere and rhizosphere from healthy and diseased pine trees naturally infected by B. xylophilus at the different stages under field conditions. In our study, 11,294 and 18,175 phyllosphere and rhizosphere bacterial OTUs, and 1272 and 3190 phyllosphere and rhizosphere fungal OTUs of healthy and diseased pines were detected. In almost all samples, the rhizosphere bacterial Chao 1 index, Pielou_e, Shannon, and Simpson index were much higher than the respective phyllosphere communities (Table 1), which was a common finding in similar studies of native and cultivated plants in different environments [8,22,102]. The differences in microbial community diversity between the two plant compartments might account for the direct influence of their surrounding environment, and their fundamental discrimination of physiology and function [28]. Mounting empirical evidences have suggested that root exudates have a strong detrimental role in selecting the growth of specific bacteria [10,88] through signal transmission of microbe-microbe and plant-microbe interactions [82], ultimately promoting the differentiation of the bacterial assemblages [8]. Additionally, phyllosphere exists generally lower bacterial richness and abundance due to the fluctuations in environmental pressures [79,87]. In regard to fungal community, phyllosphere fungal community diversity presented higher than rhizosphere fungal community diversity (Table 2), which was not in agreement with previous researches from Chen et al. [14] and Jia et al. [40]. Our results nicely demonstrated that the effects of root and leaf compartments on the α-diversity indices of fungal community were different from those of bacterial community [48,73].
What's more, phyllosphere and rhizosphere microbial community diversity between healthy and infected P. koraiensis presented obvious difference. At the early stage of the infection, rhizosphere and phyllosphere bacterial Pielou_e and Simpson index, rhizosphere fungal Shannon and Simpson index, as well as phyllosphere fungal Chao 1 index exhibited slightly higher than those of healthy P. koraiensis (Table 1; Table 2). At the last stage of the infection, rhizosphere and phyllosphere bacterial Pielou_e, Shannon, and Simpson index; rhizosphere fungal Pielou_e, Shannon, and Simpson index; and phyllosphere fungal Chao 1 index existed abundantly higher than those of healthy P. koraiensis and the early stage of the infection (Tables 1 and  2). Our findings were consistent with a previous study from Proença et al. [62] who demonstrated that the endophytic bacterial diversity of P. pinaster tree was the highest at the late stage of pine wood nematode infection, and there was no conspicuous difference in bacterial diversity at the early stage of the disease, while the research from Ma et al. [52] indicated that there were no significant differences of rhizospheric bacterial diversities between healthy and wilted pines. Besides another investigation that showed that B. xylophilus infection appeared to reduce soil bacterial diversity [69], similar findings were reported by Zhang et al. [97] who demonstrated that B. xylophilus infection likely decreased the richness and diversity of endophytic microbes. It thus appeared that the inconsistent results might be due to different tree species and the sampling period after the infection of PWD. In the present study, we revealed that PWD could increase the phyllosphere and rhizosphere bacterial and fungal diversity and the microbial community diversity differed as the disease progressed, suggesting the importance of the host microbiome in disease development. The differences might be caused by the growing abundance of the dominant microbial groups crowding out the weaker microbial groups or that microbial species unsuited towards living in infected pines disappeared.

Microbial Community Composition Response to Different Treatments
As shown by a growing body of works [14,75,90], we also observed that the microbial compositions from rhizosphere and phyllosphere samples formed distinct clusters. Collectively, these studies suggested that although the assemblies of rootassociated bacteria and fungi differ substantially from the phyllosphere microbial communities, both represent a subset of the microbe derived from soil communities and enriched in different plant-associated niches [16,32]. As previous findings indicated that the infection of plant pathogens could affect the host microbial community [52,76], we also documented that the PWD had a profound impact on the host rhizosphere bacterial and fungal community and phyllosphere fungal community, which was not complete in line with previous work that demonstrated that the community structure of healthy and diseased trees was only significantly different in the roots, and not in the needles and soil [52]. It has become evident that root exudates are the essential factor determining the structure of the rhizosphere microbial community [5,91]. The occurrence of pine wilt disease can lead to a decreased secretion of soluble sugar, total sugar, and protein in roots [66], which might have caused the observed difference in the microbial community structure in the rhizosphere.
Intriguingly, overall, the bacterial community compositions were similar (in terms of dominant groups) in all samples, different plant compartments at different stages of disease dominated mainly by Proteobacteria, followed by Actinobacteria, and this finding was consistent with several previous studies that displayed that Proteobacteria and Actinobacteria were the dominant groups in rhizosphere bacterial communities [2,41] and phyllosphere microorganisms [19,86]. In addition, these groups represent ubiquitous rhizosphere taxa were detected in various stressed environments [92], while the opposite observations from Pinus massoniana infected by B. xylophilus showed that Acidobacteria was the predominant species in infected soils [69]. Interestingly, due to B. xylophilus infection of P. koraiensis, the relative abundances of Acidobacteria, Bacteroidetes, and Proteobacteria were significantly higher in diseased pine roots, and the shifts of Proteobacteria have been observed in previous findings [52], which collectively demonstrated that Proteobacteria might be phytopathogens and parasites in plant tissues and cause a variety of diseases [44]. The root and leaf metabolism of diseased trees was weakened relative to the healthy roots and leaves, resulting in a decreased ability of the root and leaf to adapt to the environment condition and it being easily colonized by microbes. Other investigations indicated that Proteobacteria prefer to grow under nutrient-rich conditions [24], which might explain the high content of Proteobacteria in the diseased roots and leaves. The rhizospheric microbial abundance of Bryobacter, RB41, and Bradyrhizobium was richer in diseased pines. Our findings were similar to the results in rhizosphere bacterial studies on P. thunbergii where bacteria in the genus Bradyrhizobium were more abundant in soil of wilted trees than in soil of non-infected trees [52].
The abundances of the genus Massilia, Sphingoaurantiacus, Acidiphilium, Acetobacteraceae, Singulisphaera, Phascolarctobacterium, and Hymenobacter in diseased needles were significantly higher than those in healthy needles, suggesting an association of particular microbial abundances with the infection of B. xylophilus in P. koraiensis. What's more, the research from Ma et al. [52] demonstrated that Massilia was obviously higher in diseased pines, which supported our results to some extent. The genus Massilia belongs to the family Oxalobacteraceae of the class Betaproteobacteria in the phylum Proteobacteria [1]. Members of this genus are characterized as Gram-negative, aerobic, non-spore-forming bacteria [99]. Some Massilia can produce cell lysis enzymes that promote tissue lysis [55]. This may be the reason for the presence of Massilia in a high abundance in diseased needles.
In our study, Ascomycota and Basidiomycota were the dominant fungal phyla with phyllosphere and rhizosphere samples, and this result was in agreement with previous research [40]. Similar results were obtained in Taxus rhizosphere communities [35] and in tropical grasslands [49]. The majorrhizosphere fungal genera in healthy P. koraiensis were Penicillium and Trichoderma, in agreement with a study from Zhang et al. [97]. Interestingly, Trichoderma is an important genus in biocontrol of nematodes because some species produce metabolites harmful to nematodes [93]. The enriched phyllosphere fungal groups in PPKc were Phaeosphaeriaceae, Wilcoxina, Pseudocosmospora, and so on. Phaeosphaeriaceae was commonly associated with plants as pathogens, though some are also saprotrophs and parasites on powdery mildews [98]. Thus, it can be seen that plant-associated microbes could influence plant health and fitness [94], resistance to pathogens [4], and ecosystem services.  6 LEfSe analysis to identify which microbial taxa (from phylum to genus level) were major contributors to the differences in rhizosphere bacterial (A), phyllosphere bacterial (B), rhizosphere fungal (C), and phyllosphere fungal (D) community compositions among different samples. PPKa, the phyllosphere of healthy Pinus koraiensis; PPKb, the phyllosphere of P. koraiensis naturally infected by Bursaphelenchus xylophilus at the early stage; PPKc, the phyllosphere of P. koraiensis naturally infected by Bursaphelenchus xylophilus at the last stage; RPKa, the rhizosphere of healthy P. koraiensis; RPKb, the rhizosphere of P. koraiensis naturally infected by Bursaphelenchus xylophilus at the early stage; RPKc, the rhizosphere of P. koraiensis naturally infected by Bursaphelenchus xylophilus at the last stage ◂

Shifts of co-Occurrence Association Network Response to Different Treatments
In our study, the co-occurrence network of phyllosphere and rhizosphere microbial community demonstrated dynamical relationships between healthy P. koraiensis and the infection of B. xylophilus in P. koraiensis, which could provide momentous details of microbial community assembly and represent interactions among different populations that regulate ecological processes [25]. Total nodes of phyllosphere and rhizosphere microbial community association network increased after B. xylophilus infection, indicating that the populations of the ecological network increased after infection ( Fig. 7; Table S1), resulting in the microbial diversity to increase in some degree. The edges of phyllosphere and rhizosphere microbial community association network  [70]. Furthermore, the role of microbial co-occurrence networks is important in revealing the interactions (such as through parasitism, competition, and mutualism) that exist among different species [21,100]. In our study, except for phyllosphere fungi, positive links of phyllosphere bacteria, rhizosphere bacteria, and fungi decrease with the infection of PWD, and at the last stage, positive links existed the lowest ( Fig. 7; Table S1), demonstrating that most of the microbial taxa tended to be in competition rather than mutualism.

Conclusions
Overall, an increase in diversity with more severe symptomatic stage was visible. What's more, the microbial compositions from rhizosphere samples and phyllosphere samples formed distinct clusters. Rhizosphere bacterial and fungal community, and phyllosphere fungal community from PKa, PKb, and PKc formed three distinct clusters, which clearly separated along the PCoA1. These findings manifested that the phyllosphere and rhizosphere microbial community changed potentially caused by B. xylophilus infection of P. koraiensis. Furthermore, LEfSe analysis demonstrated that variations of some microbial abundances were associated with the infection of B. xylophilus in P. koraiensis, including Bradyrhizobium (rhizosphere bacteria), Massilia (phyllosphere bacteria), and Phaeosphaeriaceae (phyllosphere fungi). With the infection of PWD, most of the bacterial taxa tended to be co-excluding rather than co-occurring. Together, our results explored PWD could increase the phyllosphere and rhizosphere microbial community diversity and microbial community composition differed as the disease progressed, and these changes would correlate with microbial ability to suppress plant pathogen. This study expanded our knowledge of the ecology of plant-microbe interactions as well as the structure and assembly of microbial communities of healthy P. koraiensis and the infection of B. xylophilus in P. koraiensis, which lay the foundation for studies that aim at improving plant growth by altering the plant microbiome.