Abstract

The Po river plain lowland springs represent unique paradigms of managed environments. Their current locations used to be swamps that were drained 6–7 centuries ago, and they have been in constant use ever since. Our aims were to identify the effects of land use on the microbial communities of these soils, look for associated diversity drivers, and assess the applicability of ecology theories with respect to identified patterns. We screened the microbial diversity across a land use transect via high-throughput sequencing of partial 16S rrRNA gene amplicons. Land use had a major effect on soil properties and microbial community structures. Total organic carbon and pH were major diversity drivers for Bacteria, and pH was important for Archaea. We identified the potential contribution of soil amendments to the indigenous microbial communities, and also gained insights into potential roles of taxa in the organic carbon turnover. Verrucomicrobia coincided with the higher values of the recalcitrant organic carbon. Actinobacteria and Acidobacteria correlated with the more labile organic carbon. Finally, the higher diversity found in the soils less enzymatically active and relatively poorer in nutrients, may be explained to an extent by niche-based theories such as the resource heterogeneity hypothesis and Connell's intermediate disturbance hypothesis.

Introduction

Soil is a highly complex and important matrix encompassing immense bio-diversity and a large number of biological processes (Barrios, 2007). Many of these processes contribute to ecosystem services (e.g. nutrient cycling, soil erosion control and biological pest control) supporting human activity in both natural and managed environments. However, human input and related perturbations, from attempts to enhance ecosystem services, alter soil traits close to or even beyond points at which natural attenuation mechanisms may not longer lead to functional restoration (Shennan, 2008).

A characteristic example of human intervention in natural ecosystems is the drainage of wetlands for agricultural use. Particular environments where this strategy was applied in the past are the lowland springs and their surroundings in the northern area of the Po river plains (Kløve et al., 2011), in Italy. Such drastic change in the physico-chemical circumstances may have a strong impact on the structure and functioning of organisms living in the system. This holds in particular for microbes, which, in parallel, are the primary agents responsible for key ecosystem processes. Although the effects of land use and management on the soil microbial community structure have long been acknowledged (Guo & Gifford, 2002; Kelliher et al., 2005; Acosta-Martínez et al., 2008; Wallenius et al., 2011), relevant detailed information is scarce (Garbeva et al., 2004). The vast numbers of microbial inhabitants of complex soil environments (Schloss & Handelsman, 2006) hampered efforts for detailed diversity screening during previous decades. Recent high-throughput sequencing technologies have overcome some methodological problems (Sogin et al., 2006; Roesch et al., 2007; Lauber et al., 2009; Gloor et al., 2010; Bartram et al., 2011; Bates et al., 2011; Delmont et al., 2012; Uroz et al., 2013), therefore increasing the depth of investigation and filling associated gaps. However, more issues need to be addressed to uncover the total existing microbial diversity in targeted, polymerase chain reaction (PCR)-based approaches, e.g. with selected nucleic acid extraction strategies, the primer choices, the PCR setup and sequencing design strategies, rendering the achievement of such a goal still expensive in terms of time and effort (Berry et al., 2011; Delmont et al., 2011; Lombard et al., 2011; Vasileiadis et al., 2012). Therefore, studies using approaches based on a single primer set and/or DNA extraction method, like the presented study, should restrict their claims to the part of the existing diversity that has been uncovered. To this extent, when using the term ‘total’ in conjunction with microbial community-associated terms hereafter, we refer to values and estimates that have been uncovered.

Our aim was to study the microbial communities of soils along a land use transect in the Gaverina (Bergamo, Italy) lowland spring environment (Supporting Information, Fig. S1). Combining the screening depth provided by novel high-throughput sequencing technologies and appropriate data analysis approaches, we explored 16S rRNA gene-based microbial diversity and associations with land use traits, soil chemical and biological properties. Further, we assessed the applicability of existing ecological hypotheses in explaining the patterns of shaping microbial communities in these soils.

Materials and methods

Soil sample physical-chemical-biochemical properties

Historically, all soil environments studied here originated from the bed of a swamp owing its existence to underground water tension. Drainage, which took place 6–7 centuries ago, led to a blossoming of agricultural activity in the area located next to the river Po and generated the land use and management gradients of the soils studied here (Kassen & Rainey, 2004) (Fig. S1: Fonti della Gaverina 45°27′55.69″N, 9°38′20.05″E, elevation 97 m, Italy – sampling carried out in May 2009). The gradients include (1) the low-land spring banks (saturated with water during the high precipitation season), (2) maintained meadow zones around the springs (plant biomass is a mixture of grass, legumes and native plants, with the above-ground biomass being removed twice a year and the plant composition being supported by no tillage sowing of legume seeds when their production decreases – every 3–5 years); and (3) the surrounding fields, cultivated during the sampling period with maize (part of a rotation program with legumes, grasses and maize, receiving cow slurry every 2 years). Top-soil samples (the upper 5–10 cm) were collected in triplicate with an auger from these three environments and soil from the selected depth range was sieved through a 2-mm mesh (soil properties and management traits are provided in Table 1996). Particle-size analysis was carried out using the pipette method (Day, 1965). Three carbon fraction measurements (labile, moderately labile and recalcitrant) were based on the Walkley–Black method as modified elsewhere (Chan et al., 2001) for obtaining the different fractions. To investigate potential links between the carbon fractions and the energy flow associated with the soil microbial communities we measured widespread associated biological activities: β-glucosidase (EC 3.2.1.21) and acid phosphatase (EC 3.1.3.2) activities were determined by the p-nitrophenol method of Eivazi & Tabatabai (1988) and Margesin & Schinner (1994), with p-nitrophenyl-β-d-glucoside and p-nitrophenylphosphate as substrate, respectively. Further, nitrate reductase, being genetically widespread (Philippot, 2002), was determined using KNO3 as a substrate according to Fu & Tabatabai (1989).

Soil properties and qualitative land use and management traits

Management traits Maize Meadow Riparian  
Slurry – Tillage Plant biomass removal Seasonal soil saturation 
Soil perturbation frequency Frequent – Not frequent 
Texture Loam – Clay loam Loam Loam – Clay loam 
Chemical properties AVG SD AVG SD AVG SD F statistic 
pH*** (c) 6.4 ± 0.10 (b) 6.9 ± 0.00 (a) 7.16 ± 0.12 58.4 
CEC (mEq 100 g−11.07 ± 0.18 1.79 ± 0.25 1.1 ± 0.71 2.5 
Humidity (%) 14.57 ± 2.26 18.06 ± 2.28 24.31 ± 9.98 
N (%) 0.17 ± 0.01 0.32 ± 0.03 0.2 ± 0.24 0.9 
TOC (%)*** (b) 1.43 ± 0.06 (a) 2.98 ± 0.26 (c) 0.74 ± 0.04 158 
TOC/N 8.25 ± 0.08 9.34 ± 0.05 7.97 ± 5.57 N/A 
Labile OC proportion*** (a) 0.74 ± 0.10 (a) 0.63 ± 0.01 (b) 0.37 ± 0.03 28.3 
Moderately labile OC proportion 0.05 ± 0.02 0.11 ± 0.07 0.12 ± 0.03 2.1 
Recalcitrant OC proportion** (b) 0.21 ± 0.08 (b) 0.26 ± 0.08 (a) 0.51 ± 0.04 15.6 
Enzymatic activities 
Nitrate reductase (μg N g−1 24 h−159.43 ± 35.68 120.07 ± 156.35 159.2 ± 72.09 0.7 
β-glucosidase (μg PNP g−1 h−1)*** (a) 55.61 ± 10.35 (a) 66.99 ± 6.80 (b) 13.70 ± 5.59 38.4 
Phosphatase (μg PNP g−1 h−1)*** (b) 61.29 ± 0.73 (a) 170.07 ± 43.24 (c) 6.22 ± 2.98 33.3 
Management traits Maize Meadow Riparian  
Slurry – Tillage Plant biomass removal Seasonal soil saturation 
Soil perturbation frequency Frequent – Not frequent 
Texture Loam – Clay loam Loam Loam – Clay loam 
Chemical properties AVG SD AVG SD AVG SD F statistic 
pH*** (c) 6.4 ± 0.10 (b) 6.9 ± 0.00 (a) 7.16 ± 0.12 58.4 
CEC (mEq 100 g−11.07 ± 0.18 1.79 ± 0.25 1.1 ± 0.71 2.5 
Humidity (%) 14.57 ± 2.26 18.06 ± 2.28 24.31 ± 9.98 
N (%) 0.17 ± 0.01 0.32 ± 0.03 0.2 ± 0.24 0.9 
TOC (%)*** (b) 1.43 ± 0.06 (a) 2.98 ± 0.26 (c) 0.74 ± 0.04 158 
TOC/N 8.25 ± 0.08 9.34 ± 0.05 7.97 ± 5.57 N/A 
Labile OC proportion*** (a) 0.74 ± 0.10 (a) 0.63 ± 0.01 (b) 0.37 ± 0.03 28.3 
Moderately labile OC proportion 0.05 ± 0.02 0.11 ± 0.07 0.12 ± 0.03 2.1 
Recalcitrant OC proportion** (b) 0.21 ± 0.08 (b) 0.26 ± 0.08 (a) 0.51 ± 0.04 15.6 
Enzymatic activities 
Nitrate reductase (μg N g−1 24 h−159.43 ± 35.68 120.07 ± 156.35 159.2 ± 72.09 0.7 
β-glucosidase (μg PNP g−1 h−1)*** (a) 55.61 ± 10.35 (a) 66.99 ± 6.80 (b) 13.70 ± 5.59 38.4 
Phosphatase (μg PNP g−1 h−1)*** (b) 61.29 ± 0.73 (a) 170.07 ± 43.24 (c) 6.22 ± 2.98 33.3 

Post-hoc pairwise comparisons: Tukey HSD (α < 0.05). Significant pairwise differences are denoted with different letters in parenthesis.

N/A: conditions were not met and non-parametric test was not significant.

anova significance: *** 0.001, ** 0.01, * 0.05.

DNA isolation, PCR conditions, multiplexing and sequencing

DNA was extracted from 500 mg of each soil sample using the FastDNA® SPIN kit for soil with a FastPrep® 24 instrument (MP Biomedicals, LLC, Solon, OH) according to the manufacturer's instructions. Extracts were quantified using the Quant-iT™ HS dsDNA Assay kit (Invitrogen, Paisley, UK) and the Qubit™ fluorometer (Invitrogen). DNA purity and shearing were tested using a Biophotometer (Eppendorf, Hamburg, Germany) and 0.8% agarose gel, respectively. DNA 2 ng extracts were used for the bacterial 16S rRNA gene PCR amplifications; 20 ng was used for the archaeal ones. The targeted V5 region of the 16S rRNA gene was selected on the basis of the available sequencing technology read length screening abilities (100-bp paired-end reads, related information is provided in the end of the present section) and also conservation of candidate V region flanking loci, as indicated in a previous database search and simulation study (Vasileiadis et al., 2012). PCR 50 μL reactions were performed according to the following program: 94 °C for 5 min, 35 cycles X [94 °C for 30 s of denaturation; 50 °C (for bacterial primer-sets) and 60 °C for 30 s of primer annealing for the bacterial and the archaeal 16S rRNA gene targeting primer-sets, respectively (Table S1); 72 °C for 30 s of elongation] and 72 °C for 10 min. PCR 50 μL reactions were performed using the following mixtures: 1 × PCR buffer, 2.5 mM MgCl2, 2.5 U AmpliTaq® Taq polymerase (Applied Biosystems, Foster City, CA), 0.4 mM of each dNTP, 0.5 μM forward primer and 0.5 μM of reverse primer, 2 and 20 ng respectively of template DNA when bacterial and archaeal 16S rRNA gene targeting primers were used. DNA extracts or PCR product were stored at −20 °C until further use.

Polymerase chain reaction (PCR) products were labelled with 6-bp indices (one for each soil sample, Table S2) as previously described by Meyer et al., (2008) and were concomitantly pooled in a single sample along with other indexed samples not related to the analysed sequences. The sequencing process of the PCR products pool was carried out by Fasteris SA (Geneva, Switzerland), using the TruSeq™ SBS v5 kit (Illumina Inc., San Diego, CA) for the sequencing library generation. The sequencing reaction was performed with a HiSeq 2000 Illumina instrument (Illumina Inc.). The samples were further multiplexed with bulk DNA samples in the same lane using the Illumina barcoded adapters, to avoid Illumina technology-related amplicon sequencing issues (Bartram et al., 2011).

Dataset preparation

The software used for the base-calling and demultiplexing (according to Illumina barcodes) processes was hiseq control soft v1.1.37, rta v1.7.45 and casava v1.7. Sequences were further demultiplexed according to the sample indexes provided in Table S2 and quality-controlled using the Fastx-toolkit (http://hannonlab.cshl.edu/fastx_toolkit/). Good quality forward and reverse amplicon reads that could unambiguously re-construct the amplicons of their origin were used in follow-up analyses. Considering that most V5 region amplicons are shorter than 150 bp (Vasileiadis et al., 2012), and also the method used for sample indexing, the 100-bp length of each of the two obtained reads per amplicon was sufficient for re-constructing the vast amplicon majority. This task was performed with the ‘pandaseq’ script, previously made available by Bartram et al., (2011), with the restrictions of at least 10 bp of overlap between read pairs and no allowed mismatches, resulting in 3 969 194 assembled reads. Further quality control steps resulted in 2 989 302 high quality reads (average and minimum PhredQ values of > 28 and ≥ 20, respectively) with a mean length of c. 162 bp (including index and primer sequences) for all samples (Archaea and Bacteria). These reads were further demultiplexed according to soil sample of origin and targeted microbial group, and their unassembled versions were submitted to the Sequence Read Archive (SRA) database of the National Center for Biotechnology Information (NCBI, http://www.ncbi.nlm.nih.gov/Traces/sra) and are available under accession numbers SRS396685, SRS396687 and SRS396690-705.

Data analysis

Qualified sequence reads were used for operational taxonomic unit (OTU) and taxonomy-supervised analyses using mothur software v1.29 (Schloss et al., 2009). Sequences were pre-processed for analysis by: (1) sequence chimera control with Uchime (Edgar et al., 2011) with the exception of archaeal sequences (due to dependence of the method in existing databases which are still poor for Archaea); (2) removing non specific amplicons that were classified with the Bayesian classifier (Wang et al., 2007) for 50% bootstrap cutoff (Claesson et al., 2009) in non-targeted kingdoms, according to the available mothur version of the ribosomal database project (RDP) database; and (3) down-sampling each sample sequence set to the sample represented by the lowest sequence number per taxonomical kingdom in order to achieve equal sequence representation of samples. This choice was previously shown to reduce sequencing artifact effects on the data analysis (Schloss et al., 2011) and also biases associated with the calculation of several diversity indices (Gihring et al., 2011). All referred steps resulted in analysis of a total of 781 560 (86 840 per sample) sequences for Bacteria and 35 775 (3975 per sample) sequences for Archaea for all samples. In the OTU-based analysis, sequences were aligned with the nast algorithm (DeSantis et al., 2006; Schloss, 2009, 2010) and a kmer search approach vs. kingdom-specific SILVA pre-aligned reference databases as amended by the mothur software developers (Schloss et al., 2009; Quast et al., 2013). Pairwise distances of the aligned sequences and the average linkage algorithm (Sun et al., 2009; Schloss, 2010; Schloss & Westcott, 2011) were used for sequence clustering into OTUs encompassing sequence reads with a distance equal or < 3%. A taxonomy-supervised approach was also followed to evaluate relationships among samples and environmental variables with identified taxa. In this approach, the sequences sampled for equal representation, as described above, were classified using the Bayesian classifier (Wang et al., 2007) with 50% bootstrap cutoff value (Claesson et al., 2009). The greengenes (McDonald et al., 2012) and the silva databases (Quast et al., 2013) were used during classification of the bacterial and archaeal datasets, respectively, due to their coverage of the respective taxa.

Calculated indices related to α-diversity were the Good's diversity coverage estimate (Good, 1953), the observed richness (S), the non-parametric Shannon (npH′) (Chao & Shen, 2003) and Shannon-based Equitability [EH = H′/Hmax = H′/ln(S)] (Sheldon, 1969). Analysis of variance of means (anova) and the Tukey's honestly significant difference (HSD) pairwise comparison test (α < 0.05) were performed to assess significant differences between land use and management with respect to total microbial community diversity indices. The anova assumptions were tested using the Shapiro normality test and the Levene's test of homogeneity of variance. Where anova was not applicable, the Kruskal–Wallis non-parametric test for significant differences estimation and the Nemenyi–Damico–Wolfe–Dunn joint ranking test (for confidence intervals of 99%) were performed, along with a Tukey post-hoc analysis for pairwise comparisons. Total community correlation with the environmental variables was performed with Mantel tests, using the Bray–Curtis distance metric for microbial data and the Spearman r coefficient. Correlations between individual OTUs and between the diversity indices with the measured environmental variables were assessed with pairwise correlations using the Spearman r coefficient. The microbial data matrices were subjected to Hellinger transformation (Legendre & Gallagher, 2001) prior to the performance of the multivariate analyses. To assess total community tendencies we performed unconstrained multivariate analyses such as hierarchical clustering using the average linkage algorithm, principal components analysis, principal coordinates analysis (PCoA) on the Bray–Curtis sample distance matrix and analysis of similarity (anosim) of the Bray–Curtis distances. The effect of environmental variables was assessed with distance-based redundancy analysis (db-RDA; Legendre & Anderson, 1999) on the Bray–Curtis distance estimations of the microbial abundance matrices, and permutation tests (9999 permutations) were carried out to assess significance of results. Prior to beginning db-RDA analysis, the measured environmental variables were normalized using the Box–Cox power law transformations (Osborne, 2010) and the output data values were expressed as deviations from the mean (z-scores) to eliminate the effect of unit and distribution range differences. Significant constraints were selected by the forward selection method (Blanchet et al., 2008). Performance of the db-RDA was preferred over canonical correspondence analysis when the length of the first axis of the performed detrended correspondence analysis was lower than 2, indicating linear responses of microbes in the observed environmental gradients (Lepš & Šmilauer, 2001).

The α-diversity-associated indices were calculated with modules of mothur software v1.28.0 (Schloss et al., 2009), and the rest of the statistical analyses were performed in r v2.14.2 (R Development Core Team, 2011) using the packages Vegan (Dixon, 2003), Hmisc (http://biostat.mc.vanderbilt.edu/wiki/Main/Hmisc) and Scales (http://cran.r-project.org/web/packages/scales/index.html).

Results

Soil properties

Soil pH, TOC and its fractions were the chemical properties that showed statistically significant differences among land use types (Table 0001). Higher average values for most properties were derived from meadow (highest TOC, CEC, N content) or riparian samples (highest recalcitrant OC and pH). The exception was the labile OC proportion, where the highest value was observed for the maize samples. Statistically significant differences were shown for the measured biochemical activities of β-glucosidase and acid phosphatase, showing similar patterns between environments. Both of them had higher values in meadow, followed by maize, with β-glucosidase being significantly different between the riparian and the other two environments, whereas for acid phosphatase the differences were significant between all three tested soil environments. Nitrate reductase, on the other hand, did not show significant differences between environments.

Taxonomy-supervised analysis

About 80% of the bacterial sequences were classified to family level, whereas only slightly over 40% of the archaeal sequences was classified to class level for a 50% bootstrap cutoff with the Bayesian classifier (Fig. S2). The dominant bacterial sequences were Acidobacteria,Actinobacteria,Firmicutes,Proteobacteria and Verrucomicrobia (Fig. 0001a). Db-RDA performance with land use as constraining factor showed significantly different communities among sample groups. Highly correlated taxa within each sample group were fermenters such as Bacilli and Clostridia for the maize group, Rhizobiales for meadow samples, and Planctomycetes,Bacteroidetes and Rhodocyclales for the riparian group. TOC and pH were significantly correlated environmental variables, with the shifts observed for the bacterial communities. Dominant archaeal taxa were Thermoplasmata,Methanobacteria,Methanomicrobia,Thaumarchaeota,Halobacteria and several unclassified Euryarchaeota (Fig. 0001b). There was a significant distinction between land use types for Archaea, with the maize soils being clearly separated from the other two soil groups. Thermoplasmata were more abundant in the maize samples, while sequences in unclassified Euryarchaeota were more abundant in meadow and riparian samples. pH was the sole environmental variable significantly associated with community shifts according to db-RDA results for Archaea.

Relative participation of taxa per sample (upper), db-RDA performance with the land use as constraining factor (middle) and db-RDA performance with variable constraints according to the forward selection approach (lower), for the bacterial (a) and archaeal (b) taxonomical classification datasets. Significance of the results was based on permutation tests (9999 permutations); the percentage of variance explained by each db-RDA model and the significance is provided in the top right corner of each plot; it is analyzed among variables in the low corner; axes percentages provide the proportions of the plotted constrained variance; dash-dotted lines show the standard deviations per soil group.

Relative participation of taxa per sample (upper), db-RDA performance with the land use as constraining factor (middle) and db-RDA performance with variable constraints according to the forward selection approach (lower), for the bacterial (a) and archaeal (b) taxonomical classification datasets. Significance of the results was based on permutation tests (9999 permutations); the percentage of variance explained by each db-RDA model and the significance is provided in the top right corner of each plot; it is analyzed among variables in the low corner; axes percentages provide the proportions of the plotted constrained variance; dash-dotted lines show the standard deviations per soil group.

Correlation of OTU diversity with land use type

The achieved coverage of the identified total diversity according to the dedicated sequencing effort in the data analysis, was over 95% and 80% for the bacterial and archaeal datasets, respectively (Fig. 0002) for the 3% sequence–distance defined OTUs (referred to as OTUs hereafter). α-Diversity screening showed that bacterial communities were more diverse, OTU-rich and even for the riparian and the maize soils compared with the meadow soils. However, statistically significant differences were observed only between the riparian and meadow soils (Fig. 0002a). An inverse relationship between soil groups and indices was observed for Archaea. Riparian soils were still more OTU-rich, diverse and even, but with maize soils having the lowest values (Fig. 0002b). β-Diversity analysis of the sampled soils using unconstrained ordination showed clear differences between the microbial communities of the different land use soil groups (Fig. S3). The OTU datasets demonstrated a lower level of association of the archaeal samples from the same land use groups than was found for the bacterial samples. In both cases the meadow samples had higher within group distances compared with the other soil groups.

Boxplots of the coverage estimate, the observed richness (S), the npH′, and the EH′, for the bacterial (a) and archaeal (b) analysed datasets. Significant differences among soil groups (ma: maize, me: meadow, ri: riparian) are indicated by different letters according to the performed anova and Tukey (HSD) test (P≤ 0.05).

Boxplots of the coverage estimate, the observed richness (S), the npH′, and the EH′, for the bacterial (a) and archaeal (b) analysed datasets. Significant differences among soil groups (ma: maize, me: meadow, ri: riparian) are indicated by different letters according to the performed anova and Tukey (HSD) test (P≤ 0.05).

OTU-based assessment of diversity drivers

Sample distances according to their OTU composition were compared with differences of the measured soil property values. Of the measured environmental variables, Mantel tests showed significant correlations between the total bacterial community shifts and TOC, together with its labile and moderately labile fractions, as well as pH (Table 0002). Furthermore, shifts in the phosphatase and β-glucosidase activities were significantly correlated with the total bacterial community shifts. For the archaeal community matrices, the correlated variables were pH and β-glucosidase. Spearman r correlation coefficient tests between the total community diversity indices and measured variables, showed significantly inverse ranked relationships between the npH′ index and the Shannon Equitability (EH′) with TOC and its labile fraction for Bacteria (Fig. S4). Three measured variables were significantly correlated with three archaeal diversity indices (Fig. S5). The pH and the recalcitrant OC fraction were positively correlated with the S, the npH′ and the EH′, whereas the β-glucosidase activity was negatively correlated with these three indices.

Mantel r-values and test significance (999 permutations) for generated OTU relative abundance matrices

 Bacteria Archaea 
r Significance r Significance 
CEC 0.17  0.23 
Humidity 0.18  0.22 
pH 0.33 0.50 
TOC 0.36 ** 0.25 
OC-labile 0.32 0.20  
OC-moderately-labile 0.23 0.11  
OC-recalcitrant −0.01  −0.07  
0.29 0.27 
Nitrate reductase 0.06  0.03  
Phosphatase 0.34 ** 0.18  
β-glucosidase 0.33 0.47 
 Bacteria Archaea 
r Significance r Significance 
CEC 0.17  0.23 
Humidity 0.18  0.22 
pH 0.33 0.50 
TOC 0.36 ** 0.25 
OC-labile 0.32 0.20  
OC-moderately-labile 0.23 0.11  
OC-recalcitrant −0.01  −0.07  
0.29 0.27 
Nitrate reductase 0.06  0.03  
Phosphatase 0.34 ** 0.18  
β-glucosidase 0.33 0.47 

Alpha value: 0.1, 0.05 (*), 0.01 (**), 0.001 (***).

Individual OTU associations with the environmental variables were examined via the Spearman r coefficient rank correlation tests. For P cutoff values of 0.01, the r coefficient values were above 0.7 or below –0.7 for all correlated OTUs (Fig. S6). Sequences of these OTUs were extracted from each dataset and assigned to taxa (Fig. 0003). Overall, higher proportions of the archaeal genotypes showed correlations with chemical and enzymatic variables (up to c. 52.6%) compared with bacterial genotypes (c. 3.8%). However, correlations yielding high sequence numbers were observed in a narrower set of the considered variables for Archaea (e.g. pH and negatively correlated OTUs with moderately labile OC) than for Bacteria. Although the taxonomy levels used in this analysis were quite broad, some taxa showed, for the majority of correlated sequences, distinct responses to the measured variables.

Taxonomical assignments of sequences residing in the highly correlated OTUs (P≤ 0.01) with measured variables according to the Spearman test, for the bacterial (a) and the archaeal (b) datasets. The sizes of the stacked bars are proportional to the total sequence reads used per dataset.

Taxonomical assignments of sequences residing in the highly correlated OTUs (P≤ 0.01) with measured variables according to the Spearman test, for the bacterial (a) and the archaeal (b) datasets. The sizes of the stacked bars are proportional to the total sequence reads used per dataset.

Bacteria

Several rhizobia were positively correlated with total N. Actinobacteria were correlated with soils having less labile carbon and higher measured phosphatase activity. Acidobacteria were more associated with lower pH values (although there were genotypes with opposite responses) and increased TOC and labile OC. Firmicutes were correlated with the lower humidity, lower pH and less labile OC levels. Verrucomicrobia were positively associated with phosphatase and all carbon fractions except for the moderately labile carbon. Furthermore, several Proteobacteria were correlated with higher pH and CEC levels, whereas the lower values of these variables were dominated by other taxa (e.g. Firmicutes, Acidobacteria, Actinobacteria). The increased phosphatase and β-glucosidase activities were mostly correlated with Acidobacteria,Actinobacteria,Verrucomicrobia and some proteobacterial taxa. High levels of nitrate reductase activity coincided with some Bacteroidetes and proteobacterial taxa.

Archaea

Methanobacteria,Methanomicrobia and Thermoplasmata were abundant in soils with a lower pH. Methanobacteria were well correlated with total and labile OC and also phosphatase and β-glucosidase enzymatic activities. Thaumarchaeota were positively associated with the measured nitrate reductase activity: a higher unclassified euryarchaeotal sequence incidence was found in soil samples with higher pH and less labile OC. Sequences classified as Halobacteria showed a preference for soils with higher pH values.

Discussion

The land use differences between the three studied soil groups had profound effects on the biogeochemical properties of soil and its microbial community structures. We have used both the OTU and the taxonomy-supervised analysis in an attempt to fill methodological gaps (Schloss & Westcott, 2011; Sul et al., 2011) and exploit the information potentials of each of these approaches with respect to our objectives.

The three land use types led to the development of distinct microbial communities for Bacteria and Archaea, as supported by both the taxonomy and the OTU-based approaches (Fig. 2008; Rousk et al., 2009; Kuramae et al., 2011; Lopez-Sangil et al., 2011). Similarly, factors significantly correlated with total bacterial or archaeal community shifts were identified in these variables, but with differences concerning taxon responses. Total Bacteria were correlated with TOC and pH gradients according to the taxonomy-supervised approach (Fig. 2011). Previously, the participation of archaeal members was estimated to be between 0% and 10% of that of the prokaryotic members (Bates et al., 2011), whereas pH was indicated to affect this balance drastically by causing reduction of archaeal member abundance with lower values of pH (Bengtson et al., 2012). The latter might also contribute to the increased difference of the archaeal maize communities from the other two soil environments compared with the difference between meadow and riparian soils. The corresponding pH differences were ≥ 0.5 between maize and meadow or riparian areas, whereas the meadow–riparian difference was 0.22. In the total bacterial community, results did not demonstrate equal dependence on the pH changes, showing distinct communities between all three soil groups. Furthermore, the previously identified abundance ratios (Bates et al., 2011) can explain the lower percentage of bacterial sequences that participated in the responsive OTUs, compared with Archaea (see Fig. 0003) and also the total community distance-based correlation with variables (Fig. 0001, Table 0002).

An interesting finding of this study is the high abundance of halobacterial assignments (24.5% of total) in the archaeal sequence dataset (Fig. 2004; Oren & Ventosa, 2009; Oxley et al., 2010).

In an attempt to identify the diversity drivers in the studied soils, we performed modeling of the microbial responses to environmental variables with db-RDA and correlation tests of the variables with individual OTUs. Results showed that maize soils were rich in microbial taxa commonly found in animal excrement (Wang et al., 2010) such as Clostridia and Bacilli and the euryarchaeotal Thermoplasmata (Figs 2003; Snell-Castro et al., 2005; Mao et al., 2011). The distinction of these taxa is also associated to the depth of our strategy, as previous works using increased depth approaches of their time, such as PCR-single-strand-conformation-polymorphism analysis (Peu et al., 2006), failed to detect the associated taxa in soil shortly after amendment with manure. A high rhizobial incidence in meadow is consistent with the presence of leguminous plants among the meadow plants (Welbaum et al., 2004), whereas an increased presence of Planctomycetes and Rhodobacterales, the Alphaproteobacteria known to thrive in freshwater environments, was observed in riparian soils (Crump et al., 1999; Schlesner et al., 2004; Imhoff, 2005; Buesing et al., 2009; Fuerst & Sagulenko, 2011). The significantly increased presence of these taxa along with the periodic exposure to water flow in the riparian soils, indicate more open and dynamic ecosystems compared with the other soils. With respect to bacterial diversity drivers, potential roles have been suggested by the results of the individual OTU–environmental variable correlation tests (Fig. 2009; Rui et al., 2009; Deangelis et al., 2011; Martinez-Garcia et al., 2012) have shown associations with different OC fractions. This could suggest a niche sharing concerning organic matter decomposition between these taxa, with Acidobacteria and Actinobacteria dominating the less recalcitrant OC decomposition stages, whereas the more versatile taxon of Verrucomicrobia, some of whose members are able to decompose more complex organic matter (Martinez-Garcia et al., 2012) and participate in methylotrophy (Chistoserdova et al., 2009), are positively correlated with the recalcitrant OC and (in a lesser degree) with the labile OC.

Another interesting finding was the reduced bacterial richness, diversity and evenness found in meadow soils as compared with the other soil environments. Such effects, when identified during fallow periods in crop rotation, have been interpreted as an outcome of resources depletion (Welbaum et al., 2004). However, measured enzymatic activities related to energy flow in soil systems (β-glucosidase, acid phosphatase), had higher values in meadow samples, whereas close relationships were found in meadow soils with root-associated (Rhizobiales) or complex organic molecule decomposing taxa (Verrucomicrobia,Actinobacteria) (Figs 2006). Therefore, rhizodeposits may favor a restricted number of microbes more efficient in metabolizing these compounds. Such a phenomenon is compatible with the previously proposed hump-shaped relationship of microbial diversity in response to the increase of dominant resources according to the resource heterogeneity hypothesis (Lynch et al., 2004). Thus, a reduced variety in resources may lead to a corresponding reduced number of niches, and therefore to reduced microbial diversity, as also identified in the meadow samples.

Another potential explanation of the identified diversity patterns is related to underlying disturbance phenomena. Disturbance effects are difficult to assess in field experiments because they are usually confounded with other factors contributing to measured variables, and a more disturbance-focused study strategy may be necessary before making associated claims. However, observed diversity indices (Fig. 1978) proposing that perturbation events, when not deleterious for the majority of existing populations, reduce dominance and promote diversity. Originally, the theory described the increase in under-canopy plant diversity after an intermediate disturbance, e.g. fire incidents with reduced effects on plant propagative material, due to post-disturbance lack of competition for resources. The proof-of-principle concerning the applicability of the theory in the microbial world was provided by a laboratory microcosm setup (Buckling et al., 2000). In that study, a unimodal response of pseudomonad population diversity according to colony phenotypes was demonstrated when static broth cultures were subjected to different disturbance frequencies. Similarly, soil bacteria with more autochthonous lifestyles may become more competitive when disturbance events affect the relative abundances of dominant populations, therefore contributing to a diversity increase as in the reported cases of the maize and riparian soils when compared with the less perturbed meadow soils.

Concluding remarks

In the present study we attempted an in-depth per-sample analysis of the microbial communities under the selective pressures associated with different land uses. The chemical and microbiological patterns identified, strongly suggest that land use and management decisions can have a profound effect on chemical and biological soil properties. In the soils studied, bacterial primary producers were more affected by related interventions in soils than were their archaeal counterparts. pH was a principal diversity driver for both Bacteria and Archaea, the latter being more dependent on it. Finally, we speculate that the increased diversity found in the energy-wise poorer and the more perturbed environments may suggest that soil microbial diversity is to a significant part regulated in line with niche-based theories, such as the resource heterogeneity hypothesis (Lynch et al., 2004) or Connell's intermediate disturbance hypothesis (1978).

Acknowledgements

This study was financially supported by the European Community 7th Framework Project GENESIS (226536) on groundwater systems; the Cariplo Foundation under project No. 2011-1088 ‘SNAC- Synthetic and Natural Agrochemical Compounds: ecological impacts on the soil system and effects on plant production’; and the Doctoral School on the Agro-Food System (Agrisystem) of the Università Cattolica del Sacro Cuore (Italy). The authors are thankful for the critical evaluation of the manuscript by Dr Ioannes Stergiopoulos (Department of Plant Pathology of the University of California, Davis, CA) and Dr Eiko E. Kuramae (Microbial Ecology group of the Netherlands Institute of Ecology – NIOO-KNAW, Wageningnen, Netherlands). Moreover we would like to acknowledge the assistance during the study setup by Dr Stefano Longo and the technical support of Giulia Maria Parisi and Pierluisa Fantini.

References

Acosta-Martínez
V
Dowd
S
Sun
Y
Allen
V
(
2008
)
Tag-encoded pyrosequencing analysis of bacterial diversity in a single soil type as affected by management and land use
.
Soil Biol Biochem
 
40
:
2762
2770
.
Bais
HP
Weir
TL
Perry
LG
Gilroy
S
Vivanco
JM
(
2006
)
The role of root exudates in rhizosphere interactions with plants and other organisms
.
Annu Rev Plant Biol
 
57
:
233
266
.
Barrios
E
(
2007
)
Soil biota, ecosystem services and land productivity
.
Ecol Econ
 
64
:
269
285
.
Bartram
AK
Lynch
MDJ
Stearns
JC
Moreno-Hagelsieb
G
Neufeld
JD
(
2011
)
Generation of multi-million 16S rRNA gene libraries from complex microbial communities by assembling paired-end Illumina reads
.
Appl Environ Microbiol
 
77
:
3846
3852
.
Bates
ST
Berg-Lyons
D
Caporaso
JG
Walters
WA
Knight
R
Fierer
N
(
2011
)
Examining the global distribution of dominant archaeal populations in soil
.
ISME J
 
5
:
908
917
.
Bengtson
P
Sterngren
AE
Rousk
J
(
2012
)
Archaeal abundance across a pH gradient in an arable soil and its relationship with bacterial and fungal growth rates
.
Appl Environ Microbiol
 
78
:
5906
5911
.
Berry
D
Mahfoudh
KB
Wagner
M
Loy
A
(
2011
)
Barcoded primers used in multiplex amplicon pyrosequencing bias amplification
.
Appl Environ Microbiol
 
77
:
7846
7849
.
Blanchet
FG
Legendre
P
Borcard
D
(
2008
)
Forward selection of explanatory variables
.
Ecology
 
89
:
2623
2632
.
Buckling
A
Kassen
R
Bell
G
Rainey
PB
(
2000
)
Disturbance and diversity in experimental microcosms
.
Nature
 
408
:
961
964
.
Buesing
N
Filippini
M
Bürgmann
H
Gessner
MO
(
2009
)
Microbial communities in contrasting freshwater marsh microhabitats
.
FEMS Microbiol Ecol
 
69
:
84
97
.
Cavicchioli
R
(
2011
)
Archaea, a timeline of the third domain
.
Nat Rev Microbiol
 
9
:
51
61
.
Chan
KY
Bowman
A
Oates
A
(
2001
)
Oxidizible organic carbon fractions and soil quality changes in an oxic paleustalf under different pasture leys
.
Soil Sci
 
166
:
61
67
.
Chao
A
Shen
TJ
(
2003
)
Nonparametric estimation of Shannon's index of diversity when there are unseen species in sample
.
Environ Ecol Stat
 
10
:
429
443
.
Chistoserdova
L
Kalyuzhnaya
MG
Lidstrom
ME
(
2009
)
The expanding world of methylotrophic metabolism
.
Annu Rev Microbiol
 
63
:
477
499
.
Claesson
MJ
O'Sullivan
O
Wang
Q
et al
(
2009
)
Comparative analysis of pyrosequencing and a phylogenetic microarray for exploring microbial community structures in the human distal intestine
.
PLoS ONE
 
4
:
e6669
.
Connell
JH
(
1978
)
Diversity in tropical rain forests and coral reefs
.
Science
 
199
:
1302
1310
.
Cotta
MA
Whitehead
TR
Zeltwanger
RL
(
2003
)
Isolation, characterization and comparison of bacteria from swine faeces and manure storage pits
.
Environ Microbiol
 
5
:
737
745
.
Crump
BC
Armbrust
EV
Baross
JA
(
1999
)
Phylogenetic analysis of particle-attached and free-living bacterial communities in the Columbia River, its estuary, and the adjacent coastal ocean
.
Appl Environ Microbiol
 
65
:
3192
3204
.
Day
PR
(
1965
)
Particle fractionation and particle-size analysis
.
Methods of Soil Analysis
  (
Black
CA
, ed), pp.
545
567
.
ASA
,
Madison
.
Deangelis
KM
Allgaier
M
Chavarria
Y
et al
(
2011
)
Characterization of trapped lignin-degrading microbes in tropical forest soil
.
PLoS ONE
 
6
:
e19306
.
Delmont
TO
Robe
P
Cecillon
S
et al
(
2011
)
Accessing the soil metagenome for studies of microbial diversity
.
Appl Environ Microbiol
 
77
:
1315
1324
.
Delmont
TO
Prestat
E
Keegan
KP
et al
(
2012
)
Structure, fluctuation and magnitude of a natural grassland soil metagenome
.
ISME J
 
6
:
1677
1687
.
DeSantis
TZ
Hugenholtz
P
Keller
K
et al
(
2006
)
NAST: a multiple sequence alignment server for comparative analysis of 16S rRNA genes
.
Nucleic Acids Res
 
34
:
W394
W399
.
Dixon
P
(
2003
)
VEGAN, a package of R functions for community ecology
.
J Veg Sci
 
14
:
927
930
.
Edgar
RC
Haas
BJ
Clemente
JC
Quince
C
Knight
R
(
2011
)
UCHIME improves sensitivity and speed of chimera detection
.
Bioinformatics
 
27
:
2194
2200
.
Eivazi
F
Tabatabai
A
(
1988
)
Glucosidases and galactosidases in soils
.
Soil Biol Biochem
 
20
:
601
606
.
Fan
TWM
Bird
JA
Brodie
EL
Lane
AN
(
2009
)
C-13-Isotopomer-based metabolomics of microbial groups isolated from two forest soils
.
Metabolomics
 
5
:
108
122
.
Fu
MH
Tabatabai
MA
(
1989
)
Nitrate reductase-activity in soils – effects of trace-elements
.
Soil Biol Biochem
 
21
:
943
946
.
Fuerst
JA
Sagulenko
E
(
2011
)
Beyond the bacterium: Planctomycetes challenge our concepts of microbial structure and function
.
Nat Rev Microbiol
 
9
:
403
413
.
Garbeva
P
van Veen
JA
van Elsas
JD
(
2004
)
Microbial diversity in soil: selection of microbial populations by plant and soil type and implications for disease suppressiveness
.
Annu Rev Phytopathol
 
42
:
243
270
.
Gihring
TM
Green
SJ
Schadt
CW
(
2011
)
Massively parallel rRNA gene sequencing exacerbates the potential for biased community diversity comparisons due to variable library sizes
.
Environ Microbiol
 
14
:
285
290
.
Gloor
GB
Hummelen
R
Macklaim
JM
Dickson
RJ
Fernandes
AD
MacPhee
R
Reid
G
(
2010
)
Microbiome profiling by Illumina sequencing of combinatorial sequence-tagged PCR products
.
PLoS ONE
 
5
:
e15406
.
Good
IJ
(
1953
)
The population frequencies of species and the estimation of population parameters
.
Biometrika
 
40
:
237
264
.
Guo
LB
Gifford
RM
(
2002
)
Soil carbon stocks and land use change: a meta analysis
.
Glob Change Biol
 
8
:
345
360
.
Imhoff
J
(
2005
)
ORDER III. Rhodobacterales ord. nov
.
Bergey's Manual of Systematic Bacteriology
 ,
1
(
Brenner
DJ
Krieg
NR
Garrity
GM
et al., eds), pp.
161
167
.
Springer
,
New York
.
Kassen
R
Rainey
PB
(
2004
)
The ecology and genetics of microbial diversity
.
Annu Rev Microbiol
 
58
:
207
231
.
Kelliher
FM
Barbour
MM
Hunt
JE
(
2005
)
Sucrose application, soil microbial respiration and evolved carbon dioxide isotope enrichment under contrasting land uses
.
Plant Soil
 
268
:
233
242
.
Kløve
B
Ala-aho
P
Bertrand
G
et al
(
2011
)
Groundwater dependent ecosystems. Part I: hydroecological status and trends
.
Environ Sci Policy
 
14
:
770
781
.
Kuramae
E
Gamper
H
van Veen
J
Kowalchuk
G
(
2011
)
Soil and plant factors driving the community of soil-borne microorganisms across chronosequences of secondary succession of chalk grasslands with neutral pH
.
FEMS Microbiol Ecol
 
77
:
285
294
.
Lauber
CL
Hamady
M
Knight
R
Fierer
N
(
2009
)
Pyrosequencing-based assessment of soil pH as a predictor of soil bacterial community structure at the continental scale
.
Appl Environ Microbiol
 
75
:
5111
5120
.
Legendre
P
Anderson
MJ
(
1999
)
Distance-based redundancy analysis: testing multispecies responses in multifactorial ecological experiments
.
Ecol Monogr
 
69
:
1
24
.
Legendre
P
Gallagher
E
(
2001
)
Ecologically meaningful transformations for ordination of species data
.
Oecologia
 
129
:
271
280
.
Lepš
J
Šmilauer
P
(
2001
)
Multivariate Analysis of Ecological Data using CANOCO
 .
Cambridge Press
,
Cambridge
.
Lombard
N
Prestat
E
van Elsas
JD
Simonet
P
(
2011
)
Soil-specific limitations for access and analysis of soil microbial communities by metagenomics
.
FEMS Microbiol Ecol
 
78
:
31
49
.
Lopez-Sangil
L
Rousk
J
Wallander
H
Casals
P
(
2011
)
Microbial growth rate measurements reveal that land-use abandonment promotes a fungal dominance of SOM decomposition in grazed Mediterranean ecosystems
.
Biol Fertil Soils
 
47
:
129
138
.
Lynch
JM
Benedetti
A
Insam
H
Nuti
MP
Smalla
K
Torsvik
V
Nannipieri
P
(
2004
)
Microbial diversity in soil: ecological theories, the contribution of molecular techniques and the impact of transgenic plants and transgenic microorganisms
.
Biol Fertil Soils
 
40
:
363
385
.
Mao
SY
Yang
CF
Zhu
WY
(
2011
)
Phylogenetic analysis of methanogens in the pig feces
.
Curr Microbiol
 
62
:
1386
1389
.
Margesin
R
Schinner
F
(
1994
)
Phosphomonoesterase, phosphodiesterase, phosphotriesterase, and inorganic pyrophosphatase activities in forest soils in an alpine area – effect of pH on enzyme-activity and extractability
.
Biol Fertil Soils
 
18
:
320
326
.
Martinez-Garcia
M
Brazel
DM
Swan
BK
et al
(
2012
)
Capturing single cell genomes of active polysaccharide degraders: an unexpected contribution of Verrucomicrobia
.
PLoS ONE
 
7
:
e35314
.
McDonald
D
Price
MN
Goodrich
J
et al
(
2012
)
An improved Greengenes taxonomy with explicit ranks for ecological and evolutionary analyses of bacteria and archaea
.
ISME J
 
6
:
610
618
.
Meyer
M
Stenzel
U
Hofreiter
M
(
2008
)
Parallel tagged sequencing on the 454 platform
.
Nat Protoc
 
3
:
267
278
.
Oren
A
Ventosa
A
(
2009
)
Halobacterium salinarum strain MMD047 – a low-salt adapted member of the Halobacteriaceae?
Biotechnol Bioprocess Eng
 
14
:
869
870
.
Osborne
JW
(
2010
)
Improving your data transformations: applying the Box–Cox transformation
.
Pract Ass Res Eval
 
15
:
1
9
.
Oxley
APA
Lanfranconi
MP
Würdemann
D
et al
(
2010
)
Halophilic archaea in the human intestinal mucosa
.
Environ Microbiol
 
12
:
2398
2410
.
Peu
P
Brugere
H
Pourcher
A-M
Kerouredan
M
Godon
J-J
Delgenes
J-P
Dabert
P
(
2006
)
Dynamics of a pig slurry microbial community during anaerobic storage and management
.
Appl Environ Microbiol
 
72
:
3578
3585
.
Philippot
L
(
2002
)
Denitrifying genes in bacterial and archaeal genomes
.
BBA-Gene Struct Expr
 
1577
:
355
376
.
Purdy
KJ
Cresswell-Maynard
TD
Nedwell
DB
McGenity
TJ
Grant
WD
Timmis
KN
Embley
TM
(
2004
)
Isolation of haloarchaea that grow at low salinities
.
Environ Microbiol
 
6
:
591
595
.
Quast
C
Pruesse
E
Yilmaz
P
et al
(
2013
)
The SILVA ribosomal RNA gene database project: improved data processing and web-based tools
.
Nucleic Acids Res
 
41
:
D590
D596
.
R Development Core Team
2011
R: A language and environment for statistical computing, reference index version 2.14.1. R Foundation for Statistical Computing
 .
Roesch
LF
Fulthorpe
RR
Riva
A
et al
(
2007
)
Pyrosequencing enumerates and contrasts soil microbial diversity
.
ISME J
 
1
:
283
290
.
Rousk
J
Brookes
PC
Baath
E
(
2009
)
Contrasting soil pH effects on fungal and bacterial growth suggests functional redundancy in carbon mineralisation
.
Appl Environ Microbiol
 
75
:
1589
1596
.
Rui
J
Peng
J
Lu
Y
(
2009
)
Succession of bacterial populations during plant residue decomposition in rice field soil
.
Appl Environ Microbiol
 
75
:
4879
4886
.
Schlesner
H
Rensmann
C
Tindall
BJ
Gade
D
Rabus
R
Pfeiffer
S
Hirsch
P
(
2004
)
Taxonomic heterogeneity within the Planctomycetales as derived by DNA-DNA hybridization, description of Rhodopirellula baltica gen. nov., sp. nov., transfer of Pirellula marina to the genus Blastopirellula gen. nov. as Blastopirellula marina comb. nov. and emended description of the genus Pirellula
.
Int J Syst Evol Microbiol
 
54
:
1567
1580
.
Schloss
PD
(
2009
)
A high-throughput DNA sequence aligner for microbial ecology studies
.
PLoS ONE
 
4
:
e8230
.
Schloss
PD
(
2010
)
The effects of alignment quality, distance calculation method, sequence filtering, and region on the analysis of 16S rRNA gene-based studies
.
PLoS Comput Biol
 
6
:
e1000844
.
Schloss
PD
Handelsman
J
(
2006
)
Toward a census of bacteria in soil
.
PLoS Comput Biol
 
2
:
786
793
.
Schloss
PD
Westcott
SL
(
2011
)
Assessing and improving methods used in OTU-based approaches for 16S rRNA gene sequence analysis
.
Appl Environ Microbiol
 
77
:
3219
3226
.
Schloss
PD
Westcott
SL
Ryabin
T
et al
(
2009
)
Introducing mothur: open-source, platform-independent, community-supported software for describing and comparing microbial communities
.
Appl Environ Microbiol
 
75
:
7537
7541
.
Schloss
PD
Gevers
D
Westcott
SL
(
2011
)
Reducing the effects of PCR amplification and sequencing artifacts on 16S rRNA-based studies
.
PLoS ONE
 
6
:
e27310
.
Sheldon
AL
(
1969
)
Equitability indices: dependence on the species count
.
Ecology
 
50
:
466
467
.
Shennan
C
(
2008
)
Biotic interactions, ecological knowledge and agriculture
.
Philos Trans R Soc Lond B Biol Sci
 
363
:
717
739
.
Sinsabaugh
RL
Lauber
CL
Weintraub
MN
et al
(
2008
)
Stoichiometry of soil enzyme activity at global scale
.
Ecol Lett
 
11
:
1252
1264
.
Snell-Castro
R
Godon
J-J
Delgenès
J-P
Dabert
P
(
2005
)
Characterisation of the microbial diversity in a pig manure storage pit using small subunit rDNA sequence analysis
.
FEMS Microbiol Ecol
 
52
:
229
242
.
Sogin
ML
Morrison
HG
Huber
JA
et al
(
2006
)
Microbial diversity in the deep sea and the underexplored ‘rare biosphere’
.
P Natl Acad Sci USA
 
103
:
12115
12120
.
Sparks
DL
ed (
1996
)
Methods of Soil Analysis: Chemical Methods
 .
Soil Science Society of America
,
Madison
.
Sul
WJ
Cole
JR
Jesus E da
C
Wang
Q
Farris
RJ
Fish
JA
Tiedje
JM
(
2011
)
Bacterial community comparisons by taxonomy-supervised analysis independent of sequence alignment and clustering
.
P Natl Acad Sci USA
 
108
:
14637
14642
.
Sun
Y
Cai
Y
Liu
L
Yu
F
Farrell
ML
McKendree
W
Farmerie
W
(
2009
)
ESPRIT: estimating species richness using large collections of 16S rRNA pyrosequences
.
Nucleic Acids Res
 
37
:
e76
.
Uroz
S
Ioannidis
P
Lengelle
J
Al
Cébron
Morin
E
Buée
M
Martin
F
(
2013
)
Functional Assays and metagenomic analyses reveals differences between the microbial communities inhabiting the soil horizons of a Norway spruce plantation
.
PLoS ONE
 
8
:
e55929
.
Vasileiadis
S
Puglisi
E
Arena
M
Cappa
F
Cocconcelli
PS
Trevisan
M
(
2012
)
Soil bacterial diversity screening using single 16S rRNA gene V regions coupled with multi-million read generating sequencing technologies
.
PLoS ONE
 
7
:
e42671
.
Wallenius
K
Rita
H
Mikkonen
A
et al
(
2011
)
Effects of land use on the level, variation and spatial structure of soil enzyme activities and bacterial communities
.
Soil Biol Biochem
 
43
:
1464
1473
.
Wang
Q
Garrity
GM
Tiedje
JM
Cole
JR
(
2007
)
Naive Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy
.
Appl Environ Microbiol
 
73
:
5261
5267
.
Wang
H
Tolvanen
K
Lehtomäki
A
Puhakka
J
Rintala
J
(
2010
)
Microbial community structure in anaerobic co-digestion of grass silage and cow manure in a laboratory continuously stirred tank reactor
.
Biodegradation
 
21
:
135
146
.
Welbaum
GE
Sturz
AV
Dong
Z
Nowak
J
(
2004
)
Managing soil microorganisms to improve productivity of agro-ecosystems
.
CRC Crit Rev Plant Sci
 
23
:
175
193
.

Supporting Information

Additional Supporting Information may be found in the online version of the article:

Table S1. Primer sequences used in the present study. References correspond to Materials section reference list.

Table S2. Adapters used in the present study for sample indexing of PCR products provided in reference Meyer et al., (2008).

Fig. S1. Satellite photograph of the Gaverina area lowland spring (Source: ‘Gaverina’. 45°2757.32″-69″ North, 9°3824.21″-27.20″ East, Google Earth, March 30, 2012).

Fig. S2. Percentages of classified sequences per taxon level with the Bayesian classifier for 50% bootstrap cutoff, of the bacterial (A) and archaeal (B) sequence data.

Fig. S3. Unconstrained multivariate analyses performed on the bacterial (A) and archaeal (B) OTU data matrices.

Fig. S4. Significant Spearman r correlations (P ≤ 0.05) between bacterial dataset diversity indices and measured environmental variables.

Fig. S5. Significant Spearman r correlations (P ≤ 0.05) between archaeal dataset diversity indices and measured environmental variables.

Fig. S6. Spearman r values for significant correlations (P ≤ 0.01) between individual OTUs and measured environmental variables for the bacterial (A) and archaeal (B) datasets.

Author notes

Editor: Cindy Nakatsu