Tree rows in temperate agroforestry croplands alter the composition of soil bacterial communities

Abstract

Background

Tree-based intercropping (agroforestry) has been advocated to reduce adverse environmental impacts of conventional arable cropping. Modern agroforestry systems in the temperate zone are alley-cropping systems that combine rows of fast-growing trees with rows of arable crops. Soil microbial communities in these systems have been investigated intensively; however, molecular studies with high taxonomical resolution are scarce.

Methods

Here, we assessed the effect of temperate agroforestry on the abundance, diversity and composition of soil bacterial communities at three paired poplar-based alley cropping and conventional monoculture cropland systems using real-time PCR and Illumina sequencing of bacterial 16S rRNA genes. Two of the three systems grew summer barley (Hordeum vulgare); one system grew maize (Zea mays) in the sampling year. To capture the spatial heterogeneity induced by the tree rows, soil samples in the agroforestry systems were collected along transects spanning from the centre of the tree rows to the centre of the agroforestry crop rows.

Results

Tree rows of temperate agroforestry systems increased the abundance of soil bacteria while their alpha diversity remained largely unaffected. The composition of the bacterial communities in tree rows differed from those in arable land (crop rows of the agroforestry systems and conventional monoculture croplands). Several bacterial groups in soil showed strong association with either tree rows or arable land, revealing that the introduction of trees into arable land through agroforestry is accompanied by the introduction of a tree row-associated microbiome.

Conclusion

The presence of tree row-associated bacteria in agroforestry increases the overall microbial diversity of the system. We speculate that the increase in biodiversity is accompanied by functional diversification. Differences in plant-derived nutrients (root exudates and tree litter) and management practices (fertilization and tillage) likely account for the differences between bacterial communities of tree rows and arable land in agroforestry systems.

About the Authors:

Lukas Beule

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Validation, Visualization, Writing – original draft

* E-mail: [email protected]

Affiliation: Molecular Phytopathology and Mycotoxin Research, Faculty of Agricultural Sciences, University of Goettingen, Goettingen, Germany

ORCID logo https://orcid.org/0000-0003-1107-7540

Petr Karlovsky

Roles Conceptualization, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Validation, Visualization, Writing – review & editing

Affiliation: Molecular Phytopathology and Mycotoxin Research, Faculty of Agricultural Sciences, University of Goettingen, Goettingen, Germany

ORCID logo https://orcid.org/0000-0002-6532-5856

Introduction

Modern temperate agroforestry systems are monoculture alley-cropping systems where rows of fast-growing trees (e.g. poplar (Populus) species) are alternated with rows of annual crops (tree-based intercropping). The cultivation of trees and crops in close spatial proximity allows various interspecific interactions that can result in complementary use of the resources [1]. Partition of resources between trees and crops in agroforestry systems is regarded as the key advantage of agroforestry over conventional arable croplands [2]. Agroforestry systems can reduce nitrate leaching through nitrate uptake by deep-rooting tree roots expanding below the cropland zone (‘safety-net’ role of tree roots) [3–5] and increase soil fertility through tree-litter input [6, 7]. Furthermore, agroforestry can increase faunal and floral diversity compared to conventional monoculture cropland, as reviewed by Udawatta et al. [8]. Temperate agroforestry consistently diminished crop yield close to the trees [9–11], yet maintained food safety requirements for small-grain cereals [12]. The benefits of integrating trees in agricultural systems have been recognized and improvement of the sustainability of agriculture through agroforestry has been proposed [13, 14].

Over the last two decades, soil microbial communities in temperate agroforestry systems have been extensively investigated using traditional methods such as fumigation-extraction, enzyme activities, and substrate-induced respiration. Soil microbial biomass has been shown to be greater in agroforestry than in conventional monoculture cropland [15, 16], whereas in poplar-based agroforestry, the increase in microbial biomass was limited to the tree rows [17]. Similarly, the bacterial biomass in an adler (Alnus rubra)-maize (Zea mays) agroforestry system was greater in the vicinity of the tree rows [18]. Furthermore, numerous studies of soil microbial communities in temperate agroforestry systems found more diverse catabolic potential [19–23] and greater efficiency of substrate use under agroforestry than in conventional cropland [17, 24]. It is reasonable to assume that these differences resulted from compositional differences among soil microbial communities. In 2013, Bardhan and co-workers [25] tested this hypothesis by investigating the composition and diversity of soil bacterial communities in a maple (Acer saccharinum)-based agroforestry system using denaturing gradient gel electrophoresis (DGGE). The authors did not detect any difference between the composition or diversity of soil bacteria in the tree and crop rows, which was likely due to limited coverage of complex communities by DGGE [26]. Despite the wide application of next generation sequencing (NGS) techniques in soil microbiology, sequencing of amplified segments of ribosomal RNA genes or spacers, which provides a high sampling depth and taxonomical resolution for cultivable and non-cultivable microorganism, has rarely been used in studies of agroforestry. Recently, Banerjee et al. [27] investigated bacterial communities in Canadian agroforestry systems using both quantitative (real-time PCR (qPCR)) and qualitative (amplicon sequencing) molecular approaches. The authors concluded that temperate agroforestry did not promote soil bacterial diversity but it increased bacterial abundance [27]. Likewise, a recent molecular investigation revealed that tree rows of agroforestry systems promoted the abundance of several bacterial groups in soil [28]. Additionally, temperate agroforestry has been shown to alter the abundance of microorganisms involved in N cycling [28, 29]. Despite these recent insights, systematic investigations of soil bacterial communities in temperate agroforestry systems with high taxonomical resolution are missing.

The present study aimed to assess the effect of poplar-based temperate alley cropping (agroforestry) on bacterial communities in soil. We expected that differences in the distribution of tree litter input [11] as well as the fertilization and tillage regime [30] between agroforestry tree and crop rows and monoculture croplands would result in changes of the soil microbiome as reported previously [27–29]. Therefore, we hypothesized that the integration of rows of trees into arable land through agroforestry i) increased the abundance and ii) diversity of soil bacteria and iii) affected the composition of the soil bacterial community along transects from the centres of the tree rows to the centres of the agroforestry crop rows. Furthermore, this study aimed to establish a novel technique for library normalization using qPCR and gel densitometry.

Materials and methods

Study sites and experimental design

The land owners gave their permission to conduct soil sampling on their property. We chose three study sites with paired agroforestry and conventional monoculture cropland in Germany (Fig 1). The soil types at the three study sites were Calcaric Phaeozem soil (near Dornburg, Thuringia), Gleyic Cambisol soil (near Forst, Brandenburg), and Vertic Cambisol soil (near Wendhausen, Lower Saxony) (Table 1). The agroforestry systems were established between 2007 and 2010 by converting conventional monoculture cropland into an alley-cropping system. At every site, 12-m wide rows of poplar trees were planted from cuttings (clone Max1; Populus nigra × P. maximowiczii) using a dibble bar. The tree rows were North-South oriented and interspersed with 48-m wide rows of crops (Fig 1B). The agroforestry crop rows were managed identically as the monoculture croplands on the same sites (identical crop rotation, fertilization, and pesticide treatment). The crop rotations included maize (Zea mays), summer and winter barley (Hordeum vulgare), winter oilseed rape (Brassica napus), and winter wheat (Triticum aestivum) (Table 1). Fertilizer was applied to the crop rows of the agroforestry systems and the monoculture croplands in spring according to standard practice (Table 1). Following common temperate agroforestry practice, the trees rows of the agroforestry systems did not receive fertilizer [30].

[Figure omitted. See PDF.]

Fig 1.

Study sites in Germany (A), experimental design of paired agroforestry and conventional monoculture cropland (B), and a picture taken at the tree-crop interface of the agroforestry cropland system at Dornburg (C). The soil types at the three study sites were Calcaric Phaeozem soil (near Dornburg, Thuringia), Gleyic Cambisol soil (near Forst, Brandenburg), and Vertic Cambisol soil (near Wendhausen, Lower Saxony). Photo taken by G. Shao.

https://doi.org/10.1371/journal.pone.0246919.g001

[Figure omitted. See PDF.]

Table 1. Site characteristics and management at the three study sites of paired temperate agroforestry and monoculture cropland.

https://doi.org/10.1371/journal.pone.0246919.t001

At every site, we established four replicate plots in both the agroforestry and the monoculture cropland systems (Fig 1B). The replicate plots in the agroforestry systems were linear transects orthogonal to the North-South orientation of the tree rows. These transects spanned from the centre of the 12-m wide tree row to the centre of the 48-m wide agroforestry crop row. Soil samples were collected along the transects within the centre of the tree row (approximately 1 m from the poplar trunk) as well as in the agroforestry crop row at 1 m, 7 m, and 24 m (centre of the agroforestry crop row) distance from the trees (Fig 1B) (see https://doi.org/10.20387/BONARES-4984-ZWYR for spatial references of the sampling locations). In the adjacent monoculture systems, soil samples were collected in the centre of the plots. In total, 20 soil samples were collected at each site.

Soil collection and DNA extraction

In 2019, soil samples of the top 5-cm depth were collected on July 15 (Vertic Cambisol), July 16 (Calcaric Phaeozem), and August 6 (Gleyic Cambisol). At each of the 20 sampling locations, soil was collected using three 250 cm3 stainless steel cylinders (5 cm height) per site and immediately thoroughly homogenised in a sterile plastic bag. Still in the field, an aliquot of 20 g of fresh soil was transferred into a 15 mL Falcon tube (SARSTEDT, Nümbrecht, Germany) and frozen at -20°C. Upon arrival in the laboratory, soil samples were freeze-dried for 72 h and finely ground using a swing mill (Retsch MM400, Retsch, Haan, Germany). DNA was extracted from 50 mg soil using a cetyltrimethylammonium bromide-based (CTAB) protocol described previously by Brandfass & Karlovsky [31]. Briefly, soil was suspended in 1 mL CTAB with 1 μL proteinase K, incubated at 42°C and subsequently at 65°C, and 800 μL phenol were added. The mixture was shaken, centrifuged, and the supernatant was extracted with aliquots of chloroform/isoamyl alcohol twice. DNA was precipitated using PEG/NaCl and pelleted by centrifugation. The pellets were washed with EtOH twice, dried, and re-dissolved in 50 μL 1 × TE buffer (10 mM Tris, 1 mM ethylenediaminetetraacetic acid (EDTA), adjusted to pH 8.0 with HCl). Extracted DNA was visualized in 0.8% (w/v) agarose gels (in 1 × TAE buffer (40 mM Tris, 20 mM sodium acetate, 1 mM Na2EDTA, adjusted to pH 7.6)) stained with ethidium bromide. Gel electrophoresis was carried out at 4.6 V cm-1 for 60 min. Soil DNA extracts were stored at -20°C. The extracts were tested for the absence of PCR inhibitors as described previously [32].

Determination of absolute abundance of bacteria

Prior to library preparation, soil bacteria in the 60 DNA extracts were quantified and normalized using qPCR. Amplifications were performed in triplicates in a CFX 384 Thermocycler (Biorad, Rüdigheim, Germany) in 384-well microplates. The reaction volume was 4 μL consisting of double-distilled water (ddH2O); buffer (10 mM Tris-HCl, 50 mM KCl, 2.0 mM MgCl2, pH 8.3 at 25°C); 200 μM of each deoxynucleoside triphosphate (Bioline, Luckenwalde, Germany); 0.4 μM of each primer (341F (5’-CCTACG GGNGGC WGCAG-3’)/785R (5’-GACTAC HVGGGT ATCTAA KCC-3’) [33], targeting the V3–V4 region of the bacterial 16S rRNA); 0.1 × SYBR Green I solution (Invitrogen, Karlsruhe, Germany)); 1 μg μL-1 bovine serum albumin; 0.03u μL-1 Hot Start Taq DNA Polymerase (New England Biolabs, Beverly, Massachusetts, USA) and 1 μL of template DNA (standards or 1:50 dilutions of the DNA extracts in 0.5 × TE) or 0.5 × TE for the negative control. Standards in 0.5 × TE were obtained from environmental soil DNA extracts and amplified in triplicates. Standard rows were obtained from 1:2 serial dilutions of 1 pg PCR product. The quantification cycle (Cq) values of the standards that covered the range of the samples were plotted against their log-transformed starting quantity. Curve fitting was conducted to account for non-linear PCR efficiency over the range of the standard row. A third order polynomial function was fitted and the absolute abundance of bacteria in the extracts was determined.

Library normalization

To achieve a narrow distribution of the titre of libraries for NGS, the concentration of bacterial 16S rRNA genes in sample extracts was determined by qPCR (see Determination of absolute abundance of bacteria) and the extracts were diluted in 0.5 × TE to achieve the same concentration of rRNA genes in all samples. The success of normalization was tested by qPCR. The normalization was regarded as successful when the difference among Cq for all sample pairs was smaller than one. The samples that did not meet this criterion were re-diluted until successful normalization was achieved.

A second normalization was carried out by diluting the libraries based on their DNA concentration determined by densitometry. A 2 μL aliquot of the each of the 60 libraries was diluted 1:10 in 0.5 × TE buffer and 4 μL of the dilutions were visualized on 1.7% (w/v) agarose gels in sets of 20 samples per gel, resulting in three gels for the 60 samples. A standard row of 11 standards was obtained from 3:1 serial dilutions of a library in 0.5 × TE buffer. The standards were loaded onto a 1.7% (w/v) agarose gel in duplicates. The densitometric quantification of library DNA was performed employing ImageJ version 1.52q [34]. A second order polynomial function was fitted to obtain a calibration curve (S1 Fig). The relative library yield of each of the 20 samples per gel was determined. All samples loaded onto one gel were normalized to the sample with the lowest library yield within the respective gel and pooled, resulting in one pooled sample per gel. To ensure successful normalization across gels, 4 μL aliquots of the three pooled samples (one from each gel) were loaded on onto a 1.7% (w/v) agarose gel in duplicates. The densitometric results were used to normalize the pooled samples as described above for individual libraries.

Library preparation and high-throughput Illumina sequencing

Amplifications for library preparation were carried out in a peqSTAR 96 universal gradient thermocycler (PEQLAB, Erlangen, Germany) in 25 μL reaction volumes consisting of ddH2O; buffer (10 mM Tris-HCl, 50 mM KCl, 2.0 mM MgCl2, pH 8.3 at 25°C); 200 μM of each deoxynucleoside triphosphate (Bioline, Luckenwalde, Germany); 0.4 μM of each primer (341F/785R [33]); 1 μg μL-1 bovine serum albumin; 0.03u μL-1 Hot Start Taq DNA Polymerase (New England Biolabs, Beverly, Massachusetts, USA) and 6.25 μL of normalized template DNA or 0.5 × TE for the negative control. The primers were a set of 48 dual-indexed primer pairs that included 0–3 frameshifting bases (Ns) to improve Illumina base-calling followed by an 8-bp barcode at the 5’-end of each primer. Thermocycling conditions were as follows: initial denaturation (95°C for 2 min), 3 touch-up cycles of denaturation (95°C for 20 sec), annealing (50°C for 30 sec), and elongation (68°C for 30 sec), 25 cycles of denaturation (95°C for 20 sec), annealing (58°C for 30 sec), and elongation (68°C for 30 sec), and final elongation (68°C for 5 min). All libraries were prepared within one PCR run using the same mastermix. The libraries were visualized on 1.7% (w/v) agarose gels (in 1 × TAE buffer) stained with ethidium bromide. Gel electrophoresis was carried out at 4.6 V cm-2 for 60 min.

Based on the densitometry, the libraries were normalized (see Library normalization) and sent to LGC Genomics (Berlin, Germany) for adapter ligation using a commercial kit (Ovation® Rapid DR Multiplex System 1–96 (NuGEN, San Carlos, CA, USA)). Finally, the libraries were sequenced in one multiplex sequencing run using the Illumina MiSeq Reagent Kit v3 (2 × 300 bp) (Illumina, San Diego, CA, USA) at the facilities of LGC Genomics, Berlin, Germany. Sequencing data have been deposited at NCBI’s Sequence Read Achieve (BioProject PRJNA667193).

Processing of sequencing data

Raw paired-end data (10,736,108 reads in total) were demultiplexed using Illumina’s bcl2fast version 2.17.1.14 (Illumina, San Diego, CA, USA) and sorted by their barcodes (allowing 1 mismatch per barcode; missing, one-sided or conflicting barcode pairs were discarded). Sequencing adapters were clipped and reads shorter than 100 bp were discarded. Primer were clipped (allowing 3 mismatches per primer) and reads were imported in QIIME 2 version 2019.10 [35]. We employed DADA2 [36] for quality-filtering, merging, chimera and singleton filtering of the reads, and clustering of reads in exact amplicon sequence variants (ASVs) [37, 38]. The dereplicated ASVs were taxonomically assigned against the SILVA SSU database (release 132) [39] using VSEARCH [40]. Non-bacterial reads were removed from the obtained ASV table. After filtering, 5,890,832 bacterial counts were obtained. We normalized the library size to 12,521 counts per sample using scaling with ranked subsampling [41] using the ‘SRS’-function in the ‘SRS’ R-package version 0.1.0 [41] in the R environment version 3.6.1 [42]. Finally, our normalized dataset contained 40,708 ASVs.

Quantification of nifH gene in soil

The most frequently used marker gene for N fixation by bacteria is nifH, which encodes for a subunit of the dinitrogenase reductase enzyme [43]. The abundance of nifH in soil was determined using qPCR as described previously. Briefly, qPCR reactions were performed in 4 μl reaction volume in 384-well microplates using a CFX384 Thermocycler (Bio-Rad, Rüdigheim, Germany). The reaction volume contained 3 μl mastermix and 1 μl of a 1:50 dilution in 0.5 × TE buffer of the DNA extracts or ddH2O for negative controls. The composition of the mastermix as well as primers and the thermocycling conditions have been described in details previously [28].

Statistical analysis

Alpha diversity was quantified by calculating the Shannon index (H’) and the Pielou’s evenness index (J’) for the ASV count data. H’ was determined using the ‘diversity’-function in the ‘vegan’ R-package version 2.5–6 [44]. J’ was determined by dividing H’ by the natural logarithm of the total number of species. The effect of sampling location (tree row, 1 m, 7 m, and 24 m distance from the tree row within the agroforestry crop row and monoculture cropland) within one soil type on alpha diversity measures, bacterial 16S rRNA and nifH gene abundance was determined using one-way analysis of variance (ANOVA) with Tukey’s honestly significant difference (HSD) test or Kruskal–Wallis test with multiple comparison extension.

The Bray-Curtis index of dissimilarity for pairwise comparisons was calculated from the square root-transformed ASV count data using the ‘vegdist’-function in the ‘vegan’ R-package version 2.5–6 [44]. The Bray-Curtis dissimilarities were visualized using non-metric multidimensional scaling (NMDS) with 100 random starts (‘metaMDS’-function in the ‘vegan’ R-package version 2.5–6 [44]). Additionally, we calculated the Bray-Curtis index of similarity (1 –Bray-Curtis dissimilarity) and visualized the similarity across samples using a network approach employing the ‘igraph’ R-package version 1.2.4.2 [45]. For the network visualization, edges were only drawn between samples with a Bray-Curtis index of similarity above the 75th quantile. The network was spatially arranged applying the Fruchterman-Reingold algorithm [46]. The effect of soil type, sampling location, and soil type × sampling location on community composition was determined employing permutational multivariate analysis of variance (PERMANOVA) with 999 permutations using the ‘adonis2’-function in the ‘vegan’ R-package version 2.5–6 [44]). Complementary, we tested for the multivariate homogeneity of group dispersions (PERMDISP) on Bray-Curtis dissimilarities with 999 permutations (‘betadisper’-function in the ‘vegan’ R-package version 2.5–6 [44]).

The relative abundance of all identified genera was visualized and inspected. The effect of sampling locations within one soil type on the absolute abundance of bacterial 16S rRNA and nifH genes and the relative abundance of bacterial genera was determined using one-way ANOVA with Tukey’s HSD test or Kruskal–Wallis test with multiple comparison extension. One-way ANOVA with Tukey’s HSD was partly carried out on square root- or log-transformed data. For visualization purposes, the relative abundances of selected genera were Z-score normalized. Z-scores (Zi) for the relative abundance of each genus was calculated by subtracting the mean relative abundance of a genus (μ) from the relative abundance the genus in sample i (xi) and dividing the difference it by the standard deviation of the genus over all samples (σ):

Statistical analysis was performed in the R environment version 3.6.1 [42].

Results

Validation of two-step normalization of NGS libraries

Normalization of DNA samples prior to library preparation based on the quantification of bacterial 16S rRNA genes by qPCR substantially reduced differences in library size among samples (Fig 2A). To achieve a narrow distribution of Cq values for libraries, DNA extracts had to be diluted up to 7.5-fold. The libraries were further normalized by dilution up to 2.4-fold, based on the quantification of DNA concentration in the libraries by densitometry (Fig 2B). The number of raw sequence reads obtained from the normalized libraries ranged from 28,935 to 243,974 reads per library with an interquartile range from 137,234 to 192,129 (Fig 2C).

[Figure omitted. See PDF.]

Fig 2.

Kernel density estimation of the quantification cycle (Cq) values for non-normalised and normalised samples (prior to library preparation) (A), prepared libraries normalized by photo densitometry (B), and the distribution of the raw sequence reads per library (C). Individual Cq values are plotted as marks below the density curves for non-normalised and above for normalised samples; curves were smoothened using a bandwidth of 0.3 (A). The box indicates the q25, q50, and q75, whiskers range from q25 or q75 to 1.5 × the interquartile range, dots represent individual data points (C). n = 60 samples.

https://doi.org/10.1371/journal.pone.0246919.g002

Overall soil bacterial community

The determination of the soil bacterial community size by qPCR revealed no differences among sampling location on the Calcaric Phaeozem and Gleyic Cambisol soil, except lower bacterial abundance on the Calcaric Phaeozem soil at 7 m distance from the tree row within the agroforestry crop row compared to the monoculture system (p = 0.05) (Fig 3A). On the Vertic Cambisol, soil bacteria were more abundant in the tree row than at all distances from the tree row within the crop row of the agroforestry system as well as in the monoculture cropland system (p ≤ 0.023) (Fig 3A).

[Figure omitted. See PDF.]

Fig 3.

16S rRNA gene abundance of soil bacteria (A) and relative abundance of taxonomic groups of soil bacteria at class level (B) in paired temperate agroforestry and monoculture cropland systems in three different soil types. Samples in the agroforestry systems were collected in the tree row as well as at 1 m, 7 m, and 24 m distance from the tree row within the agroforestry crop row. Bars represent means with standard deviation (n = 4 per soil type × sampling location) (A). Different uppercase letters indicate significant differences among the sampling locations (one-way analysis of variance with Tukey’s honestly significant difference test or Kruskal–Wallis test with multiple comparison extension) (A). Soil bacterial classes with an overall relative abundance of < 0.5% were considered as rare classes and merged with unassigned classes (B).

https://doi.org/10.1371/journal.pone.0246919.g003

The 40,708 ASVs found in the normalized bacterial reads were distributed among 10 dominant (≥ 0.5% overall relative abundance across all samples) and 28 rare phyla (< 0.5% overall relative abundance across all samples). Actinobacteria (38.8%) followed by Proteobacteria (21.8%), Acidobacteria (10.7%), Chloroflexi (9.0%), and Planctomycetes (5.5%) were the most abundant phyla. Overall, 99.8% of all taxonomic groups were assigned at phylum level. Alphaproteobacteria (13.9%), Actinobacteria (19.9%), Thermoleophilia (12.7%), Subgroup 6 (5.5%), and Gammaproteobacteria (5.3%) were the most abundant bacterial classes (Fig 3B). At genus level, Nocardioides (3.2%), Microlunatus (3.0%), and Sphingomonas spp. (2.9%) were most abundant across samples.

Diversity and composition of soil bacterial communities

Shannon diversity of ASVs on the Calcaric Phaeozem soil was greater in the agroforestry crop row at 7 m distance from the tree row than in the tree row, 1 m within the crop row and the monoculture system (p ≤ 0.048) (Fig 4A). On the Vertic Cambisol and Gleyic Cambisol soil, no differences in Shannon diversity among sampling locations were present. On the Calcaric Phaeozem and Gleyic Cambisol soil, Pielou’s evenness was lower in the tree row than in the agroforestry crop row and the monoculture (p ≤ 0.038) (Fig 4B). On the Vertic Cambisol soil, no differences in evenness were found among sampling locations.

[Figure omitted. See PDF.]

Fig 4.

Alpha diversity measures (Shannon index (H’) (A) and Pielou’s evenness index (J’) (B)), non-metric multidimensional scaling (NMDS) of Bray-Curtis dissimilarities (C), and Bray-Curtis similarity network (D) of soil bacterial ASVs in paired temperate agroforestry and monoculture cropland systems in three different soil types. Samples in the agroforestry systems were collected in the tree row as well as at 1 m, 7 m, and 24 m distance from the tree row within the agroforestry crop row. Coloured vertical bars represent the standard deviation, coloured horizontal bars the mean (A, B). Coloured shapes represent individual data points (A, B). Different uppercase letters indicate statistically significant differences at p < 0.05 (one-way analysis of variance with Tukey’s honestly significant difference test or Kruskal–Wallis test with multiple comparison extension) (A, B). Solid lines in the NMDS span from the centroid of each group to the individual data points (C). Edges in the similarity network were drawn between samples with a Bray-Curtis similarity ≥ 0.75 (D). The network was spatially arranged applying the Fruchterman-Reingold algorithm (D). AF = agroforestry system, ASV = amplicon sequence variant.

https://doi.org/10.1371/journal.pone.0246919.g004

The soil bacterial community composition was strongly influenced by soil types (PERMANOVA; F = 15.021, p = 0.001), which explained 28.8% of the observed variance (S1 Table). Furthermore, the interaction of soil type × sampling location explained 16.2% of the variance and was found to influence the community composition (PERMANOVA; F = 2.110, p = 0.001). The influence of sampling locations within each soil type on community composition was reflected by the large proportion of ASVs unique to the different sampling locations (S2 Fig). Furthermore, proportions of ASVs shared between tree rows and arable land (crop rows and monocultures) generally declined with distance from the tree rows (S2 Fig). Finally, the different sampling locations across soil types explained 12.0% of the variance and affected the composition of the bacterial community (PERMANOVA; F = 3.120, p = 0.001). Multivariate homogeneity of group dispersions was given for soil type, sampling location, and soil type × sampling location (PERMDISP; p ≥ 0.057). These results were reflected by the NMDS (Fig 4C) and the similarity network (Fig 4D).

Relative abundance of taxonomic soil bacterial groups

We observed that the relative abundance of several genera of soil bacteria was enhanced by the integration of poplar rows in agricultural systems (Fig 5). For example, members of the genera Actinomycetospora (p ≤ 0.004), Bradyrhizobium (p ≤ 0.041), Flavobacterium (p ≤ 0.007), Mesorhizobium (p ≤ 0.010), Reyranella (p < 0.001), and Sporocytophaga (p ≤ 0.028) were more abundant in the tree row compared to the arable land (agroforestry crop row and the monoculture) in all three soil types. In all three soil types, the monoculture showed lower abundance of Chthoniobacter spp. as compared to the tree row (p ≤ 0.039). The genus Dinghuibacter was more abundant in the row on the Gleyic and the Vertic Cambisol soil than in the arable land (p ≤ 0.038). On the Gleyic Cambisol soil, Legionella affiliates showed greater relative abundance under the trees than in the agroforestry crop row and the monoculture (p ≤ 0.007).

[Figure omitted. See PDF.]

Fig 5.

Z-score normalized relative abundance (A) and relative abundance (B) of selected taxonomic soil bacterial genera in paired temperate agroforestry and monoculture cropland systems in three different soil types. Samples in the agroforestry systems were collected in the tree row as well as at 1 m, 7 m, and 24 m distance from the tree row within the agroforestry crop row. Horizontal bars represent the mean relative abundance, coloured shapes represent individual data points (B).

https://doi.org/10.1371/journal.pone.0246919.g005

In contrast, several groups showed a decline in their relative abundance in response to the trees (Fig 5). For example, in all three soil types, the relative abundance of Pseudarthrobacter affiliates was lower in the tree row than the arable land (p ≤ 0.021). On the Calcaric Phaeozem and Vertic Cambisol soil, the genus Paenibacillus showed lower abundance under the trees than in the arable land (p ≤ 0.015). On the Gleyic Cambisol soil, Paenibacillus spp. were lower in the tree row than in the monoculture (p ≤ 0.002). Likewise, lower relative abundance of Flavisolibacter (p ≤ 0.012) and Nitrosospira spp. (p ≤ 0.016) in the tree row versus the monoculture were observed in all three soil types. Nitrospira spp. were less abundant under the trees than in the arable land on the Calcaric Phaeozem soil (p ≤ 0.04). Additionally, Nitrospira spp. showed lower abundance under the trees than in the agroforestry crop row on the Vertic Cambisol soil (p ≤ 0.003). On the Calcaric Phaeozem, Bacillus affiliates showed lower abundance in the tree row than at 7 m and 24 m in the agroforestry crop row as well as in the monoculture (p ≤ 0.02), while on the Vertic Cambisol soil, the genus was less abundant in the tree row than in the agroforestry crop row (p ≤ 0.005).

Abundance of N fixation gene nifH in soil

The abundance of nifH gene on the Vertic Cambisol soil was greater in the tree row than in the crop row of the agroforestry system as well as in the monoculture system (p < 0.02) (Fig 6). In contrast, no differences in the abundance of nifH genes were found among sampling locations on the Calcaric Phaeozem and the Gleyic Cambisol soils.

[Figure omitted. See PDF.]

Fig 6. Abundance of N fixation gene nifH in paired temperate agroforestry and monoculture cropland systems in three different soil types.

Samples in the agroforestry systems were collected in the tree row as well as at 1 m, 7 m, and 24 m distance from the tree row within the agroforestry crop row. Bars represent means with standard deviation (n = 4 per soil type × sampling location). Different uppercase letters indicate significant differences among the sampling locations within the same soil type (one-way analysis of variance with Tukey’s honestly significant difference test or Kruskal–Wallis test with multiple comparison extension).

https://doi.org/10.1371/journal.pone.0246919.g006

Discussion

Main impacts of temperate agroforestry on soil bacterial communities

Our results revealed that poplar tree rows in temperate agroforestry increase the abundance of soil bacteria and alter the composition of the soil bacterial community, which is in line with previous findings [27, 28]. Although the integration of trees in arable land did not increase soil bacterial diversity, which agrees with the results of Banerjee et al. [27], tree row-associated soil bacteria enhanced the overall diversity of the system. Furthermore, the relative abundance of genera involved in soil N cycling was affected: N-fixing bacteria of the genera Bradyrhizobium and Mesorhizobium were promoted by the trees, whereas nitrifying bacteria of the genera Nitrosospira and Nitrospira were less abundant under the trees as compared to agroforestry crop row and the monoculture cropland. Overall, our study demonstrated that temperate agroforestry practice introduces a tree row-associated bacterial microbiome, which likely affects soil functional diversity and nutrient cycling.

Abundance and alpha diversity of soil bacterial communities in temperate agroforestry systems

The greater abundance of soil bacteria under the trees on the Vertic Cambisol soil (Fig 3A) confirmed our first hypothesis agrees with our previous findings that poplar rows in temperate agroforestry systems benefit soil microorganisms [28]. Likewise, Banerjee et al. [27] reported greater bacterial abundance in plots with trees than without trees in Canadian agroforestry systems. Also in line with the findings of Banerjee et al. [27] temperate agroforestry did not affect the Shannon diversity of the soil bacterial ASVs (Fig 4A). The evenness of soil bacterial ASVs was lower in the tree rows than in the arable land of the two agroforestry systems planted with summer barley (Calcaric Phaeozem and Gleyic Cambisol soil) (Fig 4B), revealing that the tree rows harboured a greater proportion of dominant ASVs. As reported previously [28], all our agroforestry systems harboured herbaceous vegetation of high biodiversity growing below the poplar trees. Therefore, our results are consistent with the assumption that aboveground plant alpha diversity does not predict the belowground microbial alpha diversity [47]. Since the poplar trees at all our study sites originated from the same poplar clone, the response of soil bacterial diversity to agroforestry was affected by the soil type and/or management regime (crop rotation, fertilization, and soil management).

Effect of trees in agroforestry systems on soil bacterial groups involved in soil N cycling

The relative abundance of the alpha-rhizobial genera Bradyrhizobium and Mesorhizobium were promoted by trees compared to the arable land (Fig 5). Bradyrhizobium and Mesorhizobium affiliates living in symbiosis with legumes are well-known for their capability to improve plant supplementation with N by assimilating atmospheric N2 [48, 49]. To further elucidate whether poplar-based agroforestry affects the genetic potential for N fixation, we quantified the nifH gene in soil. In a strong agreement with our previous findings [28], we observed that on the Vertic Cambisol, the tree row increased the absolute abundance of nifH compared to the arable land, whereas this was not found in the other two soil types (Fig 6). Although poplars are non-nodulating plants, plants found in the herbal vegetation below the trees on the Vertic Cambisol may be capable of nodulation and, thus, could contribute to N fixation. If confirmed, N fixation by herbal vegetation in the tree rows of temperate agroforestry systems would be a yet unaccounted benefit of temperate agroforestry practice.

The lower relative abundance of the nitrifiers comprising Nitrosospira and Nitrospira spp. under the trees as compared to the arable land (Fig 5) is congruent with our previous findings of reduced absolute abundance of ammonium-oxidizing bacteria in tree rows of temperate agroforestry systems as compared to the agroforestry crop rows and monoculture systems [28, 29]. In agreement with our results, the presence of trees reduced the relative abundance of Nitrosospira in Canadian agroforestry systems as compared to arable land [27]. The reduction of Nitrosospira and Nitrospira spp. under the trees was likely related to the absence of fertilization in the tree row [50]. Considering the increased abundance of genes involved in denitrification in the tree rows of temperate agroforestry systems [28], we conclude that agroforestry alters the genetic potential for N fixation, nitrification, and denitrification.

Effect of agroforestry on the composition of soil bacterial community

Our results demonstrated that the composition of soil bacterial community was strongly affected by the soil type (S1 Table), which was expected considering the strong effect of soil properties on the assembly of microbial populations (e.g. [51]). Altogether, our data reveals that poplar rows in temperate cropland agroforestry system harbour a different soil bacterial community than the arable land. Furthermore, the composition of the soil bacterial community in each soil type was influenced by the sampling location within the agroforestry system (tree row, 1 m, 7 m, and 24 m distance from the tree row) or in the conventional monoculture cropland without trees (S1 Table, S2 Fig), confirming our third hypothesis. In particular, soil bacterial communities collected in the tree rows were well separated from those from arable land (agroforestry crop row and monoculture cropland) as revealed by the clustering in the NMDS (Fig 4C). Therefore, the integration of trees into crop production through agroforestry diversifies the soil microbiome by introducing tree row-associated soil bacteria that enhance the overall diversity of the system, confirming our second hypothesis. These differences on DNA-level likely reflect functional differences, as reported previously for soil microbial communities of agroforestry systems [19, 20, 52]. We do not believe, however, that the functions of bacteria in soil can be accurately predicted from partial sequences of their 16S rRNA genes, considering the poor discriminatory power of the method at species level [53]. Furthermore, biological, chemical, and physical interactions in soil [54] strongly modulate the function of the soil microbiome.

There is substantial evidence that plant communities shape the composition of soil bacterial communities [47, 55–57]. Several studies showed that this effect is mediated by plant root exudates (e.g. [58]). For example, the observation that different arable crops established different bacterial communities in their rhizosphere was accounted for plant host habitat and root exudation [55]. Secondary metabolites (syn. special metabolites) were shown to play an important role in the modulation of soil microflora by plants [59, 60]. The persistent and abundant poplar biomass and the herbaceous layer in the tree rows of our agroforestry systems [28] are expected to release large quantities of root exudates with a much higher chemical diversity than the root exudates secreted into soil by in crop plants in the agroforestry system and in conventional monoculture cropland. Surprisingly, the high diversity of the herbal vegetation under the trees and their root exudates were not associated with an increase of alpha diversity of soil bacteria (Fig 4A). We assume that the tree litter and root exudation of the trees exerted a dominant effect on soil bacterial populations due to their large biomass. Therefore, the biodiversity of herbal vegetation under the trees and its exudations did not lead to the diversification of soil bacteria.

The amount of aboveground tree litter input (leaves, twigs, and branches) typically decreases exponentially with increasing distance from the trees [61]. Previous studies showed that the assembly of soil bacterial communities was affected by the quantity [62] as well as the quality [63] of tree leaf litter. In line with our results discussed above, we assume that the accumulation, decay, and incorporation of tree litter into the soil of the tree rows were among the major effects shaping belowground microbial communities in the tree rows. Because primary decomposers of tree litter are fungi, a DNA-based inventory and comparative analysis of fungal soil community would advance our understanding of soil biology in agroforestry systems.

In addition to harboring completely different plant communities, the cropland was tilled and fertilized while these treatments were absent in the tree rows. It is well established that fertilization strongly affects bacterial community composition. Both the kind of fertilizer and its amount affect bacterial communities in soil [64]. Fertilization and tree litter provide nutrients for soil microflora. Thus, nutrient input appears to be the major factor fostering the differentiation of soil microbiome between tree rows and arable land in agroforestry systems. In addition to nutrient input, tillage has repeatedly been shown to affect belowground communities [65–67], which is plausible considering the physical impact that the tillage exerts on soil as microbial habitat [68]. Therefore, differences in nutrient input and tillage in temperate agroforestry cropland systems appear to be the major factors modulating the soil microbiome.

Library normalization using qPCR

It is well established that the generation of amplicon sequencing libraries using PCR inevitably induces amplification bias [69], which increases with the number of PCR cycles (e.g. [70]). Therefore, is it generally recommended to limit the number of PCR cycles. Additionally, we suggest normalizing the starting quantity of target DNA prior to PCR to reduce differences in amplification among samples. For example, competition of PCR products with the primers for primer binding sites will slow down the amplification of dominant templates in samples with larger quantity of DNA faster than in samples with low amount of DNA. Because minor templates are less likely to be affected, differences in starting DNA concentration may affect relative amplification efficiencies of DNA species among samples, contributing to PCR bias. On the other hand, stochastic phenomena affecting rare templates during the first PCR cycles are likely to be stronger in samples with low starting DNA concentrations. Differences in the amount of extracted DNA are unavoidable with environmental samples and/or in studies in which a treatment affects the abundance of the target organism(s). Since our preliminary experiments showed that the quantities of target DNA varied among sampling locations, we normalized all our samples to a narrow range of below one Cq value using qPCR prior to library preparation (Fig 2A). Furthermore, we used a limited number of 3 touch-up cycles followed by only 25 PCR cycles for library preparation to stop the amplification before entering the plateau phase of the PCR.

Conclusion

The composition of the bacterial communities in tree rows our temperate agroforestry systems differed from those in arable land (crop rows of the agroforestry systems and conventional monoculture croplands). Therefore, the integration of trees in arable land through agroforestry diversifies the soil microbiome by introducing and promoting tree row-associated soil bacteria that enhance the overall diversity of the system. Furthermore, we suggest that the compositional alterations induced by the tree rows result in functional diversification of the soil microbiome of agroforestry systems. Differences in plant-derived nutrients (root exudates and tree litter) and management practices (fertilization and tillage) likely contributed to the observed differences between bacterial communities of tree rows and arable land.

Supporting information

S1 Fig.

Agarose gel (A) and calibration curve (B) used for densitometric quantification of library yield. The agarose gel was loaded with duplicates of a 3:1 serial dilution of a library in 0.5 × TE buffer. The agarose concentration was 1.7% (w/v) and gel electrophoresis was carried out at 4.6 V cm-2 for 60 min. Gels were stained in 0.1% (w/v) ethidium bromide solution for 10 min and de-stained in demineralized H2O for 20 min prior to visualization using UV light. For densitometry, the mean grey pixels per area unit of the library bands were determined from the agarose gel by using ImageJ version 1.52q [34].

https://doi.org/10.1371/journal.pone.0246919.s001

(DOCX)

S2 Fig. Venn diagram of the number of soil bacterial amplicon sequence variants (ASVs) in paired temperate agroforestry and monoculture cropland systems in three different soil types.

Samples in the agroforestry systems were collected in the tree row as well as at 1 m, 7 m, and 24 m distance from the tree row within the agroforestry crop row on the Calcaric Phaeozem (A), Gleyic Cambisol (B), and Vertic Cambisol soil (C) (n = 4 per soil type × sampling location).

https://doi.org/10.1371/journal.pone.0246919.s002

(DOCX)

S1 Table. Permutational multivariate analysis of variance (PERMANOVA) results.

PERMANOVA was performed with 999 permutations using amplicon sequence variant (ASV) count data. df = degrees of freedom; Sum Sq = sum of squares; R2 = coefficient of determination; F = pseudo-F ratio; p-values marked in bold indicate statistical significance at p < 0.05. a three soil types (Calcaric Phaeozem, Gleyic Cambisol, and Vertic Cambisol). b five sampling locations within each soil type, of which four were located in the agroforestry cropland (tree row, 1 m, 7 m, and 24 m distance from the tree row) and one in the adjacent monoculture cropland (Fig 1).

https://doi.org/10.1371/journal.pone.0246919.s003

(DOCX)

Acknowledgments

The authors would like to thank Guodong Shao and Jie Luo for their help with soil sampling and for providing us with photos of the study sites.

Citation: Beule L, Karlovsky P (2021) Tree rows in temperate agroforestry croplands alter the composition of soil bacterial communities. PLoS ONE 16(2): e0246919. https://doi.org/10.1371/journal.pone.0246919

References

1. Jose S, Gillespie AR, Pallardy SG. Interspecific interactions in temperate agroforestry. Agroforestry Systems. 2004;61: 237–255.

2. Cannell MGR, Noordwijk MV, Ong CK. The central agroforestry hypothesis: the trees must acquire resources that the crop would not otherwise acquire. Agroforest Syst. 1996;34: 27–31.

3. Allen SC, Jose S, Nair PKR, Brecke BJ, Nkedi-Kizza P, Ramsey CL. Safety-net role of tree roots: evidence from a pecan (Carya illinoensis K. Koch)–cotton (Gossypium hirsutum L.) alley cropping system in the southern United States. Forest Ecology and Management. 2004;192: 395–407.

4. Wang Y, Zhang B, Lin L, Zepp H. Agroforestry system reduces subsurface lateral flow and nitrate loss in Jiangxi Province, China. Agriculture, ecosystems & environment. 2011;140: 441–453.

5. Kay S, Crous-Duran J, Ferreiro-Domínguez N, García de Jalón S, Graves A, Moreno G, et al. Spatial similarities between European agroforestry systems and ecosystem services at the landscape scale. Agroforestry Systems. 2018;92: 1075–1089.

6. Torralba M, Fagerholm N, Burgess PJ, Moreno G, Plieninger T. Do European agroforestry systems enhance biodiversity and ecosystem services? A meta-analysis. Agriculture, Ecosystems & Environment. 2016;230: 150–161.

7. Pardon P, Reubens B, Reheul D, Mertens J, De Frenne P, Coussement T, et al. Trees increase soil organic carbon and nutrient availability in temperate agroforestry systems. Agriculture, Ecosystems & Environment. 2017;247: 98–111.

8. Udawatta RP, Rankoth L, Jose S. Agroforestry and Biodiversity. Sustainability. 2019;11: 2879.

9. Pardon P, Reubens B, Mertens J, Verheyen K, De Frenne P, De Smet G, et al. Effects of temperate agroforestry on yield and quality of different arable intercrops. Agricultural Systems. 2018;166: 135–151.

10. Pardon P, Mertens J, Reubens B, Reheul D, Coussement T, Elsen A, et al. Juglans regia (walnut) in temperate arable agroforestry systems: effects on soil characteristics, arthropod diversity and crop yield. Renewable Agriculture and Food Systems. 2020;35: 533–549.

11. Swieter A, Langhof M, Lamerre J, Greef JM. Long-term yields of oilseed rape and winter wheat in a short rotation alley cropping agroforestry system. Agroforest Syst. 2018;93: 1853–1864.

12. Beule L, Lehtsaar E, Rathgeb A, Karlovsky P. Crop Diseases and Mycotoxin Accumulation in Temperate Agroforestry Systems. Sustainability. 2019;11: 2925.

13. Smith J, Pearce BD, Wolfe MS. A European perspective for developing modern multifunctional agroforestry systems for sustainable intensification. Renewable Agriculture and Food Systems. 2012;27: 323–332.

14. Smith J, Pearce BD, Wolfe MS. Reconciling productivity with protection of the environment: Is temperate agroforestry the answer? Renewable Agriculture and Food Systems. 2013;28: 80–92.

15. Lee K-H, Jose S. Soil respiration and microbial biomass in a pecan—cotton alley cropping system in Southern USA. Agroforestry Systems. 2003;58: 45–54. pmid:11099962

16. Mao R, Zeng D-H, Li L-J, Hu Y-L. Changes in labile soil organic matter fractions following land use change from monocropping to poplar-based agroforestry systems in a semiarid region of Northeast China. Environ Monit Assess. 2012;184: 6845–6853. pmid:22124586

17. Beuschel R, Piepho H-P, Joergensen RG, Wachendorf C. Similar spatial patterns of soil quality indicators in three poplar-based silvo-arable alley cropping systems in Germany. Biol Fertil Soils. 2019;55: 1–14.

18. Seiter S, Ingham ER, William RD. Dynamics of soil fungal and bacterial biomass in a temperate climate alley cropping system. Applied Soil Ecology. 1999;12: 139–147.

19. Mungai NW, Motavalli PP, Kremer RJ, Nelson KA. Spatial variation of soil enzyme activities and microbial functional diversity in temperate alley cropping systems. Biology and Fertility of Soils. 2005;42: 129–136.

20. Udawatta RP, Kremer RJ, Adamson BW, Anderson SH. Variations in soil aggregate stability and enzyme activities in a temperate agroforestry practice. Applied Soil Ecology. 2008;39: 153–160.

21. Udawatta RP, Kremer RJ, Garrett HE, Anderson SH. Soil enzyme activities and physical properties in a watershed managed under agroforestry and row-crop systems. Agriculture, Ecosystems & Environment. 2009;131: 98–104.

22. Unger IM, Goyne KW, Kremer RJ, Kennedy AC. Microbial community diversity in agroforestry and grass vegetative filter strips. Agroforest Syst. 2013;87: 395–402.

23. Weerasekara C, Udawatta RP, Jose S, Kremer RJ, Weerasekara C. Soil quality differences in a row-crop watershed with agroforestry and grass buffers. Agroforest Syst. 2016;90: 829–838.

24. Rivest D, Lorente M, Olivier A, Messier C. Soil biochemical properties and microbial resilience in agroforestry systems: effects on wheat growth under controlled drought and flooding conditions. Sci Total Environ. 2013;463–464: 51–60. pmid:24300458

25. Bardhan S, Jose S, Udawatta RP, Fritschi F. Microbial community diversity in a 21-year-old temperate alley cropping system. Agroforest Syst. 2013;87: 1031–1041.

26. Gloor GB, Hummelen R, Macklaim JM, Dickson RJ, Fernandes AD, MacPhee R, et al. Microbiome Profiling by Illumina Sequencing of Combinatorial Sequence-Tagged PCR Products. PLoS One. 2010;5: e15406. pmid:21048977

27. Banerjee S, Baah-Acheamfour M, Carlyle CN, Bissett A, Richardson AE, Siddique T, et al. Determinants of bacterial communities in Canadian agroforestry systems: Co-occurrence patterns of soil bacterial communities. Environmental Microbiology. 2016;18: 1805–1816. pmid:26184386

28. Beule L, Lehtsaar E, Corre MD, Schmidt M, Veldkamp E, Karlovsky P. Poplar Rows in Temperate Agroforestry Croplands Promote Bacteria, Fungi, and Denitrification Genes in Soils. Front Microbiol. 2020;10: 3108. pmid:32038551

29. Beule L, Corre MD, Schmidt M, Göbel L, Veldkamp E, Karlovsky P. Conversion of monoculture cropland and open grassland to agroforestry alters the abundance of soil bacteria, fungi and soil-N-cycling genes. PLOS ONE. 2019;14: e0218779. pmid:31246995

30. Tsonkova P, Böhm C, Quinkenstein A, Freese D. Ecological benefits provided by alley cropping systems for production of woody biomass in the temperate region: a review. Agroforest Syst. 2012;85: 133–152.

31. Brandfass C, Karlovsky P. Upscaled CTAB-Based DNA Extraction and Real-Time PCR Assays for Fusarium culmorum and F. graminearum DNA in Plant Material with Reduced Sampling Error. Int J Mol Sci. 2008;9: 2306–2321. pmid:19330077

32. Guerra V, Beule L, Lehtsaar E, Liao H-L, Karlovsky P. Improved Protocol for DNA Extraction from Subsoils Using Phosphate Lysis Buffer. Microorganisms. 2020;8: 532. pmid:32272709

33. Herlemann DP, Labrenz M, Jürgens K, Bertilsson S, Waniek JJ, Andersson AF. Transitions in bacterial communities along the 2000 km salinity gradient of the Baltic Sea. ISME J. 2011;5: 1571–1579. pmid:21472016

34. Schneider CA, Rasband WS, Eliceiri KW. NIH Image to ImageJ: 25 years of image analysis. Nat Methods. 2012;9: 671–675. pmid:22930834

35. Bolyen E, Rideout JR, Dillon MR, Bokulich NA, Abnet CC, Al-Ghalith GA, et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat Biotechnol. 2019;37: 852–857. pmid:31341288

36. Callahan BJ, McMurdie PJ, Rosen MJ, Han AW, Johnson AJA, Holmes SP. DADA2: High-resolution sample inference from Illumina amplicon data. Nature Methods. 2016;13: 581–583. pmid:27214047

37. Callahan BJ, McMurdie PJ, Holmes SP. Exact sequence variants should replace operational taxonomic units in marker-gene data analysis. ISME J. 2017;11: 2639–2643. pmid:28731476

38. Edgar RC. Updating the 97% identity threshold for 16S ribosomal RNA OTUs. Bioinformatics. 2018;34: 2371–2375. pmid:29506021

39. Quast C, Pruesse E, Yilmaz P, Gerken J, Schweer T, Yarza P, et al. The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Res. 2013;41: D590–D596. pmid:23193283

40. Rognes T, Flouri T, Nichols B, Quince C, Mahé F. VSEARCH: a versatile open source tool for metagenomics. PeerJ. 2016;4: e2584. pmid:27781170

41. Beule L, Karlovsky P. Improved normalization of species count data in ecology by scaling with ranked subsampling (SRS): application to microbial communities. PeerJ. 2020;8: e9593. pmid:32832266

42. R Core Team. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. 2017.

43. Gaby JC, Buckley DH. A comprehensive aligned nifH gene database: a multipurpose tool for studies of nitrogen-fixing bacteria. Database (Oxford). 2014;2014. pmid:24501396

44. Oksanen J, Kindt R, Legendre P, O’Hara B, Simpson G, Solymos P, et al. The vegan package. 2019. Available: https://cran.r-project.org/web/packages/vegan/vegan.pdf

45. Csardi G, Nepusz T. The igraph software package for complex network research. International Journal of Complex Systems. 2006;1695: 1–9.

46. Fruchterman TMJ, Reingold EM. Graph drawing by force-directed placement. Software: Practice and Experience. 1991;21: 1129–1164.

47. Prober SM, Leff JW, Bates ST, Borer ET, Firn J, Harpole WS, et al. Plant diversity predicts beta but not alpha diversity of soil microbes across grasslands worldwide. Ecol Lett. 2015;18: 85–95. pmid:25430889

48. Phillips DA. Efficiency of symbiotic nitrogen fixation in legumes. Annual Review of Plant Physiology. 1980;31: 29–49.

49. Laranjo M, Alexandre A, Oliveira S. Legume growth-promoting rhizobia: An overview on the Mesorhizobium genus. Microbiological Research. 2014;169: 2–17. pmid:24157054

50. Liu S, Coyne MS, Grove JH. Long-term tillage and nitrogen fertilization: Consequences for nitrifier density and activity. Applied Soil Ecology. 2017;120: 121–127.

51. Rousk J, Bååth E, Brookes PC, Lauber CL, Lozupone C, Caporaso JG, et al. Soil bacterial and fungal communities across a pH gradient in an arable soil. ISME J. 2010;4: 1340–1351. pmid:20445636

52. Gomez E, Bisaro V, Conti M. Potential C-source utilization patterns of bacterial communities as influenced by clearing and land use in a vertic soil of Argentina. Applied Soil Ecology. 2000;15: 273–281.

53. Johnson JS, Spakowicz DJ, Hong B-Y, Petersen LM, Demkowicz P, Chen L, et al. Evaluation of 16S rRNA gene sequencing for species and strain-level microbiome analysis. Nat Commun. 2019;10: 1–11. pmid:30602773

54. Young IM, Crawford JW. Interactions and Self-Organization in the Soil-Microbe Complex. Science. 2004;304: 1634–1637. pmid:15192219

55. Haichar F el Z, Marol C, Berge O, Rangel-Castro JI, Prosser JI, Balesdent J, et al. Plant host habitat and root exudates shape soil bacterial community structure. The ISME Journal. 2008;2: 1221–1230. pmid:18754043

56. Schlatter DC, Bakker MG, Bradeen JM, Kinkel LL. Plant community richness and microbial interactions structure bacterial communities in soil. Ecology. 2015;96: 134–142. pmid:26236898

57. Lange M, Eisenhauer N, Sierra CA, Bessler H, Engels C, Griffiths RI, et al. Plant diversity increases soil microbial activity and soil carbon storage. Nature Communications. 2015;6: 6707. pmid:25848862

58. Zhalnina K, Louie KB, Hao Z, Mansoori N, da Rocha UN, Shi S, et al. Dynamic root exudate chemistry and microbial substrate preferences drive patterns in rhizosphere microbial community assembly. Nat Microbiol. 2018;3: 470–480. pmid:29556109

59. Veach AM, Morris R, Yip DZ, Yang ZK, Engle NL, Cregger MA, et al. Rhizosphere microbiomes diverge among Populus trichocarpa plant-host genotypes and chemotypes, but it depends on soil origin. Microbiome. 2019;7: 76. pmid:31103040

60. Voges MJEEE, Bai Y, Schulze-Lefert P, Sattely ES. Plant-derived coumarins shape the composition of an Arabidopsis synthetic root microbiome. PNAS. 2019;116: 12558–12565. pmid:31152139

61. Oelbermann M, Paul Voroney R, Gordon AM. Carbon sequestration in tropical and temperate agroforestry systems: a review with examples from Costa Rica and southern Canada. Agriculture, Ecosystems & Environment. 2004;104: 359–377.

62. Zeng Q, Liu Y, An S. Impact of litter quantity on the soil bacteria community during the decomposition of Quercus wutaishanica litter. PeerJ. 2017;5: e3777. pmid:28894648

63. Fanin N, Hättenschwiler S, Fromin N. Litter fingerprint on microbial biomass, activity, and community structure in the underlying soil. Plant Soil. 2014;379: 79–91.

64. Ramirez KS, Lauber CL, Knight R, Bradford MA, Fierer N. Consistent effects of nitrogen fertilization on soil bacterial communities in contrasting systems. Ecology. 2010;91: 3463–3470. pmid:21302816

65. Mathew RP, Feng Y, Githinji L, Ankumah R, Balkcom KS. Impact of No-Tillage and Conventional Tillage Systems on Soil Microbial Communities. In: Applied and Environmental Soil Science [Internet]. 2012 [cited 9 Apr 2019].

66. Smith CR, Blair PL, Boyd C, Cody B, Hazel A, Hedrick A, et al. Microbial community responses to soil tillage and crop rotation in a corn/soybean agroecosystem. Ecol Evol. 2016;6: 8075–8084. pmid:27878079

67. Yin C, Mueth N, Hulbert S, Schlatter D, Paulitz TC, Schroeder K, et al. Bacterial Communities on Wheat Grown Under Long-Term Conventional Tillage and No-Till in the Pacific Northwest of the United States. Phytobiomes Journal. 2017;1: 83–90.

68. Hill RL. Long-Term Conventional and No-Tillage Effects on Selected Soil Physical Properties. Soil Science Society of America Journal. 1990;54: 161–166.

69. Aird D, Ross MG, Chen W-S, Danielsson M, Fennell T, Russ C, et al. Analyzing and minimizing PCR amplification bias in Illumina sequencing libraries. Genome Biol. 2011;12: R18. pmid:21338519

70. Krehenwinkel H, Wolf M, Lim JY, Rominger AJ, Simison WB, Gillespie RG. Estimating and mitigating amplification bias in qualitative and quantitative arthropod metabarcoding. Sci Rep. 2017;7. pmid:29247210