The epidemicity of facultative microsymbionts in faba bean rhizosphere soils

The epidemicity of facultative microsymbionts in faba bean rhizosphere soils

Soil Biology & Biochemistry 115 (2017) 243e252 Contents lists available at ScienceDirect Soil Biology & Biochemistry journal homepage: www.elsevier...

3MB Sizes 7 Downloads 23 Views

Soil Biology & Biochemistry 115 (2017) 243e252

Contents lists available at ScienceDirect

Soil Biology & Biochemistry journal homepage: www.elsevier.com/locate/soilbio

The epidemicity of facultative microsymbionts in faba bean rhizosphere soils Hui Yang Xiong a, Xing Xing Zhang a, Hui Juan Guo a, Yuan Yuan Ji a, Ying Li a, Xiao Lin Wang a, Wei Zhao a, Fei Yu Mo a, Jin Cheng Chen a, Tao Yang b, Xuxiao Zong b, Wen Xin Chen a, Chang Fu Tian a, * a State Key Laboratory of Agrobiotechnology, MOA Key Laboratory of Soil Microbiology, Rhizobium Research Center, College of Biological Sciences, China Agricultural University, 100193 Beijing, China b Institute of Crop Sciences/The National Key Facility for Crop Gene Resources and Genetic Improvement, Chinese Academy of Agricultural Sciences, 100081 Beijing, China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 28 March 2017 Received in revised form 26 August 2017 Accepted 29 August 2017 Available online 4 September 2017

The epidemicity of bacteria facultatively associated with eukaryotes, involves not only housekeeping genes but also genes linked with the pathogenesis and symbiosis. Here, by characterizing both housekeeping (rpoB and 16S rRNA gene) and nodulation (nodD) genes, we explore processes shaping epidemic patterns of facultative microsymbionts from rhizosphere soils of faba bean in three ecoregions. Although total bacterial communities in rhizosphere were not significantly differentiated across ecoregions, rpoB amplicon sequencing uncovered that Rhizobium laguerreae and Rhizobium anhuiense were predominant in different samples with contrasting pH or salt content. However, R. anhuiense can outcompete R. laguerreae in certain sterilized soils where R. laguerreae originally dominated, and viceversa. Contrasting bacterial taxa associated with either R. laguerreae or R. anhuiense in soils. The biogeographical pattern of nodD was more clear than that of rhizobial species in both rhizosphere soils and nodules. Competitive nodulation experiments demonstrated a hierarchical selection on nodD genotypes and their genomic backgrounds by faba bean cultivars. Taken together, abiotic and biotic factors in soils and the selection by legume hosts are either indirectly or directly involved in shaping rhizobial species-level taxonomic biogeography, which however cannot be used to infer spatial patterns of nodulation gene. © 2017 Elsevier Ltd. All rights reserved.

Keywords: Biogeography Epidemicity Legume Rhizobia Symbiosis Vicia faba

1. Introduction Understanding the epidemicity of mutualistic and pathogenic bacteria associated with eukaryotes is crucial for the further practical management of their positive or negative impacts on host fitness (Hussa and Goodrich-Blair, 2013; Martínez and Baquero, 2002; Nemergut et al., 2013). Specific interactions between these bacteria and their hosts are usually at the species or even intraspecific level as demonstrated using pure cultures isolated from infected host tissues. However our knowledge on the diversity of these bacteria in the environment is limited due to both the lacking of selective medium and the low resolution (above the genus level) of the 16S rRNA gene survey, which is widely used in the cultureindependent studies (Nemergut et al., 2013). Alternative

* Corresponding author. E-mail address: [email protected] (C.F. Tian). http://dx.doi.org/10.1016/j.soilbio.2017.08.032 0038-0717/© 2017 Elsevier Ltd. All rights reserved.

housekeeping genes such as rpoB encoding the beta subunit of RNA polymerase has been proposed as a good candidate due to its intraspecific resolution and its robustness as a taxonomic and phylogenetic marker (Adekambi et al., 2008; Drancourt et al., 2009; Vos et al., 2012). The rhizobium-legume symbiosis has been widely studied as one of the model mutualistic systems (more than 7000 papers in the past 5 years) and is an essential component of sustainable agriculture due to its high efficiency in nitrogen fixation. The symbiotic interaction is initiated through the specific sensing of host symbiotic signal by rhizobial NodD. The activated form of NodD further induces the transcription of nodulation genes involved in the biosynthesis of Nod factors that are specifically recognized by legume receptors before the intracellular infection of nodules by rhizobia (Remigi et al., 2016). These key nodulation genes including nodD are organized in a genomic island in either the symbiosis plasmid or the chromosome of different rhizobial

244

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

species and usually exhibit a distinct evolutionary history compared to the chromosomal background due to their horizontal transfer potential (Guo et al., 2014; Ling et al., 2016; Tian et al., 2010; Zhang et al., 2014). As typical facultative microsymbionts, rhizobia live saprophytically in soils and can enter the nitrogen-fixing symbiosis with legumes under appropriate conditions. It is hypothesized that rhizobia in soils experience a transient bottleneck mediated by the host selection of the most compatible strain during the symbiotic interactions. Indeed, by examining the nodule occupancy of signature-tagged mutants in an inoculant mixture, it has been demonstrated that most nodules are usually occupied by a single mutant clone (Pobigaylo et al., 2008; Wang et al., 2016). Competition by indigenous strains in rhizosphere soils can affect symbiotic interaction between an inoculant strain and the host plant, and this has been the major obstacle for improving the inoculation efficiency of rhizobia (Brockwell and Bottomley, 1995). However, rhizobial relative abundance in rhizosphere soils and factors affecting its variations have not been directly addressed. Faba bean (Vicia faba L.) was domesticated at least 7000 years ago in southwestern Asia (Tanno and Willcox, 2006; Zeder, 2009). It is now being cultivated either in the season of winter or spring around the world, with China and Europe as the main producers (Duc et al., 2010). Recent studies have revealed a considerable genetic differentiation in faba bean germplasms, generally corresponding to the geographical origin (Bao et al., 2006; Duc et al., 2010; Wang et al., 2012; Zong et al., 2009). Rhizobium leguminosarum bv. viciae (Rlv), R. fabae, R. laguerreae and R. anhuiense have been isolated from faba bean nodules in the field (Mutch and Young, 2004; Saïdi et al., 2014; Tian et al., 2008; Zhang et al., 2015), and these isolates can be collectively grouped into the symbiovar of viciae reflecting bacterial adaptation to legumes (Kumar et al., 2015; Remigi et al., 2016; Rogel et al., 2011). In this study, we aimed to investigate the epidemicity of rhizobia in rhizosphere soils of faba bean across three ecoregions in China and the potential factors shaping the epidemic patterns. To this end, the diversity of rhizobial genospecies and nodulation gene in rhizosphere soils was studied by using the high-throughput sequencing analyses of rpoB and nodD. Given the horizontal transfer ability of symbiotic functions across bacterial species, theoretically it is challenging to distinguish all rhizobial species compatible with faba bean from the huge pool of bacteria in soils. To determine a conservative list of faba bean microsymbionts in rhizosphere soils, we also studied the diversity of nodule isolates from the same faba bean plants, from the rhizosphere of which soil samples were collected. The geographical patterns of rpoB genospecies and nodD genotypes in rhizosphere soils and nodule isolates were analyzed, and abiotic or biotic factors potentially shaping these patterns were investigated. The significance of our findings was further discussed in the context of the ecology and evolution of bacteria, which form either pathogenic or mutualistic associations with eukaryotic hosts. 2. Material and methods 2.1. Sampling, DNA extraction and soil characterization Rhizosphere soils were collected from 44 sites in the northwest (Gansu and Qinghai provinces), east (Jiangxi, Anhui, Henan and Zhejiang), and southwest (Yunnan) growing regions of faba bean in China (Table S1 in Supporting Information). The distance between the sites ranges from 3.3 km to 2264.4 km. During the blooming stage of faba bean, soil attached to roots of 5e10 random plants within 100 m2 were pooled together after gentle shaking by hand, and defined as one “thick” rhizosphere soil sample (~500 g) as

described earlier (Sørensen et al., 2009; Thrane et al., 2000). DNA was extracted from 0.5 g samples using the Fast DNAtm SPIN Kit for soil (MP Biomedicals, Cleveland, OH, USA). Root nodules were collected from the same faba bean plants. Isolation of rhizobia and DNA extraction were done as described earlier (Tian et al., 2007). Soil chemical properties (Table S1) were analyzed, using the standard procedures described earlier (Lu, 2000), at the Plant Nutrient and Resource Research Institute, Beijing Academy of Agriculture and Forestry Sciences, China. Briefly, pH was determined for a suspension of 1:5 soil:water (w/v) using a pH meter; total N was measured by the Kjeldhl method; available phosphorus (P) was extracted by sodium bicarbonate (0.5 M, pH 8.5) and determined by Mo-Sb colorimetric method; available potassium (K) was extracted with ammonium acetate (1 M, pH 7.0) and determined using a flame photometer; organic carbon was measured by potassium dichromate oxidation and subsequent titration with ammonium ferrous sulfate (0.2 M), and organic matter (OM) was determined by multiplying the organic carbon content by 1.724 assuming that 58% of OM is carbon; concentrations of soil Kþ, Ca2þ, Naþ, Mg2þ, and SO2 4 extracted using deionized water (1:5 w/v) were determined by inductively coupled plasma atomic emission spectrometry; the level of HCO 3 in a suspension of 1:5 soil:water was measured using an automatic potentiometer and hydrochloric acid (20 mM); Cl content in the soil:water suspension (1:5 w/v) was titrated using silver nitrate (20 mM); electrical conductivity (EC) was measured in deionized water (1:5 w/v). To determine the total salt content, the aqueous extract from a suspension of 1:5 soil:water was evaporated in a water bath; organic matter within the resultant product was then removed by adding 15% H2O2 and by heating in a water bath; and the salt content was finally determined by gravimetry of the residue evaporated in a drying oven. 2.2. DNA fingerprinting and analyses of rpoB and nodD for nodule isolates Total DNA from 1456 purified nodule isolates was used as the template for repetitive extragenic palindromic PCR (REP PCR) with the BOX-A1R primer as described earlier (Koeuth et al., 1995; Tian et al., 2007). Three hundred ninety-seven representative strains of 85 BOX fingerprint types were further characterized by amplification and sequencing of rpoB using primers rpoB-F/rpoB-R (TCGCAGTTCATGGACCAGG/GTAGCCGTTCCAGGGCATG) as described previously (Zhang et al., 2017). Since all published genomes and nod genes of rhizobia belonging to the symbiovar viciae have a conserved single copy of nodD (Kumar et al., 2015; MarekKozaczuk et al., 2013; Tian et al., 2010; Young et al., 2006), the primers nodDviciaeF88/Y6 (TGCAGAGACGGGAGCTARTTC/CGCA WCCANATRTTYCCNGGRTC) were used to amplify nodD gene of all  ze  nodule isolates as described earlier (Macdonald et al., 2011; Ze et al., 2001). The nodD PCR products were analyzed by restriction fragment length polymorphism (RFLP) analysis using endonucleases HaeIII, MseI and BstNI. The restriction fragments were separated by electrophoresis and visualized under UV light as previously described (Tian et al., 2007). The nodD PCR products of 366 representative strains were purified and sequenced on an ABI 3730xl sequencer. 2.3. High-throughput sequencing of rpoB, 16S rRNA and nodD from rhizosphere soils The rpoB, 16S rRNA gene and nodD were amplified using the primers rpoB1479-F/rpoB1831-R (rpoB) (Zhang et al., 2017), 515f/ 806r (the V4/V5 region of the 16S rRNA gene) (Caporaso et al., 2012) and nodDviciaeF88/nodDviciaeR443 (nodD) (Macdonald et al., 2011), respectively. The nodD fragment obtained using this pair of

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

primers is within the larger fragment amplified using nodDviciaeF88/Y6, which were used in the analysis of nodule isolates as described above. The first round of PCR amplifications was conducted with the primers containing the transposase sequences 50 TCG TCG GCA GCG TCA GAT GTG TAT AAG AGA CAG-30 and 50 -GTC TCG TGG GCT CGG AGA TGT GTA TAA GAG ACA G-30 at the 50 end of the forward and reverse primers, respectively. PCR products were purified with AMPure XP Beads for the second round of amplification with P5 þ index i5 (50 -AAT GAT ACG GCG ACC ACC GAG ATC TAC AC [i5] TCG TCG GCA GCG TC-30 ) and P7 þ index i7 (50 -CAA GCA GAA GAC GGC ATA CGA GAT [i7] GTC TCG TGG GCT CGG-30 ) primers. The PCR products with P5/P7 adapters and double indexes were purified with double volume of AMPure XP Beads, and quantified using Qubit 2.0 Fluorometer and Qubit dsDNA HS assay kit (Invitrogen, United States of America). Multiple samples were pooled together in equimolar concentrations. Then 10 ng of pooled PCR products were sequenced on the Illumina MISeq platform at the National Human Genome Center of China at Shanghai, China. 2.4. Processing of raw data from high-throughput sequencing The MOTHUR v.1.35.1 software was used to process raw data produced by MISeq according to MISeq SOP (Schloss et al., 2009), using the corresponding reference databases (Table S2). Briefly, reads were combined and screened for reducing sequencing errors. The reads were aligned according to MISeq SOP (Schloss et al., 2009) based on the reference for 16S rRNA gene (the SILVA-based alignment release 123), rpoB (rpoB.reference.fas; Appendix S1 in Supporting information), and nodD (nodD.reference.fas; Appendix S2). The aligned reads were de-noised using the PyroNoise method (Quince et al., 2011) and chimeras were removed using the UCHIME algorithm (Edgar et al., 2011). OTUs for 16S rRNA gene and rpoB were defined at the 97% and 97.7% (Zhang et al., 2017) sequence similarity levels, respectively. The taxonomic assignments for rpoB OTUs were determined according to the taxonomic information of 113 rhizobial rpoB sequences in rpoB.reference.fas (Appendix S1). The genotypes of nodD sequences were determined based on nodD.genotypes.fas (Appendix S3), which harbors five nodD genotypes of rhizobia isolated from faba bean nodules.

245

microsymbionts could be potentially linked to the distribution patterns of other bacteria in rhizosphere soils, the relative abundance of bacterial genera defined using 16S rRNA gene was regressed on the PCoA ordination diagram of Rhizobium species using the vegan package. Alpha-diversity was measured by calculating Pielou's evenness, Shannon index and Simpson's index of diversity with the vegan package, and Faith's phylogenetic diversity with the MOTHUR software (Schloss et al., 2009). The neighborjoining phylogenetic tree of nodD was constructed using MEGA6 (Tamura et al., 2013) and representative sequences obtained in this study and reference sequences described earlier (Mutch and Young, 2004; Tian et al., 2010). 2.6. Competition in soils and nodulation To compare the survival ability of R. anhuiense and R. laguerreae in soils, nine samples of rhizosphere soils exhibiting a linear gradient of pH were sterilized by gamma rays with a minimum absorption of 15.0 kGy. The competition experiment was initiated by equal inoculation of two strains at the level of 5  104 cells per gram soils. Three replicates were performed and inoculated soils were kept at 28  C. The colony-forming units (CFUs) were recorded every three days until reaching a plateau. All colonies presented at a same dilution level in three replicates (ranging from 29 to 98 colonies in different competition experiments) were characterized using colony PCR with primers rpoB-F/rpoB-R and subsequent RFLP analysis of PCR products with endonuclease AluI and KpnI. The nodule occupancy of Rhizobium strains was tested on two modern cultivars of faba bean (Qidou2 and Fengdou1). Faba bean seeds were rinsed in 95% ethanol for 1 min and then surface sterilized in 15% (vol/vol) NaClO solution. Rhizobium strains were mixed in an equal quantity (OD600, 0.2; 1 ml volume) before inoculating germinated seeds. At 30 days post inoculation, nodules were surface sterilized and rhizobia inside were isolated as described previously (Tian et al., 2007). Colony PCR with either rpoB-F/rpoB-R or nodF88/Y6 as primers was used to amplify rpoB and nodD fragments, which were further subject to RFLP analysis as described above to distinguish different test strains. 2.7. Accession numbers

2.5. Statistical and phylogenetic methods All data manipulation was performed in the vegan package (Oksanen et al., 2015) of R statistic environment (R Development Core Team, http://www.R-project.org) unless otherwise indicated. Rare OTUs containing less than 0.01% of total reads of 16S rRNA or rpoB were removed before further analyses to avoid inflated estimates of diversity due to the inclusion of erroneous reads (Bokulich et al., 2013). The site-to-site variability in OTU composition was measured by Bray-Curtis dissimilarity. The dissimilarities were visualized by either principal coordinates analysis (PCoA) or the geospatial biodiversity analysis with GENGIS 2.0.0 (Parks et al., 2009). The input tree used in GENGIS was constructed by the unweighted pair group method with arithmetic mean (UPGMA) of Bray-Curtis dissimilarity. Monte Carlo permutation test (1000 iterations) was used to test whether the fit of ordered leaf nodes to geographical locations in the GENGIS result is significantly better than expected by chance alone. Edaphic factors with strong cocorrelation (Spearman coefficient above 0.70, P-value < 0.05) were removed before the forward selection (P-value < 0.001) of those factors explaining significant amounts of variation. The relationship between significant edaphic factors and the BrayCurtis dissimilarity was defined by the constrained analysis of principle coordinates (CAP). To investigate whether the betadiversity of Rhizobium species identified as faba bean

Raw data reported in this study are publicly deposited in the Sequence Read Archive database (accession no. SAMN06198893SAMN06198936 in the BioProject PRJNA359745). The sequences of rpoB (KY424542-KY424951) and nodD (KY424952-KY425351) are deposited in GenBank. 3. Results 3.1. Biased geographical distribution of facultative microsymbionts of faba bean The rpoB high-throughput sequencing produced 281,847 reads after quality trimming for 44 soil samples across three faba beangrowing ecoregions (Tables S1 and S2). At the 97.7% identity of rpoB sequences, there were 3265 OTUs in total belonging to: Rhizobium (497 OTUs), Sinorhizobium (186 OTUs), Mesorhizobium (1107 OTUs) and Bradyrhizobium (1475 OTUs). 1606 OTUs containing less than 0.01% of total reads were removed before further analyses. Based on Bray-Curtis dissimilarities, the east and southwest samples were highly intermingled with each other and the whole dataset did not exhibit a statistically significant biogeographical pattern (Fig. S1a, Monte Carlo permutation test of geographical pattern, P-value ¼ 0.488), although the northwest samples formed a distinct cluster. Forty-four OTUs were identified

246

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

Fig. 1. Beta-diversity of Rhizobium from rhizosphere and nodules of faba bean. (a) Geographic distribution of Rhizobium species isolated from nodules (left; based on rpoB genospecies assignment of each BOX fingerprint group) or identified in rhizosphere of faba bean (right; based on rpoB amplicon sequencing). (b) Bray-Curtis dissimilarity dendrogram of Rhizobium species in nodules (left) and rhizosphere (right). P-values indicate the significant level at which the fit of ordered leaf nodes to geographical locations is better than expected by chance alone (Monte Carlo permutation test, 1000 iterations). Branches are colored according to the corresponding regions, and branches encompassing two or more regions are indicated by white color. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

as documented legume microsymbionts (Fig. S1b). A dendrogram was further constructed based on Bray-Curtis dissimilarities of this subset of rhizobial genospecies and showed a highly intermingled beta-diversity pattern across three ecoregions (Fig. S1b; Monte Carlo permutation test of geography, P-value ¼ 0.768). To determine which genospecies in rhizosphere soils could be assigned as faba bean microsymbionts, we characterized the taxonomic information of 1456 nodule isolates of the same plants, from the rhizosphere of which soil samples were collected. There were 391, 600, and 465 isolates from the east, northwest, and southwest regions, respectively. These isolates belonged to 85 BOX fingerprint groups (Table S3). The rpoB sequence analysis of representative strains identified R. anhuiense, R. laguerreae, Rlv, R. pisi/R. fabae, R. sp. I and R. sp. II among these nodule isolates (Fig. 1a). R. pisi and R. fabae could not be differentiated herein using rpoB sequences. This is in line with earlier evidence that they could be the same species (Tian et al., 2010; Zhang et al., 2012), although the type strains of R. fabae and R. pisi harbor distinct nodD genotypes reflecting their original hosts faba bean and Pisum sativum L. (as discussed below). These six genospecies were all identified in rhizosphere samples (Fig. 1a). The relative abundance values of these Rhizobium genospecies were similar between rhizosphere and nodules from the northwest and east regions (Fig. 1a). The southwest samples harbored more genospecies than those from the other two regions. In the southwest region, different genospecies dominated in rhizosphere (R. sp. I followed by R. anhuiense) and nodules (R. laguerreae followed by R. anhuiense) (Fig. 1a). Dendrograms based on Bray-Curtis dissimilarities of these Rhizobium genospecies in rhizosphere and nodules (Fig. 1b) suggested that a significant biogeographical pattern existed for nodule isolates (P < 0.01) but not for rhizosphere samples

(P ¼ 0.105). For nodule isolates, east and northwest regions preserved two distinct pools of faba bean microsymbionts generally homogeneous within each ecoregions whereas these two pools were all intermingled with that of the southwest region (Fig. 1b). 3.2. Environmental factors shaping geographical distribution of R. anhuiense and R. laguerreae R. anhuiense and R. laguerreae were the major faba bean rhizobia present in rhizosphere in three ecoregions and dominated in the east and northwest regions, respectively (Fig. 2a). The CAP analysis based on Bray-Curtis dissimilarities (Fig. 2b) revealed that the abundance of R. laguerreae correlated positively with pH, total potassium and salt content in soils, while R. anhuiense abundance correlated negatively with these edaphic factors. To test whether contrasting edaphic conditions would directly lead to different survival fates of R. anhuiense and R. laguerreae, nine samples of rhizosphere soils exhibiting a linear gradient of pH were chosen (Fig. 2c). Representative strains of R. anhuiense and R. laguerreae were equally inoculated into the sterilized soils, and the relative abundance of two species was examined 21 days post inoculation. R. anhuiense survived in all samples and showed a better competitiveness in samples Y14, Y6, Y1, Y3, S7 and S12. Viable R. laguerreae cells were found in all soil samples except Y3 and R. laguerreae dominated in A2. Notably Y1 and S12 with the highest level of salt content and pH, respectively, were originally dominated by R. laguerreae, and A2 was dominated by R. anhuiense in the field. Therefore, it is possible that other undetermined abiotic and biotic factors might also affect the relative abundance of two Rhizobium species in soils. The beta-diversity of bacterial community across 44 rhizosphere

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

247

Fig. 2. Biased distribution of Rhizobium species from rhizosphere is not just due to survival ability to edaphic conditions. (a) Ternary plot of the relative abundance of Rhizobium species from rhizosphere of faba bean (based on rpoB amplicon sequencing). (b) Constrained analysis of principal coordinates for samples from rhizosphere of faba bean. Samples are colored by their corresponding regions (Blue: east; red, northwest; green, southwest). (c) pH and salt content of representative rhizosphere soils. (d) Relative abundance of Rhizobium species in intact rhizosphere soils (upper sector) and in the competition experiment using sterilized rhizosphere soils (bottom). R. anhuiense and R. laguerreae inoculated in a ratio of 1:1; data are mean ± s.e.m. based on three replicates. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

samples was analyzed using 16S rRNA gene amplicon sequencing. A total of 289,277 reads were obtained after quality trimming (Table S2). At the 97% identity of 16S rRNA gene sequences, 33,201 OTUs were identified. 1239 OTUs containing more than 0.01% of total reads were kept for the analyses. No significant geographical patterns were observed across ecoregions (Fig. 3a and Fig. S2; Pvalues > 0.05) indicating a global intermingled bacterial community in faba bean rhizosphere. However, further regression analysis uncovered that the relative abundances of 67 genera significantly correlated to the beta-diversity of R. anhuiense and R. laguerreae (Fig. 3b and Table S4). Noteworthy, there were more bacterial genera whose abundances correlated positively with the abundance of R. anhuiense than with that of R. laguerreae (Fig. 3b and Table S4). This is in line with the significant higher values of Shannon index, Pielou's evenness and Simpons's index of diversity for rhizosphere samples with R. anhuiense as the predominant microsymbionts of faba bean than those samples with the abundant R. laguerreae (Fig. 3c; T-test, P-values < 0.01), though no significant difference was observed in their bacterial phylogenetic diversity (T-test, P-value > 0.05).

3.3. Strong biogeographical patterns of symbiosis gene nodD in both nodules and rhizosphere Five genotypes of nodD were found among nodule isolates using RFLP and they belonged to five phylogenetic classes based on nodD sequences from 366 representative strains (Fig. 4a). The nodD amplicon sequencing of rhizosphere samples produced totally 417,246 reads (96.4% of total reads after quality trimming) belonging to these five nodD genotypes (Fig. S3). The relative abundances of five nodD genotypes were similar between nodule isolates and rhizosphere samples (Fig. S3). These nodD sequences exhibited a significant geographical distribution pattern across three ecoregions in both nodules and rhizosphere (Fig. 4b and c; P-values < 0.05). The east regionwas dominated by nod-1 and nod-2, the southwest region was dominated by nod-4 and nod-2, and the northwest region harbored nod-1 and nod-3 as the main nodD genotypes (Fig. 4d and Fig. S3). The dominant R. anhuiense and R. laguerreae can harbor all nodD genotypes (Fig. 4d), suggesting that different nodD genotypes may help these two species to adapt to host plants of diverse genotypes, which remain largely unexplored.

248

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

Fig. 3. R. anhuiense and R. laguerreae are accompanied by contrasting bacterial communities in rhizosphere. (a) Bray-Curtis dissimilarity dendrogram of bacterial communities based on genera determined using 16S rRNA gene amplicon sequencing of rhizosphere samples. P-value indicates the significant level at which the fit of ordered leaf nodes to geographical locations is better than expected by chance alone (Monte Carlo permutation test, 1000 iterations). Branches are colored according to the corresponding regions, and branches encompassing two or more regions are indicated by white color. (b) Regressed bacterial abundance (arrows) that significantly (P-value < 0.05) correlated with the betadiversity (Bray-Curtis dissimilarity of data obtained using rpoB amplicon sequencing of rhizosphere samples) of R. anhuiense and R. laguerreae presented in a PCoA plot. Names for representative bacterial taxa are shown. Samples are colored by their corresponding regions (Blue: south; red, northwest; green, southwest). (c) Alpha-diversity of bacteria in rhizosphere samples dominated either by R. anhuiense or R. laguerreae. Data are mean ± s.e.m.; ** (P-value < 0.01, T-test); NS, nonsignificant. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

3.4. Selection of symbiotic partners by faba bean cultivars The presence of dominant Rhizobium species and nodD genotypes in each ecoregion implied a potential selection of symbiotic partners by faba bean cultivars. To test this hypothesis, four Rhizobium species (R. laguerreae, R. sp. I, Rlv and R. anhuiense) harboring the same nodD genotype nod-4 were equally inoculated on two faba bean cultivars (Qidou2 and Fengdou1). Rlv and R. laguerreae occupied more nodules than the other two species on Qidou2 (Fig. 5a). Rlv was preferred by Fengdou1 and showed significantly higher nodule occupancy than R. anhuiense and R. sp. I (Fig. 5b). When Qidou2 was inoculated with equal amounts of five R. laguerreae strains harboring different nodD genotypes (Fig. 5c), R. laguerreae with nod-2 was the major occupant of nodules. However, no significant individual dominant strain in Qidou2 nodules was found among R. anhuiense strains containing different nodD genotypes (Fig. 5d). R. laguerreae strain with nod-2 was also the major microsymbiont of Fengdou1 nodules (Fig. 5e). R. anhuiense strains with nod-3 and nod-5 occupied 60% and 39% of Fengdou1 nodules, respectively (Fig. 5f). When ten strains of R. anhuiense and R. laguerreae were equally inoculated on Qidou2 (Fig. 5g), R. laguerreae with nod-5 was the most abundant in nodules, though R. laguerreae harboring nod-2/nod-3 and R. anhuiense

harboring nod-2 also occupied a considerable ratio of nodules. Among the same ten strains, Fengdou1 preferred R. anhuiense with nod-3 (Fig. 5h). 4. Discussion Similar to diverse bacterial pathogens associated with eukaryotes (Moran et al., 2008), rhizobia have a typical facultative symbiotic life cycle. To effectively manage the competitive nodulation by indigenous rhizobia, it is particular important to reveal the epidemic pattern of these facultative symbionts both within host tissues and in the surrounding environment. However, the 16S rRNA gene survey of low taxonomic resolution is unable to address these species- and strain-level questions (Nemergut et al., 2013). Using high-throughput sequencing analyses of 16S rRNA gene, rpoB and a key symbiosis gene nodD, this study not only unveiled the pattern and processes of the biogeography of rhizobial species but also highlighted an independent epidemic pattern of the symbiosis gene. The general theory of community ecology outlines that speciation and dispersal can add species into communities, and the relative abundances of these species are then structured by drift and selection (Vellend, 2010). Faba bean was domesticated in

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

249

Fig. 4. Beta-diversity of rhizobial nodD from rhizosphere and nodules of faba bean. (a) Neighbor-joining phylogenetic tree of NodD. Bootstrap values above 50% are shown. Representative sequences of five nodD genotypes obtained from nodule isolates in this study are colored by red. (b-c) Bray-Curtis dissimilarity dendrogram of Rhizobium nodD from nodule isolates (b) and nodD amplicon sequencing of rhizosphere samples (c). P-values indicate the significant level at which the fit of ordered leaf nodes to geographical locations is better than expected by chance alone (Monte Carlo permutation test, 1000 iterations). Branches are colored according to the corresponding regions, and branches encompassing two or more regions are indicated by white color. (d) Relative abundance of nodD genotypes identified among nodule isolates in three ecoregions. Rhizobium species harboring these nodD genotypes are indicated. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

southwest Asia and then spread to all inhabited continents (Duc et al., 2010) and it has been reported that legume seeds naturally rez-Ramírez et al., 1998). Bacteria carry rhizobia on their testa (Pe could be cotransferred with faba bean seeds. Indeed, despite the long geographical distances, the total bacterial communities (OTUs defined by 16S rRNA gene) or rhizobial communities (OTUs defined by rpoB) in faba bean rhizosphere can be similar across ecoregions. The same BOX fingerprint groups including those of low abundance can be isolated from faba bean nodules in different ecoregions. Therefore, the contrasting distribution patterns of R. laguerreae and R. anhuiense should not be due to the dispersal limitation. Faba bean was introduced to the northwest region of China from southwest Asia 2100 years ago via traders along the Silk Road (Lang and Ying, 1997). Noteworthy, R. laguerreae, the predominant microsymbiont of faba bean in the northwest region, has been reported as effective microsymbiont of faba bean in Europe, Africa  and America (Alvarez-Martínez et al., 2009; Saïdi et al., 2014; Santillana et al., 2008). However, in the east region R. anhuiense was the dominant partner of faba bean. Faba bean belongs to the cross-nodulation group containing species of Vicia, Pisum, Lens and Lathyrus (Mutch and Young, 2004), indicating that introduced faba bean, may also recruit R. anhuiense from other wild and cultivated

legumes in this region. The abundance of R. laguerreae correlated positively with soil pH and salt content whereas R. anhuiense correlated negatively with these two edaphic factors. Noteworthy, the statistical correlation between the abundance of certain bacterial phyla and edaphic conditions has been widely reported (Hanson et al., 2012), but rarely experimentally tested. In this study, the competitive ability of R. anhuiense and R. laguerreae was further investigated. Unexpectedly, both species can survive and even dominate in sterilized soil samples, which were originally dominated by the other one of these two species. Recently R. anhuiense was reported as the major microsymbiont of a wild Lathyrus species on a seashore line with typical saline-alkaline edaphic conditions (Li et al., 2016). Thus pH and salt content are not the direct cause of the contrasting distribution patterns of R. laguerreae and R. anhuiense. In addition to abiotic factors, the fitness in one species could also depend either directly or indirectly on the density of other species (Vellend, 2010). Bacteria such as those belonging to Acidobacteria Gp6 and Gp3 identified using 16S rRNA amplicon-sequencing, exhibited significant correlations with the beta-diversity pattern of R. laguerreae and R. anhuiense in rhizosphere samples. Acidobacteria are ubiquitous and among the most abundant phyla in soil

250

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

Fig. 5. Selection of symbiotic partners by faba bean cultivars. (a-b) Nodule occupancy on faba bean cultivars Qidou2 (a) and Fengdou1 (b) by four Rhizobium species harboring the same nodD genotype (nod-4). (cef) Nodule occupancy on Qidou2 (ced) and Fengdou1 (eef) by five R. laguerreae strains (c and e) or five R. anhuiense strains (d and f) harboring different nodD genotypes. (g-h) Nodule occupancy on Qidou2 (g) and Fengdou1 (h) by ten Rhizobium strains harboring different nodD genotypes. All strains were inoculated in equal amounts. Data are mean ± s.e.m. based on at least from 100 nodules of three replicated experiments. Statistics: Kruskal-Wallis test, different letters indicate significant differences (False discovery rate < 0.05).

and different subgroups of Acidobacteria respond differentially to pH (Kielak et al., 2009). The abundance of Acidobacteria Gp6 correlated positively with soil pH and the abundance of R. laguerreae, while Acidobacteria Gp3 correlated negatively with pH but positively with the abundance of R. anhuiense. However, Acidobacteria members have been difficult to culture and the adaptation mechanisms under different pH conditions remain unknown. Noteworthy, Bradyrhizobium correlated positively with the abundance of R. anhuiense. Bradyrhizobium genomes are among the largest ones in the bacterial world and are enriched in lipid and secondary metabolisms (Tian et al., 2012). Further experiments should be performed to investigate mechanisms underlying these positive or negative relationships between Rhizobium species and other bacteria. Similar to pathogenesis genes of many pathogenic bacteria, key nodulation genes of rhizobia are located in a genomic island named symbiosis island, which can be found in either the chromosome or the symbiosis plasmid of different species (Remigi et al., 2016). With their specific role in symbiotic interactions with legumes, these nodulation genes have been strongly selected and can be horizontally transferred between species nodulating the same legume (Guo et al., 2014; Ling et al., 2016; Sullivan et al., 1995; Sullivan and Ronson, 1998; Young et al., 2006; Zhang et al., 2014). The nodD amplicon sequencing uncovered a geographical pattern of nodD genotypes stronger than that of rpoB genospecies,

indicating a potential selection of rhizobial nodD genotypes by faba bean plants. Alternatively, the geographical pattern of nodD could be governed by a stochastic process during the dispersal history of rhizobia associated with faba bean seeds, and we would expect similar nodule occupancy level between different nodD genotypes. However, competitive nodulation experiments demonstrated that different strains of the same rpoB genospecies harboring distinct nodD genotypes could be differentially selected by faba bean. Considering the significant genetic differentiations of faba bean landraces from different ecoregions (Bao et al., 2006; Duc et al., 2010; Wang et al., 2012; Zong et al., 2009), we hypothesized that host germplasms have contributed to the beta-diversity patterns of nodD genotypes. The bottleneck effect on the rhizobial populations during the infection process (Pobigaylo et al., 2008; Wang et al., 2016) may further enhance the geographical patterns. European faba beans preferred a distinct dominant nodD genotype, corresponding to nod-3 in this study (Mutch and Young, 2004; Tian et al., 2010). The nod-3 genotype together with nod-1 and nod-5 were also present in the faba bean isolates from Jordan, which is near the domestication center of faba bean (Duc et al., 2010; Tian et al., 2010). The dominant nodD genotypes nod-1, -2 and -3 in Chinese samples belong to a monophyletic clade containing sequences from faba bean isolates around the world, such as the type strains of R. fabae, R. laguerreae, and R. anhuiense (Saïdi et al., 2014; Tian et al., 2008; Zhang et al., 2015). Similarly, the

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

nod-5 genotype corresponds to another distinct but smaller monophyletic clade containing only faba bean isolates (Tian et al., 2010). By contrast, the nod-4 genotype belongs to a clade containing sequences from nodule isolates of diverse legumes within the cross-nodulation group of Vicia-Pisum-Lens-Lathyrus (Mutch and Young, 2004; Tian et al., 2010), such as the type strains of Rlv and R. pisi isolated from P. sativum (Ramírez-Bahena et al., 2008). The nod-4 genotype, regardless of which species harbored it, was less competitive than those characteristic nodD genotypes associated with faba bean rhizobia (the clade containing nod-1, -2 and -3; and the clade of nod-5) in the competitive nodulation experiment using two test cultivars. The relatively higher abundance of nod-4 in the southwest region than in the other two regions could be due to distinct germplasms of faba bean landraces in this region (Bao et al., 2006; Duc et al., 2010; Zong et al., 2009). Different species harboring the same nodD genotype can be also subject to differential selection by faba beans, for example R. anhuiense carrying nod-4 was less competitive in nodulation compared to R. laguerreae, while R. anhuiense with nod-3 was preferred by the faba bean cultivar Fengdou1 in the presence of R. laguerreae carrying nod-3. These findings imply a potential contribution of lineage-specific accessory genes or the involvement of diversified core genes in symbiosis. Although the integration of key symbiosis genes with genomic background is still largely unexplored, the importance of lineage-specific genes in symbiotic adaptation has been highlighted by recent studies of rhizobial comparative genomics and experimental evolution (Amadou et al., 2008; Capela et al., 2017; Remigi et al., 2016; Tian et al., 2012). The epidemic patterns of both housekeeping and symbiosis genes, as revealed in this study, could allow us to identify inoculant candidates, which showed better adaptive ability to legume cultivars within an ecoregion. These abundant elite strains might have been selected due to their better integration of key symbiosis genes with genome background.

5. Conclusions The epidemicity of bacteria, which are facultatively associated with eukaryotic hosts, involves not only housekeeping genes but also accessory functions directly participating in the pathogenesis or symbiosis. Using the system of Rhizobium-faba bean symbiosis, we revealed that different dominant Rhizobium species can be found from both nodules and rhizosphere soils in three ecoregions, but bacterial communities, including rhizobia, in rhizosphere soils were not significantly differentiated across ecoregions. However, nodD genotypes in rhizosphere soils exhibited a significant biogeographical pattern, which was similar to that of nodD genotypes of nodule isolates. In line with the hierarchical selection of nodD genotypes and their genomic backgrounds by faba bean cultivars, the taxonomic biogeography of nodule isolates cannot be directly used to infer the stronger spatial patterns of nodulation gene. Taken together, everything can be everywhere, but selection regimes particularly those acting on functions essential for interacting with hosts have shed substantial effects on the biogeography, ecology and evolution of facultative microsymbionts, which undergo a population bottleneck during the infection process.

Conflict of interest We declare that we do not have any commercial or associative interest that connected with the work submitted.

251

Acknowledgements We thank two anonymous reviewers for their valuable comments. This work was supported by Program for New Century Excellent Talents in University (NCET-13-0561), National Basic Research Program of China (973 program 2015CB158300) and Chinese Universities Scientific Fund (2017SY003). Appendix A. Supplementary data Supplementary data related to this article can be found at http:// dx.doi.org/10.1016/j.soilbio.2017.08.032. References Adekambi, T., Shinnick, T.M., Raoult, D., Drancourt, M., 2008. Complete rpoB gene sequencing as a suitable supplement to DNA-DNA hybridization for bacterial species and genus delineation. International Journal of Systematic and Evolutionary Microbiology 58, 1807e1814. http://dx.doi.org/10.1099/ijs.0.65440-0.   Ramírez-Bahena, M.H., García-Fraile, P., Alvarez-Martínez, E.R., Valverde, A., ~ iga, D., Peix, A., Vel Tejedor, C., Mateos, P.F., Santillana, N., Zún azquez, E., 2009. The analysis of core and symbiotic genes of rhizobia nodulating Vicia from different continents reveals their common phylogenetic origin and suggests the distribution of Rhizobium leguminosarum strains together with Vicia seeds. Archives of Microbiology 191, 659e668. http://dx.doi.org/10.1007/s00203-0090495-6. Amadou, C., Pascal, G., Mangenot, S., Glew, M., Bontemps, C., Capela, D., Carrere, S., Cruveiller, S., Dossat, C., Lajus, A., Marchetti, M., Poinsot, V., Rouy, Z., Servin, B., Saad, M., Schenowitz, C., Barbe, V., Batut, J., Medigue, C., Masson-Boivin, C., re, S., Me digue, C., 2008. Genome sequence of the beta-rhizobium Carre Cupriavidus taiwanensis and comparative genomics of rhizobia. Genome Research 18, 1472e1483. http://dx.doi.org/10.1101/gr.076448.108. Bao, S., Wang, L., Lu, M., He, Y., 2006. Analysis for genetic diversity of landraces of Vicia faba in Yunnan, China. In: Avila, C.M., Cubero, J.I., Moreno, M.T., Suso, M.J., Torres, A.M. (Eds.), International Workshop on Faba Bean Breeding and Agronomy. Junta de Andalucia, Cordoba, Spain, pp. 196e200. Bokulich, N.A., Subramanian, S., Faith, J.J., Gevers, D., Gordon, J.I., Knight, R., Mills, D.A., Caporaso, J.G., 2013. Quality-filtering vastly improves diversity estimates from Illumina amplicon sequencing. Nature Methods 10, 57e59. http:// dx.doi.org/10.1038/nmeth.2276. Brockwell, J., Bottomley, P.J., 1995. Recent advances in inoculant technology and prospects for the future. Soil Biology & Biochemistry 27, 683e697. http:// dx.doi.org/10.1016/0038-0717(95)98649-9. rissi, C., Perrier, A., Guetta, D., Gris, C., Valls, M., Capela, D., Marchetti, M., Cle Jauneau, A., Cruveiller, S., Rocha, E.P.C., Masson-Boivin, C., 2017. Recruitment of a lineage-specific virulence regulatory pathway promotes intracellular infection by a plant pathogen experimentally evolved into a legume symbiont. Molecular Biology and Evolution 3. http://dx.doi.org/10.1093/molbev/msx165. Caporaso, J.G., Lauber, C.L., Walters, W.a., Berg-Lyons, D., Huntley, J., Fierer, N., Owens, S.M., Betley, J., Fraser, L., Bauer, M., Gormley, N., Gilbert, J.A., Smith, G., Knight, R., 2012. Ultra-high-throughput microbial community analysis on the Illumina HiSeq and MiSeq platforms. The ISME Journal 6, 1621e1624. http:// dx.doi.org/10.1038/ismej.2012.8. kambi, T., 2009. The rpoB gene as a tool for clinical Drancourt, M., Raoult, D., Ade microbiologists. Trends in Microbiology 17, 37e45. http://dx.doi.org/10.1016/ j.tim.2008.09.008. Duc, G., Bao, S., Baum, M., Redden, B., Sadiki, M., Suso, M.J., Vishniakova, M., Zong, X., 2010. Diversity maintenance and use of Vicia faba L. genetic resources. Field Crops Research 115, 270e278. http://dx.doi.org/10.1016/j.fcr.2008.10.003. Edgar, R.C., Haas, B.J., Clemente, J.C., Quince, C., Knight, R., 2011. UCHIME improves sensitivity and speed of chimera detection. Bioinformatics 27, 2194e2200. http://dx.doi.org/10.1093/bioinformatics/btr381. Guo, H.J., Wang, E.T., Zhang, X.X., Li, Q.Q., Zhang, Y.M., Tian, C.F., Chen, W.X., 2014. Replicon-dependent differentiation of symbiosis-related genes in Sinorhizobium nodulating Glycine max. Applied and Environmental Microbiology 80, 1245e1255. http://dx.doi.org/10.1128/AEM.03037-13. Hanson, C.A., Fuhrman, J.A., Horner-Devine, M.C., Martiny, J.B.H., 2012. Beyond biogeographic patterns: processes shaping the microbial landscape. Nature Reviews Microbiology 10, 497e506. http://dx.doi.org/10.1038/nrmicro2795. Hussa, E.A., Goodrich-Blair, H., 2013. It takes a village: ecological and fitness impacts of multipartite mutualism. Annual Review of Microbiology 67, 161e178. http:// dx.doi.org/10.1146/annurev-micro-092412-155723. Kielak, A., Pijl, A.S., van Veen, J.A., Kowalchuk, G.A., 2009. Phylogenetic diversity of Acidobacteria in a former agricultural soil. The ISME Journal 3, 378e382. http:// dx.doi.org/10.1038/ismej.2008.113. Koeuth, T., Versalovic, J., Lupski, J.R., 1995. Differential subsequence conservation of interspersed repetitive Streptococcus pneumoniae box elements in diverse bacteria. Genome Research 5, 408e418. http://dx.doi.org/10.1101/gr.5.4.408. Kumar, N., Lad, G., Giuntini, E., Kaye, M.E., Udomwong, P., Shamsani, N.J., Young, J.P.W., Bailly, X., 2015. Bacterial genospecies that are not ecologically

252

H.Y. Xiong et al. / Soil Biology & Biochemistry 115 (2017) 243e252

coherent: population genomics of Rhizobium leguminosarum. Open Biology 5, 140133. http://dx.doi.org/10.1098/rsob.140133. Lang, L., Ying, H., 1997. Faba bean. In: Zheng, Z., Wang, S., Zong, X. (Eds.), Food Legume Crops in China. China Agricultural Press, Beijing, China, pp. 53e92 (in Chinese). Li, Y., Wang, E.T., Liu, Y., Li, X., Yu, B., Ren, C., Liu, W., Li, Y., Xie, Z., 2016. Rhizobium anhuiense as the predominant microsymbionts of Lathyrus maritimus along the Shandong Peninsula seashore line. Systematic and Applied Microbiology 39, 384e390. http://dx.doi.org/10.1016/j.syapm.2016.07.001. Ling, J., Wang, H., Wu, P., Li, T., Tang, Y., Naseer, N., Zheng, H., Masson-Boivin, C., 2016. Plant nodulation inducers enhance horizontal gene transfer of Azorhizobium caulinodans symbiosis island. Proceedings of the National Academy of Sciences U. S. A. 113, 13875e13880. http://dx.doi.org/10.1073/pnas.1615121113. Lu, R.K., 2000. Methods for Soil Agrochemistry Analysis. China Agricultural Science and Technology Press, Beijing, China (in Chinese). Macdonald, C.A., Clark, I.M., Hirsch, P.R., Zhao, F.J., McGrath, S.P., 2011. Development of a real-time PCR assay for detection and quantification of Rhizobium leguminosarum bacteria and discrimination between different biovars in zinccontaminated soil. Applied and Environmental Microbiology 77, 4626e4633. http://dx.doi.org/10.1128/AEM.02232-10.  bel, S., Skorupska, A., 2013. Marek-Kozaczuk, M., Leszcz, A., Wielbo, J., Wdowiak-Wro Rhizobium pisi sv. trifolii K3.22 harboring nod genes of the Rhizobium leguminosarum sv. trifolii cluster. Systematic and Applied Microbiology 36, 252e258. http://dx.doi.org/10.1016/j.syapm.2013.01.005. Martínez, L., Baquero, F., 2002. Interactions among strategies associated with bacterial Infection: pathogenicity, epidemicity, and antibiotic resistance. Clinical Microbiology Reviews 15, 647e679. http://dx.doi.org/10.1128/CMR.15.4.647. Moran, N.A., McCutcheon, J.P., Nakabachi, A., 2008. Genomics and evolution of heritable bacterial symbionts. Annual Review of Genetics 42, 165e190. http:// dx.doi.org/10.1146/annurev.genet.41.110306.130119. Mutch, L.A., Young, J.P.W., 2004. Diversity and specificity of Rhizobium leguminosarum biovar viciae on wild and cultivated legumes. Molecular Ecology 13, 2435e2444. http://dx.doi.org/10.1111/j.1365-294X.2004.02259.x. Nemergut, D.R., Schmidt, S.K., Fukami, T., O'Neill, S.P., Bilinski, T.M., Stanish, L.F., Knelman, J.E., Darcy, J.L., Lynch, R.C., Wickey, P., Ferrenberg, S., 2013. Patterns and processes of microbial community assembly. Microbiology and Molecular Biology Reviews 77, 342e356. http://dx.doi.org/10.1128/MMBR.00051-12. Oksanen, J., Blanchet, F.G., Kindt, R., Legendre, P., Minchin, P.R., O'Hara, R.B., Simpson, G.L., Solymos, P., Stevens, M.H.H., Wagner, H., 2015. Vegan: Community Ecology Package. R package version 2.3-0. https://CRAN.R-project.org/ package¼vegan. Parks, D.H., Porter, M., Churcher, S., Wang, S., Blouin, C., Whalley, J., Brooks, S., Beiko, R.G., 2009. GenGIS: a geospatial information system for genomic data. Genome Research 19, 1896e1904. http://dx.doi.org/10.1101/gr.095612.109. rez-Ramírez, N.O., Rogel, M.A., Wang, E., Castellanos, J.Z., Martírez-Romero, E., Pe 1998. Seeds of Phaseolus vulgaris bean carry Rhizobium etli. FEMS Microbiology Ecology 26, 289e296. http://dx.doi.org/10.1111/j.1574-6941.1998.tb00513.x. Pobigaylo, N., Szymczak, S., Nattkemper, T.W., Becker, A., 2008. Identification of genes relevant to symbiosis and competitiveness in Sinorhizobium meliloti using signature-tagged mutants. Molecular Plant Microbe Interactions 21, 219e231. http://dx.doi.org/10.1094/MPMI-21-2-0219. Quince, C., Lanzen, A., Davenport, R.J., Turnbaugh, P.J., 2011. Removing noise from pyrosequenced amplicons. BMC Bioinformatics 12, 38. http://dx.doi.org/ 10.1186/1471-2105-12-38. Ramírez-Bahena, M.H., García-Fraile, P., Peix, A., Valverde, A., Rivas, R., Igual, J.M., zquez, E., 2008. Revision of the taxonomic Mateos, P.F., Martínez-Molina, E., Vela status of the species 1889 AL, Rhizobium phaseoli Dangeard 1926 AL and later synonym of R. leguminosarum. Reclassification of the strain R. leguminosarum DSM 30132 (¼ NCIMB 11478) as Rhizobium pisi. International Journal of Systematic and Evolutionary Microbiology 58, 2484e2490. http://dx.doi.org/ 10.1099/ijs.0.65621-0. Remigi, P., Zhu, J., Young, J.P.W., Masson-Boivin, C., 2016. Symbiosis within symbiosis: evolving nitrogen-fixing legume symbionts. Trends in Microbiology 24, 63e75. http://dx.doi.org/10.1016/j.tim.2015.10.007. ~ o-Orrillo, E., Martinez Romero, E., 2011. Symbiovars in rhizobia Rogel, M.A., Ormen reflect bacterial adaptation to legumes. Systematic and Applied Microbiology 34, 96e104. http://dx.doi.org/10.1016/j.syapm.2010.11.015. € ~ iga, D., Alvarez-Martínez, Saïdi, S., Ramírez-Bahena, M.-H., Santillana, N., Zún E., zquez, E., 2014. Rhizobium laguerreae sp. nov. nodPeix, A., Mhamdi, R., Vela ulates Vicia faba on several continents. International Journal of Systematic and Evolutionary Microbiology 64, 242e247. http://dx.doi.org/10.1099/ijs.0.0521910. zquez, E., Zún ~ iga, D., Santillana, N., Ramírez-Bahena, M.H., García-Fraile, P., Vela 2008. Phylogenetic diversity based on rrs, atpD, recA genes and 16S-23S intergenic sequence analyses of rhizobial strains isolated from Vicia faba and Pisum sativum in Peru. Archives of Microbiology 189, 239e247. http://dx.doi.org/ 10.1007/s00203-007-0313-y. Schloss, P.D., Westcott, S.L., Ryabin, T., Hall, J.R., Hartmann, M., Hollister, E.B., Lesniewski, R.A., Oakley, B.B., Parks, D.H., Robinson, C.J., Sahl, J.W., Stres, B., Thallinger, G.G., Van Horn, D.J., Weber, C.F., 2009. Introducing mothur: opensource, platform-independent, community-supported software for describing and comparing microbial communities. Applied and Environmental Microbiology 75, 7537e7541. http://dx.doi.org/10.1128/AEM.01541-09. Sørensen, J., Nicolaisen, M.H., Ron, E., Simonet, P., 2009. Molecular tools in rhizosphere microbiology d from single-cell to whole-community analysis. Plant

and Soil 321, 483e512. http://dx.doi.org/10.1007/s11104-009-9946-8. Sullivan, J.T., Ronson, C.W., 1998. Evolution of rhizobia by acquisition of a 500-kb symbiosis island that integrates into a phe-tRNA gene. Proceedings of the National Academy of Sciences U. S. A. 95, 5145e5149. Sullivan, J.T., Patrick, H.N., Lowther, W.L., Scott, D.B., Ronson, C.W., 1995. Nodulating strains of Rhizobium loti arise through chromosomal symbiotic gene transfer in the environment. Proceedings of the National Academy of Sciences U. S. A. 92, 8985e8989. http://dx.doi.org/10.1073/pnas.92.19.8985. Tamura, K., Stecher, G., Peterson, D., Filipski, A., Kumar, S., 2013. MEGA6: molecular evolutionary genetics analysis version 6.0. Molecular Biology and Evolution 30, 2725e2729. http://dx.doi.org/10.1093/molbev/mst197. Tanno, K., Willcox, G., 2006. The origins of cultivation of Cicer arietinum L. and Vicia faba L.: early finds from Tell el-Kerkh, north-west Syria, late 10th millennium BP. Vegetation History and Archaeobotany 15, 197e204. http://dx.doi.org/ 10.1007/s00334-005-0027-5. €rensen, J., Olsson, S., 2000. ViscosinamideThrane, C., Nielsen, T.H., Nielsen, M.N., So producing Pseudomonas fluorescens DR54 exerts a biocontrol effect on Pythium ultimum in sugar beet rhizosphere. FEMS Microbiology Ecololgy 33, 139e146. http://dx.doi.org/10.1111/j.1574-6941.2000.tb00736.x. Tian, C.F., Wang, E.T., Han, T.X., Sui, X.H., Chen, W.X., 2007. Genetic diversity of rhizobia associated with Vicia faba in three ecological regions of China. Archives of Microbiology 188, 273e282. http://dx.doi.org/10.1007/s00203-007-0245-6. Tian, C.F., Wang, E.T., Wu, L.J., Han, T.X., Chen, W.F., Gu, C.T., Gu, J.G., Chen, W.X., 2008. Rhizobium fabae sp nov., a bacterium that nodulates Vicia faba. International Journal of Systematic and Evolutionary Microbiology 58, 2871e2875. http://dx.doi.org/10.1099/ijs.0.2008/000703-0. Tian, C.F., Young, J.P.W., Wang, E.T., Tamimi, S.M., Chen, W.X., 2010. Population mixing of Rhizobium leguminosarum bv. viciae nodulating Vicia faba: the role of recombination and lateral gene transfer. FEMS Microbiology Ecology 73, 563e576. http://dx.doi.org/10.1111/j.1574-6941.2010.00909.x. Tian, C.F., Zhou, Y.J., Zhang, Y.M., Li, Q.Q., Zhang, Y.Z., Li, D.F., Wang, S., Wang, J., Gilbert, L.B., Li, Y.R., Chen, W.X., 2012. Comparative genomics of rhizobia nodulating soybean suggests extensive recruitment of lineage-specific genes in adaptations. Proceedings of the National Academy of Sciences U. S. A. 109, 8629e8634. http://dx.doi.org/10.1073/pnas.1120436109. Vellend, M., 2010. Conceptual synthesis in community ecology. The Quarterly Review of Biology 85, 183e206. Vos, M., Quince, C., Pijl, A.S., de Hollander, M., Kowalchuk, G.A., 2012. A Comparison of rpoB and 16S rRNA as markers in pyrosequencing studies of bacterial diversity. PLoS One 7, e30600. http://dx.doi.org/10.1371/journal.pone.0030600. Wang, H.-F., Zong, X.-X., Guan, J.-P., Yang, T., Sun, X.-L., Ma, Y., Redden, R., 2012. Genetic diversity and relationship of global faba bean (Vicia faba L.) germplasm revealed by ISSR markers. Theoretical and Applied Genetics 124, 789e797. http://dx.doi.org/10.1007/s00122-011-1750-1. Wang, D., Wang, Y.C., Wu, L.J., Liu, J.X., Zhang, P., Jiao, J., Yan, H., Liu, T., Tian, C.F., Chen, W.X., 2016. Construction and pilot screening of a signature-tagged mutant library of Sinorhizobium fredii. Archives of Microbiology 198, 91e99. http://dx.doi.org/10.1007/s00203-015-1161-9. Young, J.P.W., Crossman, L.C., Johnston, A.W., Thomson, N.R., Ghazoui, Z.F., Hull, K.H., Wexler, M., Curson, A.R., Todd, J.D., Poole, P.S., Mauchline, T.H., East, A.K., Quail, M.A., Churcher, C., Arrowsmith, C., Cherevach, I., Chillingworth, T., Clarke, K., Cronin, A., Davis, P., Fraser, A., Hance, Z., Hauser, H., Jagels, K., Moule, S., Mungall, K., Norbertczak, H., Rabbinowitsch, E., Sanders, M., Simmonds, M., Whitehead, S., Parkhill, J., 2006. The genome of Rhizobium leguminosarum has recognizable core and accessory components. Genome Biology 7, R34. http://dx.doi.org/10.1186/gb-2006-7-4-r34. Zeder, M.A., 2009. The neolithic macro-(r)evolution: macroevolutionary theory and the study of culture change. Journal of Archaeological Science 17, 1e63. http:// dx.doi.org/10.1007/s10814-008-9025-3. ze , A., Mutch, L.A., Young, J.P.W., 2001. Direct amplification of nodD from comZe munity DNA reveals the genetic diversity of Rhizobium leguminosarum in soil. Environmental Microbiology 3, 363e370. http://dx.doi.org/10.1046/j.14622920.2001.00202.x. Zhang, Y.M., Tian, C.F., Sui, X.H., Chen, W.F., Chen, W.X., 2012. Robust markers reflecting phylogeny and taxonomy of rhizobia. PLoS One 7, e44936. http:// dx.doi.org/10.1371/journal.pone.0044936. Zhang, X.X., Guo, H.J., Wang, R., Sui, X.H., Zhang, Y.M., Wang, E.T., Tian, C.F., Chen, W.X., 2014. Genetic divergence of Bradyrhizobium nodulating soybeans as revealed by multilocus sequence analysis of genes inside and outside the symbiosis island. Applied and Environmental Microbiology 80, 3181e3190. http://dx.doi.org/10.1128/AEM.00044-14. Zhang, Y.J., Zheng, W.T., Everall, I., Young, J.P.W., Zhang, X.X., Tian, C.F., Sui, X.H., Wang, E.T., Chen, W.X., 2015. Rhizobium anhuiense sp.nov., isolated from effective nodules of Vicia faba and Pisum sativum. International Journal of Systematic and Evolutionary Microbiology 65, 2960e2967. http://dx.doi.org/ 10.1099/ijs.0.000365. Zhang, X.X., Guo, H.J., Jiao, J., Zhang, P., Xiong, H.Y., Chen, W.X., Tian, C.F., 2017. Pyrosequencing of rpoB uncovers a significant biogeographical pattern of rhizobial species in soybean rhizosphere. Journal of Biogeography 44, 1491e1499. http://dx.doi.org/10.1111/jbi.12891. Zong, X., Liu, X., Guan, J., Wang, S., Liu, Q., Paull, J.G., Redden, R., 2009. Molecular variation among Chinese and global winter faba bean germplasm. Theoretical and Applied Genetics 118, 971e978. http://dx.doi.org/10.1007/s00122-0080954-5.