Aquatic Macrophytes and Local Factors Drive Bacterial Community Distribution and Interactions in a Riparian Zone of Lake Taihu

Aquatic macrophytes rhizosphere are biogeochemical cycling hotspots in freshwater ecosystems. However, little is known regarding the effect of aquatic macrophytes on bacterial community and interactions in the riparian zones. We investigated the bacterial community composition and network structures along a gradient of the riparian zone as follows: The supralittoral and eulittoral zones with Phragmites australis, the eulittoral and infralittoral zones without P. australi. The bacterial communities in the four zones differed significantly based on taxonomic dissimilarity, but the two zones with P. australis exhibited phylogenetic closeness of the bacterial communities. The characteristics of the bacterial networks, such as connectivity, modularity, and topological roles of OTUs, were totally different between the P. australis and non-P. australis zones. Some bacterial phyla enriched in the P. australis zones were found to be putative keystone taxa in the networks, which might be involved in the regulation of bacterial interactions and plant growth. Moreover, the hydrological regime and particle size were shown to be determinants of the bacterial community and network structures in the riparian zones. In summary, our results show that the role of P. australis and local factors are crucial for constructing bacterial community and interactions in the riparian zones of lakes.


Introduction
In the context of lake hydrological changes, the riparian zone, the interface between the terrestrial and aquatic domains, can provide a heterogeneous environment and diverse habitats for variety of flora and fauna. The riparian zones of lakes can also maintain the stability of lake ecosystems by receiving and transferring the majority of terrestrial allochthonous input [1]. Aquatic macrophytes are an important component in complex ecosystems of the riparian zones of lakes [2]. Meanwhile, the rhizosphere of aquatic plant has also become widely recognized biogeochemical hotspots because of interactions between roots and microorganisms [3,4], which effectively enhances the function of wetlands in purifying environmental contaminants [5]. The function of the riparian zones of lakes is inseparable from the diversity and structure of riparian zone-associated microorganisms. Aquatic macrophytes and local factors are crucial for the bacterial community distribution and assembly along the gradient of riparian zones. Although numerous studies have shown that the environment factors of riparian zones are crucial for microbial community composition and diversity [6][7][8], we lack an understanding of the effect of aquatic macrophytes on riparian bacterial community composition and interactions in the riparian zones of freshwater lakes.
Plant roots interact with surrounding soil/sediment, resulting in plant species-specific rhizosphere microorganisms, which has consequences for variation in microbial community composition and interactions [9]. Aquatic plant root exudates, such as organic acids, sugars, and other secondary metabolites, have been regarded as a driving forcing in determining microbial community composition and activity [10], especially promoting nitrogen cycling processes [11]. With the development of high-throughput sequencing technology, 16S rRNA gene sequencing has provided a deeper understanding of aquatic plant rhizosphere microbiota. Thus, more taxa of rhizosphere-associated microbiota have been shown to contribute to the processes of pollutant degradation and plant fitness in wetlands [4,12]. Some studies have also shown that the species and lineages of aquatic plants play a crucial role in the distribution and composition of microbial communities [4,13,14]. Moreover, interactions among members of a microbial community are an indispensable part the construction of a microbial community and impacts its functions [15]. Network analysis is an effective method to predict correlations between microbes and the dynamics of microbial interactions [16,17]. Recently, an increasing number of studies have investigated the influence of the rhizosphere on bacterial networks, with the results showing that the characteristics of rhizosphere microbial community networks are completely different from those of the surroundings, including properties such as connectivity, potential links, and complexity [18][19][20]. The putative key species of rhizosphere microbial networks have been recognized to be involved in plant growth and health [21,22].
Microbial communities are governed by abiotic (environment factors) and biotic (species and diversity of plant) factors in different ecosystems [9]. The characteristics of riparian microbial communities are typically considered to be highly sensitive to environment factors, such as nitrogen content, organic matter, and pH [8,23,24]. Dynamic littoral hydrology can also predominately affect riparian bacterial community assembly and carbon cycling [6,7,25]. In addition, aquatic macrophytes are adapted to a changing environment by being able to regulate their physiological characteristics. For example, the relatively high evapotranspiration rates of P. australis can help protect plants and affects nutrient uptake under flooding conditions [2,26]. In addition, aquatic macrophytes gradually modify the physiochemical properties of soil/sediment in a variety of ways, such as nutrient uptake by roots, plant residue decomposition, and root exudates [27][28][29]. These results suggest that interactions between aquatic plants and local factors are crucial to the distribution and interactions of bacterial communities in the riparian zones of lakes.
P. australis is widely distributed in land-lake ecotones, and larger roots systems promote its effect over a wide range [2]. In particular, P. australis can inhibit the growth of other aquatic macrophytes [30], which is useful for investigating the effects of single aquatic macrophytes on bacterial communities along the gradient of the riparian zones of lakes. P. australis grow in different hydrology conditions along the riparian zone. Meanwhile, gradient environment conditions of riparian zone can be considered as a good model to evaluate the relative contribution of local factors and P. australis on bacterial community characteristics of riparian zone. In this study, 16S rRNA gene amplicon sequencing and co-occurrence network analysis [31] were performed to explore bacterial community and network structures in the riparian zones of a lake. We hypothesized that the P. australis rhizosphere selection exerts a strong effect on bacterial community and network structures. The objectives of this study were (i) to understand the P. australis rhizosphere bacterial community composition and core species of bacterial networks; and (ii) to elucidate the correlations among P. australis, local factors, and bacterial communities in the riparian zone of freshwater lake.

Site Description and Sampling
Lake Taihu, a shallow eutrophic lake in a region of the Yangtze Delta in China, has a mean depth of 1.9 m and a surface area of 2338 m 2 . The study area is located in a riparian zone in the southeastern portion of Lake Taihu (Figure 1), which is covered with lush P. australis. Given the environmental heterogeneity of regions within the riparian zone of Lake Taihu with respect to hydrological characteristics and the presence of P. australis, we divided the area into two parts along the gradient of riparian zone (30 m long), where one zone with P. australis is the rhizosphere zone (RZ) and another zone without aquatic macrophytes is the nonrhizosphere zone (NZ). In addition, the two zones were subdivided into four sites, and six replicated plots (1 × 1 m) were placed along each transect (Figure 1), as follows: The supralittoral zone with P. australis (RZa) in the nonflooded zone; the eulittoral zone with P. australis (RZb), and the eulittoral zone without any aquatic macrophytes (NZc), which are closely adjacent to each other and located at the middle areas of the riparian zones of the lake, where the zones are irregularly inundated and regulated by the hydrological characteristics of Lake Taihu; and the infralittoral zone without any aquatic macrophytes (NZd), which is located at the lowest region and is inundated. Twenty-four composite samples were collected in May 2016. At each plot, three parallel samples were randomly collected and evenly combined into a single composite sample. The RZ samples were collected by lightly shaking roots [32]. The NZ samples (5-10 cm, similar depth as rhizosphere samples) were collected using a polymer glass tube (6 cm in diameter and 60 cm in length). All the samples were stored on ice and transported to the laboratory. For each sample, one portion was stored at −80 • C for molecular analysis, while another portion was stored at 4 • C to assess the chemical properties.

Sample Chemical Properties
For all the samples, the sample water content (SWC) was measured as gravimetric weight loss after drying at 105 • C for 24 h and reweighing, while the organic matter content (OM) was estimated by the loss-on-ignition method after combustion at 550 • C in a muffle furnace for 4 h. The pH values of the samples were determined using a soil to water ratio of 1:5 (v/v) with a pH meter. Total carbon (TC) and total nitrogen (TN) were analyzed with an elemental analyzer (Euro Vector EA3000) and then used to calculate the TC/TN. The dissolved organic carbon (DOC) was measured for the sample extracts (5 g freeze-dried samples, 50 mL pure water) in a TOC analyzer. Particle size distributions (silt, clay, and sand) were determined using a laser particle size analyzer. To assess ammonium nitrogen (NH 4 + ) and nitrate nitrogen (NO 3 − ) contents, samples were extracted with 2 M KCl (1:10 sample: Extraction ratio) and the extract was analyzed using a continuous flowing analyzer.

S rRNA Gene Pyrosequencing and Bioinformatics Analysis
Total DNA was extracted from 0.5 g of dry soil/sediment using a FastDNA SPIN for Soil kit (MP Biomedicals, Santa Ana, CA, USA) according to the manufacturer's instructions. The V4-V5 hypervariable regions of the bacterial 16S rRNA gene was amplified using specific primers with the barcodes (515F GTGCCAGCMGCCGCGGTAA and 907R CCGTCAATTCCTTTGAGTTT). All PCR-amplification reactions were in a 30 µL volume with 15 µL of Phusion ® High-Fidelity PCR Master Mix (New England Biolabs, Ipswich, MA, USA), 0.2 µM forward and reverse primers, and approximately 10 ng of template DNA. The thermocycling conditions were as follows: 98 • C for 1 min, followed by 30 cycles of denaturation at 98 • C for 10 s, annealing at 50 • C for 30 s, and elongation at 72 • C for 30 s, with a final extension at 72 • C for 5 min. The PCR products were assessed by electrophoresis on a 2% agarose gel, and a mixture of PCR products was purified with a GeneJET Gel Extraction Kit (Thermo Scientific, Waltham, MA, USA). Sequencing libraries were generated using an Illumina TruSeq DNA PCR-Free Library Preparation kit (Illumina, San Diego, CA, USA) following the manufacturer's recommendations, and index codes were added. The library quality was assessed Water 2020, 12, 432 4 of 16 using a Qubit@ 2.0 Fluorometer (Thermo Scientific) and an Agilent Bioanalyzer 2100 system. Finally, the library was sequenced on an Illumina HiSeq platform, and 250 bp paired-end reads were generated. All of the sequences were deposited in to the NCBI Sequence Read Archive (SRA) database (Accession Number: SRP235199).

S rRNA Gene Pyrosequencing and Bioinformatics Analysis
Total DNA was extracted from 0.5 g of dry soil/sediment using a FastDNA SPIN for Soil kit (MP Biomedicals, Santa Ana, CA, USA) according to the manufacturer's instructions. The V4-V5 hypervariable regions of the bacterial 16S rRNA gene was amplified using specific primers with the barcodes (515F GTGCCAGCMGCCGCGGTAA and 907R CCGTCAATTCCTTTGAGTTT). All PCRamplification reactions were in a 30 μL volume with 15 μL of Phusion ® High-Fidelity PCR Master Mix (New England Biolabs, Ipswich, MA, USA), 0.2 μM forward and reverse primers, and approximately 10 ng of template DNA. The thermocycling conditions were as follows: 98 °C for 1 min, followed by 30 cycles of denaturation at 98 °C for 10 s, annealing at 50 °C for 30 s, and elongation at 72 °C for 30 s, with a final extension at 72 °C for 5 min. The PCR products were assessed by electrophoresis on a 2% agarose gel, and a mixture of PCR products was purified with a GeneJET Gel Extraction Kit (Thermo Scientific, Waltham, MA, USA). Sequencing libraries were generated using an Illumina TruSeq DNA PCR-Free Library Preparation kit (Illumina, San Diego, CA, USA) following the manufacturer's recommendations, and index codes were added. The library quality was assessed using a Qubit@ 2.0 Fluorometer (Thermo Scientific) and an Agilent Bioanalyzer 2100 system. Finally, the library was sequenced on an Illumina HiSeq platform, and 250 bp paired-end reads were generated. All of the sequences were deposited in to the NCBI Sequence Read Archive (SRA) database (Accession Number: SRP235199).

Network Construction and Analysis
Networks were constructed for the RZ and NZ bacterial communities based on the relative abundances of OTUs. A random matrix theory (RMT)-based network analysis was performed using the molecular ecological network analysis pipeline, and additional information regarding the structure and theories of the pipeline are described elsewhere [31]. Only OTUs detected in six samples of 12 replicates were selected to construct the networks. The primary procedures performed in this study were as follows. First, the OTU table, environmental factors data, and OTU taxonomy file were prepared as guided in the pipeline. Second, a similarity matrix was constructed based on the Spearman's Rho between pairwise OTUs. A cutoff value (similarity threshold, St) for the similarity matrix was generated automatically following the default settings. When the correlation between the abundance of OTUs is larger than the St, a link is assigned between the pair of OTUs. Third, the "Global Network Properties", "Individual Nodes' Centrality", and "Module Separation and Modularity" functions were performed to calculate some network properties (e.g., average path, average degree, and average clustering coefficient). A module is a group of nodes that is more highly connected to each other than to nodes outside the group. Modularity (M) is a value that measures how well a network is divided into a module (M > 0.4 as the threshold) [33]. The network hub analysis was based on the connectivity of each node, which was determined by the within-module connectivity (Zi) and among-module connectivity (Pi). The node topologies are divided into four categories: Network hubs (Zi > 2.5 and Pi > 0.62), module hubs (Zi > 2.5 and Pi < 0.62), connectors (Zi < 2.5 and Pi > 0.62), and peripherals (Zi < 2.5 and Pi < 0.62). Fourth, "Output for Cytoscape visualization" was run, and three files were generated for network graph visualization using Cytoscape. Fifth, "Randomize the network structure and calculate network properties" was performed to compare the properties of the empirical and random-generated networks. Sixth, "Gene/OTU significances (GS)" with environmental traits were calculated, and a Mantel test was performed to assess the correlations between the GS and network connectivity. Finally, "Module-EigenGene analyses" were performed to identify relationships between local factors and network modules.

Statistical Analysis
One-way ANOVA with Tukey's post-hoc tests were used to assess the effects of the riparian zones on environmental characteristics, alpha diversity indices, and the relative abundances of phyla (significance level of 0.05). If the distribution of dates were non-normal, the Kruskal-Wallis ANOVA was used. OTUs with relative abundances of >0.2% were selected to study the differences in dominant OTUs in four sampling sites and were visualized using the "ggplot2" package. The characterization of microbial features differentiated along riparian sampling sites at different taxonomic levels was performed using the linear discriminant analysis (LDA) effect size (LEfSe) method for biomarker discovery, and we used an LDA > 3.0 to indicate significant differences [34]. We used the principal coordinate analysis (PCoA) to study the shifts in bacterial community structure based on Bray-Curtis/Weighted UniFrac distances permutational multivariate analysis of similarity (ANOSIM) and permutational multivariate analysis of variance (ADONIS) to test whether the bacterial communities in the four sites were significantly different. Redundancy analysis (RDA) was performed using R to explore the linkages between bacterial community structure and local factors (i.e., NH 4 + , NO 3 − , pH, OM, SWC, DOC, and silt). To elucidate the relationships between the local factors distance and bacterial community dissimilarities (Bray-Curtis/Weighted UniFrac distances), Mantle tests were performed in R. The analyses described above were performed in R.

Local Factors in Sampling Zones
The results revealed that the NO 3 − and DOC concentrations of sediments exhibited significant differences between the RZ and NZ in the following order: NO 3 − (RZ > NZ) and DOC (RZ < NZ) (Table S1). The NH 4 + concentrations in the NZ samples ranged from 21.87 to 63.36 mg·kg −1 , which was much higher than those observed in the RZ samples (5.85 to 5.89 mg·kg −1 ) (p < 0.05) (Table S1). However, the SWC was significantly different among four zones and increased in the order of RZa, RZb, NZd, and NZc (p < 0.05) (Table S1). For the pH values, the lowest concentrations were observed in the RZb samples (5.28 ± 0.30). Particle size, including clay, silt, and sand differed among the four zones (p < 0.05). From the supralittoral to infralittoral regions, the clay content increased, whereas the silt content notably decreased.

Bacterial Community Composition and Diversity in the Riparian Zone of Lake Taihu
The obtained OTUs were affiliated with 16 phyla and subphyla (relative abundance >1%) ( Figure S1). The abundances of Alphaproteobacteria, Acidobcteria, and Gemmatimonadetes were significantly higher in the RZa and RZb samples compared to the NZc and NZd samples, whereas Chloroflexi and Firmicutes were enriched in the NZc and NZd samples. Bacterroidetes was dominant in the NZc samples compared with the other zones. The richness and diversity of the bacterial communities did not exhibit remarkable differences among the four zones based on the Chao, Ace, Shannon, and phylogenetic diversity (PD) indices (p > 0.05) (Table S2). Furthermore, the 30 most abundant bacterial OTUs were visualized by hierarchical clustering analysis ( Figure 2). The results revealed that the top 30 OTUs in the RZa and RZb samples clearly grouped together, whereas the top 30 OTUs in the NZc and NZd samples exhibited obvious differences, suggesting that the P. australis rhizosphere had similar abundant bacterial OTUs in different riparian zones. appeared to be the primary discriminant clades.
Principal coordinate analysis of the taxonomic (Bray-Curtis distance matrix) and phylogenetic (Weighted UniFrac distance matrix) dissimilarity were performed to investigate patterns of separations among the bacterial communities in the four zones. We observed that RZ and NZ samples separated across the first principal coordinate explained 31.23% and 49.35% of the variability in the bacterial communities with respect to taxonomic-and phylogenetic-based dissimilarity, respectively ( Figure 4). Moreover, the ANOSIM and ADONIS results demonstrated that the structures of the bacterial communities among the four zones were significantly different based on taxonomic dissimilarity (p < 0.05), whereas the structure of the bacterial communities in the RZa and RZb samples were similar based on phylogenetic dissimilarity (p > 0.05) ( Table S3), indicating that P. australis roots may enrich phylogenetic closeness bacterial communities in different riparian zones.
To further explore the significance of discriminatively specific bacterial taxa, the LEfSe analysis was performed from the phylum to genus level among the four sampling zones using the default logarithmic (LDA) value of 3 ( Figure 3). The clade graph showed that 11 groups of bacteria were enriched in the RZa samples, namely, Gemmatimonadetes (from the phylum to genus levels), Rhizobiales (Pseudolabrys), and Caulobacterales (Phenylobacterium). The discriminative clade of RZb was affiliated with Alphaproteobacteria, Acidobacteria_Gp1, and Myxococcales. We observed that the RZa and RZb samples had similar biomarkers, such as those affiliated with Alphaproteobacteria. In contrast, Bacteroidetes (Cryomorphaceae) and Clostridium_sensu_stricto were biomarkers of the NZc samples. In the NZd samples, Chloroflexi (from the phylum to family levels), Actinobacteria, and Clostridium_XI appeared to be the primary discriminant clades.
Principal coordinate analysis of the taxonomic (Bray-Curtis distance matrix) and phylogenetic (Weighted UniFrac distance matrix) dissimilarity were performed to investigate patterns of separations among the bacterial communities in the four zones. We observed that RZ and NZ samples separated across the first principal coordinate explained 31.23% and 49.35% of the variability in the bacterial communities with respect to taxonomic-and phylogenetic-based dissimilarity, respectively ( Figure 4). Moreover, the ANOSIM and ADONIS results demonstrated that the structures of the bacterial communities among the four zones were significantly different based on taxonomic dissimilarity (p < 0.05), whereas the structure of the bacterial communities in the RZa and RZb samples were similar based on phylogenetic dissimilarity (p > 0.05) (Table S3), indicating that P. australis roots may enrich phylogenetic closeness bacterial communities in different riparian zones.

Relationships between Bacterial Communities and Local Factors in the Riparian Zone of Lake Taihu
To investigate the associations of bacterial community dissimilarity with local factors in the riparian zone of Lake Taihu, Mantel tests were performed ( Table 1). The results showed that the taxonomic and phylogenetic dissimilarity in the NZ samples was significantly correlated with local factors (p < 0.05). The local factors also exhibited a significant effect on the RZ-associated bacterial community taxonomic dissimilarity (r = 0.37, p < 0.05), whereas they were uncorrelated with the RZassociated bacterial community phylogenetic dissimilarity (r = 0.13, p > 0.05) ( Table 1). These results suggested that the selection of P. australis roots strongly influenced the bacterial community at the phylogenetic level. In addition, Mantel tests showed that the NZ-associated bacterial community dissimilarity was strongly correlated with NH4 + , DOC, OM, N%, and C%. We also observed that the SWC and particle size was associated with the RZ-and NZ-associated bacterial community dissimilarity (Table 1). Furthermore, RDA were performed to assess associations between local factors and bacterial community structure (Figure 4c). The RDA results showed that two axes (RDA1 and RDA2) explained 33.91% of the bacterial community variance. The results of Monte Carlo permutations indicated that the local factors NH4 + , NO3 − , SWC, DOC, OM, pH, and silt were significantly correlated with the bacterial community structures (p < 0.05) (Table S4).

Relationships between Bacterial Communities and Local Factors in the Riparian Zone of Lake Taihu
To investigate the associations of bacterial community dissimilarity with local factors in the riparian zone of Lake Taihu, Mantel tests were performed ( Table 1). The results showed that the taxonomic and phylogenetic dissimilarity in the NZ samples was significantly correlated with local factors (p < 0.05). The local factors also exhibited a significant effect on the RZ-associated bacterial Water 2020, 12, 432 8 of 16 community taxonomic dissimilarity (r = 0.37, p < 0.05), whereas they were uncorrelated with the RZ-associated bacterial community phylogenetic dissimilarity (r = 0.13, p > 0.05) ( Table 1). These results suggested that the selection of P. australis roots strongly influenced the bacterial community at the phylogenetic level. In addition, Mantel tests showed that the NZ-associated bacterial community dissimilarity was strongly correlated with NH 4 + , DOC, OM, N%, and C%. We also observed that the SWC and particle size was associated with the RZ-and NZ-associated bacterial community dissimilarity (Table 1). Furthermore, RDA were performed to assess associations between local factors and bacterial community structure (Figure 4c). The RDA results showed that two axes (RDA1 and RDA2) explained 33.91% of the bacterial community variance. The results of Monte Carlo permutations indicated that the local factors NH 4 + , NO 3 − , SWC, DOC, OM, pH, and silt were significantly correlated with the bacterial community structures (p < 0.05) (Table S4).

Bc Wu Bc Wu
All

Co-Occurrence Bacterial Networks in the Riparian Zone of Lake Taihu
To explore the bacterial interactions in the RZ and NZ samples, we constructed bacterial co-occurrence networks using the random matrix theory (RMT) [31] (Figure 5). As shown in Table 2, the similarity threshold (St) values were 0.86 and 0.89 for the RZ and NZ network topologies, respectively. The modularity of each co-occurrence network was higher than 0.4, suggesting that the two networks were modular [33]. Moreover, the average path (GD), average degree (avgK), and average clustering coefficient (avgCC) values of the two networks were significantly higher than the values of random networks, indicating the occurrence of small-word behavior [35] (p < 0.05). Overall, the co-occurrence networks were markedly different between the RZ and NZ samples, revealing that the RZ network had more nodes but fewer links than the NZ network. Furthermore, the NZ network had a significantly higher avgK value (5.646) than the RZ network (2.95) (p < 0.05), and the avgCC value of the NZ network (0.311) was also notably higher than that observed for the RZ network (0.161) (p < 0.05), suggesting that the NZ-associated bacterial community had a greater potential for interactions.   Table S5. Abbreviations: P. australi zone (RZ); non-P. australis zone (NZ).  The analysis of modular topological roles is important to identify keystone species that play distinct topological roles in a network. The nodes of the network module were classified into four categories based on their within-module connectivity (Zi) and among-module connectivity (Pi) values ( Figure 6). Most of the nodes in each network were peripherals that had most of their links inside their modules. There were one module hub and four connectors in the RZ network, whereas the NZ network had three module hubs and two connectors. The module hub and connectors of the RZ network were affiliated with the genus Sphingomonas (Alphaproteobacteria), the order Gp3 (Acidobacteria), the class Betaproteobacteria (Proteobacteria), and the order Acidimicrobiales (Actinobacteria), whereas the module hubs and connectors of the NZ network belonged to the family Anaerolineaceae (Chloroflexi), the genus Ignavibacterium (Ignavibacteriae), and the genus Pelosinus (Firmicutes), suggesting that these taxa were crucial in the construction of co-occurrence networks (detailed taxonomic information is provided in Table S5). Furthermore, the connector of the RZ network appeared at the peripheral of the NZ network (OTU_166), whereas the module hub of the NZ network occurred at the peripheral of the RZ network (OTU_115) ( Table S5), suggesting a role shift of the nodes between the two networks.   Table S5. Abbreviations: P. australi zone (RZ); non-P. australis zone (NZ).  Table S5. Abbreviations: P. australi zone (RZ); non-P. australis zone (NZ).

Relationships between Network Connectivity and Local Factors in the Riparian Zone of Lake Taihu
Overall, the Mantel test results revealed that the connectivity of the RZ network was uncorrelated with local factors, whereas that of the NZ network was significantly correlated with local factors (Table S6). Partial Mantel test results showed that pH was significantly correlated with the RZ network, whereas the variables NH 4 + and DOC affected the NZ network. Notably, the silt and sand contents were both correlated with the two networks. Furthermore, network module eigengene analysis was performed to analyze the correlations between modules and local factors (Table S7). For the RZ and NZ networks, the eigengenes of the primary modules could explain 54-77% and 63-83%, respectively. pH was only correlated with the modules of RZ, whereas NH 4 + , NO 3 − , TC, OM, and DOC were only correlated with the NZ modules. Moreover, the variables TN, SWC, clay, silt, and sand affected some modules in both the RZ and NZ networks.

Roots Enrich and Filter the Bacterial Community Composition in the Riparian Zone of Lake Taihu
Similar to the studies on the P. australis rhizosphere [12,14,36], Alphaproteobacteria and Acidobcteria were substantially enriched in the P. australis zones. Alphaproteobacteria exhibit fast growth by utilizing a broad range of root-derived carbon substrates [37], with members of this phylum being widely distributed in different types of plant rhizosphere soils, such as rice [38], lettuce [39], rabbiteye blueberry [20], and wheat [19]. Acidobacteria play an important role in the carbon cycle by performing functions such as cellulose and lignin degradation [40], and Acidobacteria has also been shown to be a primary phylum in rhizosphere soil [41,42]. Moreover, Rhizobiales and Gemmatimonadetes were shown to be biomarkers of the P. australis zones (Figure 3), members of which are associated with the metabolism of nitrogen and phosphorus [43,44], suggesting that the P. australis rhizosphere microbiota may be involved in promoting plant nutritional requirements. In contrast, Chloroflexi and Firmicutes were primarily dominated in non-P. australis zones (Figure 3 and Figure S1), consistent with the results of previous studies on anaerobic sediments, such as in Lake Taihu [45] and freshwater wetlands [19].

Co-Occurrence Network Structures in the Riparian Zone of Lake Taihu
Network analysis is crucial to understanding the community dynamics of bacterial interactions and niches [15,16,46], providing a new dimension to our understanding of plant rhizosphere bacterial community ecology [18]. For the RZ and NZ networks, positive links accounted for 78.8% and 87.7% of the total links, respectively, suggesting that the nodes (i.e., OTUs) of the two networks represent cooperative interactions (interdependency) or a response to a common local factor [16,18]. Similar to the results of a study on wheat [19], the bacterial network of the P. australis zones exhibited lower connectivity than the surroundings. In contrast, rhizosphere bacterial networks have shown to be of greater complexity and connectivity compared with nonrhizosphere bacterial networks as plants grow, where the flowering stage of plants is particularly more complex than the vegetative stage [18]. These contradictory results indicate that plant species, growth stages, and habitats can affect the characteristics of bacterial networks. Thus, we need to pay greater attention to the bacterial networks of aquatic macrophyte rhizospheres at different growth stages. In addition, our results showed that bacterial network of the P. australis zones had more nodes and higher modularity than the zones without P. australis. Higher modularity was shown beneficial for increasing the stability of microbial networks and the resistance profiles of microbial communities during environment changes [47].
With respect to topology, numerous studies have focused on putative keystone species that play distinct roles in co-occurrence networks [16,18]. We showed that the module hubs and connectors of the networks were completely different between the P. australis and non-P. australis zones, suggesting that P. australis roots can affect ecological functions and structuring bacterial networks in the riparian zone of Lake Taihu. Acidimicrobiales (Actinobacteria) and Gp3 (Acidobacteria) were module connectors in the bacterial network of the P. australis zones. Similarly, Actinobacteria and Acidobcteria have emerged as keystone taxa in blueberry rhizosphere soil [20]. Actinobacteria have the ability to produce diverse antibiotics [48], potentially allowing members of this phylum to take part in protecting plants from pathogenic bacteria. Our results showed that Sphingomonas (Alphaproteobacteria) was the module hub in the P. australis zones network, which could be correlated with signaling molecule regulation (e.g., quorum sensing) [49]. Previous studies showed that microorganisms can use quorum sensing (QS) to control competitive and cooperative behaviors in the plant rhizosphere [18,[50][51][52]. In contrast, the module hubs and connectors of bacterial network in the non-P. australis zones were assigned to Anaerolineaceae (Chloroflexi), Ignavibacterium (Ignavibacteriae), and Pelosinus (Firmicutes). These members have been reported to be involved in biogeochemical cycling, such as hydrocarbon degradation, sulfur cycling, and Fe (III) reduction [53][54][55]. In addition, we also observed that the role shifting of keystone species appeared between the two networks, which was the result of changes in keystone species as the surrounding environment changed [18].

Links among Aquatic Macrophytes, Local Factors, and Bacterial Communities in the Riparian Zone of Lake Taihu
We showed that taxonomic-based dissimilarities in bacterial communities were different in the four assayed zones, which exhibited significant correlations with local factors. Environmental heterogeneity has been suggested to be an important factor driving the distribution of bacterial communities in the riparian zones of lakes [56]. However, P. australis zones were shown to exhibit phylogenetic closeness of the bacterial communities under different hydrologic conditions, indicating the existence of common selective pressures shaping bacterial communities in the rhizosphere [20]. Moreover, the local factors had a lower correlation with bacterial community structures and networks of P. australis zones compared with non-P. australis zones (Tables 1 and S6). A number of mechanisms could explain why P. australis exerted a strong selection on riparian bacterial communities. Root exudates of the P. australis, such as oxalic, acetic, and succinic acid act as endogenous carbon sources for heterotrophic bacteria [11]. The properties of the soil/sediment close to roots are modified by the plant roots. For example, the roots of P. australis secretes oxygen to oxygenate the anaerobic sediments [57], which can affect the rhizosphere bacterial community composition and diversity. Furthermore, the assembly of rhizosphere bacterial community is based on functional activities related to the nutrient metabolism of plants [58,59]. Similarly, microbes with functional genes encoding transporters were shown to be crucial to the assembly of bacterial rhizosphere communities [60].
Our results showed that the bacterial community structure in the nonrhizosphere zones was correlated with NH 4 + , DOC, and OM, similar to the results of other studies in different freshwater ecosystems [6,61]. Notably, SWC and particle size correlated with the bacterial community and network structures in both the P. australis and non-P. australis zones. Although P. australis exhibited the ability to select and enrich riparian bacterial community composition and interactions, local factors can also directly or indirectly drive the variation in bacterial communities. Previous studies showed that dynamic riparian hydrology can influence the diffusion of oxygen and organic matter decomposition in the sediment, which can affect bacterial community structure and activity in riparian zones [6,62]. In the context of hydrologic variation, the anaerobic conditions of sediments can prompt the replacement of the rhizomes of P. australis by adventitious roots [2,63], and these changes can indirectly affect rhizosphere bacterial community composition and network structures. Furthermore, the existing literature has shown that differences in bacterial assemblages are associated with different particle sizes in various ecosystems [64][65][66]. Soil particle size can provide different surface properties and microhabitats that strongly contribute to the activity and composition of microbial communities [67,68]. For aquatic macrophytes, the riparian particle size can also affect their growth and root exudates, indirectly regulating the characteristics of rhizosphere bacterial communities.

Conclusions
Our results revealed that P. australis roots exerted effects on bacterial community composition and network structures in the riparian zone of Lake Taihu. Although the riparian bacterial communities exhibited differences based on taxonomic dissimilarity in the four assayed riparian zones, the P. australis zones enriched the phylogenetic closeness of the bacterial community. The results showed that two of the bacterial networks were nonrandom networks, and completely different topological properties were observed between the P. australis and non-P. australis zones. The putative keystone taxa of the bacterial network of the P. australis zones were affiliated with Sphingomonas (Alphaproteobacteria), Acidimicrobiales (Actinobacteria), and Gp3 (Acidobcteria), which might be involved in the regulation of bacterial interactions and plant growth. Furthermore, the hydrological regime and particle size of the riparian zone of Lake Taihu were correlated with the bacterial community and network structures. Overall, these results suggest that the distribution of P. australis and local factors are crucial for shaping bacterial communities and interactions in the riparian zone of Lake Taihu.