High-throughput metabarcoding reveals the effect of physicochemical soil properties on soil and litter biodiversity and community turnover across Amazonia

Sampling design and localities

We sampled four localities along the Amazon River (Fig. 1) following the sampling design of Tedersoo et al. (2014). Detailed locality descriptions are available in Ritter et al. (2018). Briefly, we sampled all depths of the litter layer above the mineral soil (all organic matter, including leaves, roots, and animal debris) and the top 5 cm of the mineral soil in a total of 39 circular plots, each with a radius of 28 m. We chose 20 random trees inside each plot and collected litter and soil on both sides of each tree. We then pooled the samples by substrate to obtain one soil and one litter sample per plot. Each sample was stored in a plastic bag with the same weight than soil/litter samples of sterilized white silica gel of 1–4 mm grain size, pre-treated for two minutes of microwave heating (800 W) and 15 min of UV light. The bags were stored at room temperature (around 30 °C) in a dark box to avoid exposure to light. Once they arrived in Sweden, a period between 7–30 days, the samples were frozen (−20 °C). We sampled in different habitat types, which can be summarized as terra-firme, várzeas, igapós, and campinas. These are four of the commonly recognised main Amazonian environments. Terra-firme is characterised by not being inundated during the annual flood season, and terra-firme forests generally have tall stature and complex canopy structure (IBGE, 2004). In contrast, várzeas and igapós are seasonally flooded and remain submerged during parts of the year. Várzeas are flooded by white-water rivers, which carry a large load of suspended sediments, whereas igapós are flooded by black-water rivers, which bring a small load of suspended sediments but a high concentration of organic acids (Junk et al., 2011). Finally, campinas have nutrient-impoverished sandy soils and forests with a reduced stature and relatively simple canopy structure (Prance, 1996; IBGE, 2004).

Our sampling was carried out in four areas: Benjamin Constant (9 plots covering terra-firme, várzea and igapó), our westernmost locality, approximately 1,100 km west of Manaus in the upper Amazonas River (4.383°S, 70.017°W; Fig. 1A); Jaú national park (6 plots covering terra-firme and igapó; 1.850°S, 61.616°W; Fig. 1B) and Novo Airão (3 plots covering campinas; 2.620°S, 60.944°W; Fig. 1C), on the west side of the Negro River; Reserva do Cuieras (six plots covering terra-firme and igapó; 2.609°S, 60.217°W; Fig. 1D) and Reserva da Campina (three plots covering campinas; 2.592°S, 60.030°W; Fig. 1E), on the east side of the Negro River; and Caxiaunã (12 plots covering campinas, terra-firme, várzea, and igapó), a national forest located 350 km west of Belém (1.7352°S, 51,463°W; Fig. 1F), which constitutes our easternmost locality. The sample collection was authorized by Brazilian authorities: ICMBio (registration number 48185-2) and IBAMA (registration number 127341).

Physicochemical soil analyses

We determined the physicochemical soil properties of each plot from three soil samples per plot, totalling 117 samples. The pH was measured in water (soil: water ratio 1:2.5). The exchangeable concentrations were measured for sodium (Na), potassium (K), and phosphorus (P) using Mehlich-1 extraction (unit mg/dm³) and for calcium and magnesium (Ca, Mg) using KCl (1 mol/L) extraction (unit cmol_c/dm³). The sum of all exchangeable bases (SB; which comprises K⁺, Ca²⁺, Mg²⁺, and Na⁺; unit cmol_c/dm³) was then calculated. We also estimated exchangeable aluminium (Al and H+Al unit cmol_c/dm³) extracted with calcium acetate (0.5 mol/L at pH 7.0), aluminium saturation index (m; unit %), and Base Saturation Index (V; unit %). The effective cation exchange capacity (t) as well as the cation exchange capacity (T) were measured at pH 7.0 (unit cmol_c/dm³). The organic matter (O.M) was quantified (unit g/kg) and the organic carbon (C) was estimated from the organic matter as: C = O.M/1.724—Walkley-Black (unit g/kg). Soil texture was characterized by quantifying the fractions of clay (<0.002 mm), silt (0.002–0.05 mm), fine sand (0.05–0.2 mm), coarse sand (0.2–2 mm), and total sand (0.05–2 mm) (unit % of soil weight). We did not quantify nitrogen levels due to the highly volatile nature of nitrogen; its concentration changes quickly during sample storage due to the activity of soil bacteria, and freezing the samples in our remote sampling localities was not feasible. All analyses were commissioned from EMBRAPA Ocidental (Brazil), following the protocol described in Donagema et al. (2011). Afterwards, we used the mean of the three soil samples from the same plot to obtain a representative value for the measurement of each variable for each plot.

DNA extraction, amplification, and sequencing

The detailed laboratory procedures are described in Ritter et al. (2018). Briefly, we extracted soil and litter using the PowerMax® Soil DNA Isolation Kit (MO BIO Laboratories, Carlsbad, CA, USA) following the manufacturer’s instructions. The amplification of 16S was performed by Macrogen (Republic of Korea) following standard protocols, and sequencing was performed using the Illumina MiSeq 2×300 platform. For metabarcoding of the 18S gene, sequencing preparation was performed at the laboratory of the University of Gothenburg as described in Ritter et al. (2018) and the amplicons were sequenced at SciLifeLab (Stockholm, Sweden) using an Illumina MiSeq 2×250 machine.

Sequence analyses

We used the USEARCH/UPARSE v9.0.2132 Illumina paired reads pipeline (Edgar, 2013) to filter out poor-quality sequences, de-replicate and sort reads by abundance, and remove singletons. We inferred operational taxonomic units (OTU) at the 97% sequence similarity level as usually used for OTU clustering (meaning that sequences differing by more than 3% are considered to belong to different OTUs; Stackebrandt and Goebel, 1994; Blaxter et al., 2005). We used the SINA v1.2.10 for ARB SVN (revision 21,008; Pruesse, Peplies & Glöckner, 2012) taxonomic reference dataset for both markers and used SILVAngs 1.3 for taxonomic assignments (Quast et al., 2012).

Construction of corrected OTU tables

Statistical analyses

OTU diversity and turnover in relation to soil properties (research 1 and 2)

In the physical soil data PCA, large values on the first PC were associated with coarse texture (coarse sand fraction loading 0.52, total sand fraction loading 0.55) and small values with fine texture (silt loading −0.45, clay loading −0.38; Table S2). The flooded forests (igapós and várzeas) generally had fine-textured soils (negative values of PC1). The unflooded forests (terra-firme and campinas) were more widely distributed along PC1, with some plots having similar values with várzeas and igapós (Fig. 2A). In the chemical soil data PCA, large values on the first PC were associated with poor soils. The most negative loading was −0.35 for the sum of exchangeable bases (SB), and the largest positive loading was 0.29 for aluminium saturation index (Table S3). The habitat types were not well separated along the first axis of the chemical PCA, as most plots of all habitat types had poor soils (large values of PC1) and just a few scattered várzea and igapó plots had more cation-rich soils (small values of PC1; Fig. 2B).

In general, organic carbon and pH had the strongest effects on OTU diversity. This was the case both for prokaryotes and eukaryotes, for richness and effective number of OTUs and for soil and litter (Table 1). In addition, PC1 of chemical soil properties was an important predictor for prokaryotic OTU richness in the soil and litter, with stronger effect in soil than in litter, and for prokaryotic effective number of OTUs just in soil. For eukaryotes, soil texture had an important effect on OTU diversity, albeit different on each substrate: positive for soil and negative for litter (Table 1). Overall, soil properties had strong effects on OTU diversity in litter.

OTU community turnover was significantly associated with soil properties, especially with organic carbon and pH, which were significant for all communities. The pH effect was strong for all prokaryote (16S) datasets and for eukaryotes (18S) in soil when relative abundance data were used (Table 2). Organic carbon had the strongest effect for eukaryotes in soil when presence/absence data were used and for eukaryotes in litter with both presence/absence and relative abundance data. Chemical PC1 was significant for prokaryotes in soil (both presence/absence and relative abundance) and for eukaryotes in litter when presence/absence data were used. Texture PC1 was significant only for eukaryotes in litter (both presence/absence and relative abundance; Table 2). Geographical distance was a significant explanatory factor for all datasets, but as closer places usually are more environmentally similar, we cannot separate the effect of spatial correlation from soil property effects.

A moderate percentage of the variation in Jaccard dissimilarities was explained by soil physicochemical properties in the presence/absence data, for prokaryotes (31% in soil and 35% in litter; Fig. 3). For eukaryotes, the total explanatory power of soil physicochemical properties was smaller (12% in soil and 16% in litter). For prokaryotes and eukaryotes, the litter communities were more structured by soil characteristics than were the soil communities. All variables explained small but significant proportions of the variance in all communities and showed some weak but significant interactions considering presence/absence matrices (Fig. 3) and a similar, strong proportion of the variance in abundance data (Fig. S2). Organic carbon had the strongest effect in all substrates and for both organism groups (ranging from 0.03 for eukaryotes in soil through 0.05 for eukaryotes in litter to 0.08 for prokaryotes in both soil and litter).

Similarities in OTU diversity and turnover patterns between litter and soil (research 3–4)

We found a weak positive regression between OTU richness of prokaryotes in litter and in soil (adj. R² = 0.25, p < 0.001; Fig. 4A) and between the effective number of prokaryote OTUs in litter and in soil (adj. R² = 0.1, p = 0.03; Fig. 4B). For eukaryotes, the corresponding correlations were not significant (Figs. 4C and 4D). The plot “CXNCAMP3” had very low soil OTU richness, and excluding this data point strengthened the correlation of OTU richness between soil and litter for prokaryotes (to adj. R² = 0.46, p < 0.001; Fig. S3A), but not for eukaryotes (Fig. S3B).

The OTU communities in litter and in soil tended to be separated in the NMDS ordination space, although there was some overlap especially for the igapó plots (Fig. 5). The PERMANOVA test indicated weak but significant effects (all p < 0.001) of substrate type on compositional dissimilarities of both prokaryotes (R² = 0.06, F = 5.83, for presence/absence data and R² = 0.07, F = 6.7, for abundance data) and eukaryotes (R² = 0.03, F = 2.8 for presence/absence data and R² = 0.04, F = 3.74, for abundance data). Habitat type had an even stronger effect on the compositional dissimilarities of both prokaryotes (R² = 0.17, F = 5.48, for presence/absence data and R² = 0.18, F = 5.68, for abundance data) and eukaryotes (R² = 0.1, F = 2.8, for presence/absence data and R² = 0.1, F = 2.98, for abundance data). Taxonomic composition at the phylum and kingdom (for fungi) level was similar in litter and in soil both for prokaryotes and for eukaryotes (Fig. 6). However, in prokaryotes Actinobacteria is the second most abundant phylum in litter when taking into account the relative abundances (Fig. 6).

Soil predictors of OTU diversity and community turnover

In this study, we tested the impact of physicochemical soil properties on the OTU diversity (richness and effective number of OTUs) and community turnover of prokaryotes and eukaryotes in soil and litter across Brazilian Amazonia. We found that that the soil properties we quantified had variable effects on OTU diversity and community turnover for litter and soil, and the effect varied between prokaryotic and eukaryotic organisms. The variable with the highest explanatory power was overall organic carbon for both prokaryotes and eukaryotes. OTU diversity and community turnover were better explained by soil properties in litter than in soil.

Considering the results from the linear models, in general organic carbon and pH were the strongest factors in explaining soil prokaryotic and litter and soil eukaryotic diversity. Our results show a positive correlation between soil pH and OTU diversity, which is expected since much of the soils in Amazonia are acidic. For instance, for soil samples, Lauber et al. (2009) found pH to be the main factor in explaining bacterial phylogenetic diversity and phylogenetic composition, where soils with pH between 4.5 and 8 had the highest bacterial diversity. Tropical forests with high macro-organismic diversity had soil with pH < 4.5 and had the lowest bacterial diversity (Lauber et al., 2009). In our samples, pH was overall low and its variation was moderate, from 3.65 to 5.14, thereby in less acid soil we found highest OTU diversity considering both richness and effective number of OTUs.

We also found that variation in pH was significant for all community turnovers, and that it was the strongest variable in explaining community turnover for prokaryotes (in both soil and litter) and for eukaryotes in soil (Table 2). However, we found that pH had no effect on prokaryote richness in litter (Table 1). The consistent effect of pH for prokaryotes and eukaryotes in the soil, with significant effect in community turnover and diversity, but inconsistent effect in litter support the findings of Gregorich et al. (2017), who found no correlation of soil properties and litter decomposition. This points to the independence of the environmental factors regulating each substrate.

We found significant effects of variation in organic carbon on all community turnovers with the strongest effect for eukaryotes in litter (Table 2). Additionally, we found a negative correlation of soil organic carbon with OTU diversity for all groups (Table 1). Soil biodiversity has previously been found to have an effect on carbon sequestration (Wagg et al., 2014). However, the relationship between soil biodiversity and carbon has varied across studies (Nielsen et al., 2011). Furthermore, Fierer et al. (2012) and De Lima Brossi et al. (2014) found that soil organic matter was related to microbial community composition in several different vegetation types. The negative correlation between soil organic carbon content and OTU diversity reported here might be related to high nutrient turnover in high-diversity soil/litter, keeping the carbon stock locked in aboveground biomass. Our results support the findings of Wall et al. (2008), who found a positive influence of the richness of soil biota on decomposition rates in wet tropical environments. Along the same line, Wagg et al. (2014) found that soil diversity and soil community composition are related through nutrient cycling. Decreases in soil diversity and the related changes in community composition alter the communities’ capacity to break down organic matter and recycle nutrients, slowing down the return of nutrients to the above-ground communities (Wardle et al., 2004). These findings stress the complex nature of carbon-diversity dynamics and the plant–soil feedback loop mediated by soil biota (Mangan et al., 2010). They furthermore highlight a connection between decomposition rates and biodiversity in Amazonia that should be better explored.

We found that prokaryotic community turnover is more strongly related with environmental distance than eukaryotic community turnover which is mostly dominated by fungi (fungi correspond to 35% of OTU richness and 50% of OTU relative abundance in our eukaryotic data; Figs. 6C and 6D). This is in agreement with the results from a global bacterial and fungi soil sampling (Bahram et al., 2018). In our results, environmental distances explained only a limited percentage of OTU community turnover (31–35% for prokaryotes and 12–16% for eukaryotes; Fig. 3). This suggests that other factors, such as precipitation and bacterial-fungal antagonistic interactions (Bahram et al., 2018), also need to be considered to better understand the community turnover these organisms.

Biotic and abiotic interactions jointly determine soil properties, making it important to consider environmental and biological interactions between variables. Indeed, our variance analysis reveals several co-variances between soil properties, such as pH and organic carbon and physical and chemical properties. Although these interactions were weak, this analysis is important for providing a better understanding of the study system. Considering physicochemical soil properties, we had a partial separation of the major environmental types by the properties of their soils. It is in agreement with previous studies, which report an association of soil types and habitats in Amazonia (e.g., Falesi, 1984; Prance, 1996). The soil texture (first axis of the physical PCA) was well separated by the habitat types of flooded forests (igapós and várzeas), whereas terra-firme and campinas were more spread in physical properties. On the other hand, the first axis of chemical PCA was less well separated for flooded areas (igapós and várzeas). This result was expected since there is a variation within habitat types, especially flooded forests (Kalliola et al., 1993; Tuomisto, Ruokolainen & Yli-Halla, 2003; Tuomisto, Zuquim & Cárdenas, 2014; Tuomisto et al., 2016). Furthermore, terra-firme forests have been reported to vary in soil nutrients (e.g., Tuomisto, Ruokolainen & Yli-Halla, 2003; Tuomisto, Zuquim & Cárdenas, 2014; Tuomisto et al., 2016; Fine et al., 2005), consistent with the variation observed in our plots. However, due to the limited sampling in our studies, the variation we detected was small (Fig. 2B). Our finding that soil texture is similar among the flooded environments (várzeas and igapós) and that soil texture was an important factor for eukaryote diversity (both soil and litter) is consistent with the previously reported community similarity among these environmental types based on the data from the same samples (Ritter et al., 2018).

Contrasting litter and soil diversity

The correlations between soil and litter OTU diversity (richness and the effective number of OTUs) were significant for prokaryotes but not for eukaryotes. This is congruent with previous reports that showed the independence of litter accumulation from properties of the underlying soils (Gregorich et al., 2017).

We expected a difference in taxonomic composition between litter and soil communities, with microbes dominating the soil (e.g., Bates et al., 2013; Mahé et al., 2017) and plant and nematode OTUs dominating the litter due to it mainly being composed of leaves and roots. However, we found the highest plant (Chloroplastida) richness and abundance in the soil samples. Furthermore, unlike Porazinska et al. (2012) who found a dominance of nematodes in the litter of tropical forests, we found very similar proportions of nematode OTUs in soil and litter, with the highest richness and abundance in the soil (Figs. 6C and 6D). This pattern of no clear taxonomic differences between soil and litter layers is consistent with respect to all dominant groups at the phylum and kingdom levels (Fig. 6). This suggests that, on the Amazon basin scale, the taxonomic composition at higher levels (phylum and kingdom) is consistent between litter and soil with some variation in abundance for some groups, such as Arthropoda and Chloroplastida for eukaryotes and Actinobacteria, Bacterioides and Chloroflexi for prokaryotes (Fig. 6).

Interestingly, the prokaryotic phyla that dominated our samples were only partly the same as those found dominant in a large global dataset (Figs. 6A and 6B; Delgado-Baquerizo et al., 2018). While we also found Proteobacteria to be the most frequent phylum considering both presence/absence and relative abundance, the second most frequent phylum was Chloroflexi in presence/absence data for soil and litter and abundance for soil in our samples, while this Chloroflexi was only the 5th most abundant in the global dataset. Actinobacteria, the second most abundant phylum in the global database, was in our data the second most abundant phylum just for abundance in litter samples. Moreover, the rank-abundance distribution of the most dominant phyla was more even in our tropical sample than in the global sample, with Proteobacteria accounting for just over 20% of all reads (versus almost 40% in the global dataset) and eight phyla representing more than 5% of relative frequency each (>70% of relative frequency) versus only four phyla in Delgado-Baquerizo et al. (2018). Taken together, these differences highlight the need for more studies across the Amazon basin to better characterize the taxonomic composition.

The OTU community compositions of both prokaryotes and eukaryotes were better explained by habitat type (terra firme, várzea, igapó, campina) than they were by substrate type (soil, litter), which was expected since both substrates should share a large number of organisms. The substrate types were weakly differentiated at the OTU level, but we could not observe any difference at the phylum or kingdom levels for presence-absence and only a small difference for abundance data (Fig. 6). For instance, fungi usually dominate eukaryotic soil communities in any environment, including tropical forests (Tedersoo et al., 2017), but the dominant fungal taxa (OTU) may vary considerably even at local and sub-local scales (Urbanová, Šnajdr & Baldrian, 2015). In a study conducted in the western parts of the Czech Republic, similar results for bacteria and fungi were found: the phylum level indicated the same taxonomic groups as dominant in soils and litter, but there were striking differences at the OTU level in these substrates (Urbanová, Šnajdr & Baldrian, 2015).