Consequences of agroindustrial sugarcane production to freshwater biodiversity

The environmental benefits of a broad‐scale adoption of biofuels are critically contingent on what current land uses will be converted for feedstock expansion and how converted land will be managed. We assessed the consequences of land use and land management for the agroindustrial production of sugarcane to the physical, chemical, and biological properties of freshwater systems. We surveyed 16 environmental variables and algae, tadpoles, predatory invertebrates, and fish in lentic water bodies distributed across a gradient in land‐use intensity ranging from seasonal Atlantic Forest and cerrado to pastures to sugarcane plantations in SE Brazil, the most important sugarcane‐producing region in the world. The gradient in land‐use intensity was not only an axis of native habitat loss but also of ecosystem productivity, as indicated by increased conductivity, turbidity, and phytoplankton biomass. Land use had a clear signal on community and metacommunity organization, with converted land being impoverished in biodiversity relative to native habitats. However, frequency of occurrence, density, biomass, and alpha diversity of tadpoles and their predators were not affected by land use. These results suggest that sugarcane fields function as habitat to a fraction of aquatic biodiversity. Within sugarcane fields, larger wetlands surrounded by buffer strips as required by law appeared comparatively buffered against land management practices and housed a disproportional fraction of animal biomass, likely acting as sources of migrants to other water bodies in the landscape. Conversion of pastures to sugarcane fields, suggested as a strategy to reduce competition for land with food production and biodiversity conservation, does not appear to have strong consequences to lentic freshwater systems, provided that wetlands and surrounding buffer strips are preserved. These observations emphasize the importance of enforcement of legislation regulating land use (i.e. the ‘Forest Code’) and certification systems verifying compliance and rewarding the voluntary adoption of better land management practices.


Introduction
Global energy demand is predicted to rise by 37% from 2012 to 2040 and, despite the continuing predominance of fossil fuels, consumption of biofuels for transportation is to be increased more than threefold (IEA, 2014). Biofuels are expected to alleviate human reliance on nonrenewable fossil fuels while contributing to a reduction in greenhouse gas emissions and to the promotion of rural development. However, the environmental benefits of biofuels are critically contingent on which, where, and how biofuel feedstocks are produced. This occurs because feedstocks vary broadly in energetic efficiency, because their expansion may promote the conversion of lands of high conservation value and release large quantities of stored carbon, and because intensive agriculture is linked to various forms of environmental degradation (OECD, 2001(OECD, , 2006Clay, 2004;Fargione et al., 2008;Lapola et al., 2010).
Sugarcane is the most energetically efficient first-generation source of ethanol, yielding over eight units of biofuel energy output per-unit energy input when compared with two for beet, wheat, or corn (WWI, 2006). Echoing a global expansion of biofuel crops, between 2002 and 2012 sugarcane land cover in Brazil almost doubled to reach a total of 9.7 million hectares, and it has to expand another 1.4-4.3 million hectares if the country is to reach its 2020 ethanol production goals (UNICA, 2008;Lapola et al., 2010). Brazil currently responds to 35% of all sugarcane produced in the world, half of which concentrated in southeastern state of São Paulo (UNICA, 2014).
Despite significant recent advances, sugarcane production in Brazil has been historically associated with serious environmental impacts. These include deforestation of the Atlantic Forest and to a lesser extent the cerrado, large-scale preharvest burns, extensive soil erosion, high water consumption in mills, and deliberate release of vinasse in streams and rivers causing eutrophication, algal blooms, anoxia, and fish kills (Dean, 1997;Ballester et al., 1999;Gunkel et al., 2007;Martinelli & Filoso, 2008). In addition, sugarcane production has a significant potential for environmental contamination from the moderately high consumption of inorganic and organic fertilizers (Cantarella & Rossetto, 2010;UNICA, 2011), and the moderate consumption of a broad range of pesticides (225 formulations containing 62 active ingredients; Schiesari & Grillitsch, 2011).
Perhaps surprisingly then, explicit assessments of the impacts of sugarcane production on biodiversity are rare. For example, a recent review of the impacts of global biofuel crop production on biodiversity failed to find a single article on sugarcane (Immerzeel et al., 2014; but see, e.g. Corbi et al., 2008;Verdade et al., 2012).
Freshwater systems constitute an appropriate system to assess the effects of land use and land management on biodiversity because are intimately connected to surrounding terrestrial systems by the transfer of energy and matter. Furthermore, water bodies usually sustain diverse communities comprising species with varying degrees of dependence on the aquatic environment. Whereas many species are obligate aquatic, many others have complex life cycles and as such integrate the cumulative effects of freshwater and terrestrial degradation over the course of their life history. Finally, there is a solid tradition in freshwater theoretical and applied ecology that can provide a basis for interpreting landuse effects on traits, individuals, populations, and communities (e.g. Wetzel, 2001).
This study aimed at assessing the consequences of land use and land management for the production of sugarcane on the physical, chemical, and biological properties of freshwater systems, with an emphasis on the diversity, composition, and structure of their communities.
This assessment was based on sampling surveys in replicated lentic (i.e. still) water bodies distributed across a gradient in land-use intensity ranging from seasonal Atlantic Forest, wooded cerrado ('cerradão') and cerrado to pastures to sugarcane plantations in the state of São Paulo, and employing amphibian larvae and their predators as the study system. The sampling design thus captures three types of land conversion associated with the historical expansion of sugarcane: the conversion of native habitats to pasture, the conversion of native habitats to sugarcane, and the conversion of pasture to sugarcane (Dean, 1997). The latter, in particular, is central in the debate regarding the current and future expansion of biofuel crops in the country. A recently approved agroecological zoning for the sustainable production of sugarcane identified as suitable 37 Mha currently covered by low-profitability pastures and concluded that sugarcane expansion could occur over pastures with no impacts on land covered by native vegetation or devoted to food production (Manzatto et al., 2009;Federal decree 6961/2009). Other studies arrived at a similar recommendation while targeting at a reduction in carbon emissions (Lapola et al., 2010). Therefore, considering that industrial scale sugarcane production follows an intensive model of agriculture, this raises the question not only of the environmental impacts of sugarcane production relative to native habitats, but also relative to pastures, which may support moderate levels of biodiversity.

Study site
Field work was conducted in the Estac ßão Ecol ogica do Jata ı (EEJ) and surrounding farms in the municipality of Luis Antonio, São Paulo, southeastern Brazil (21°34 0 11 00 S, 47°44 0 07 00 W). The EEJ represents one of the very few preserved (mostly second growth) remnants of cerrado, wooded cerrado ('cerradão'), and seasonal Atlantic Forest vegetation in the state of São Paulo, and is embedded in a matrix largely dominated by intensive sugarcane cultivation, with a smaller occurrence of pastures and silviculture. Around the EEJ, we were granted access to two farms cultivating sugarcane in industrial scale and one small ranch raising cattle and sheep. According to local agronomists, land management for sugarcane production in these farms are typical for the broader region of Ribeirão Preto, which is the most traditional sugarcane-producing region in the country.

Selection of water bodies
We selected 15 lentic water bodies in native habitats, 15 in pastures, and 17 in sugarcane fields based on the following criteria: (i) that the water body was amenable to pipe sampling, that is, that a considerable fraction of the water body was less than 60 cm deep; (ii) that the water body in a given land-use type was at least 100 m away from another land-use type; (iii) that the water body was independent from and at least 50 m away from other sampled water bodies; (iv) that each land-use type included temporary, semipermanent, and permanent water bodies, given the importance of hydroperiod in influencing freshwater community composition (Wellborn et al., 1996). Water bodies comprised puddles, ditches, ponds, farm ponds ('cacimbas'), reservoirs, swamps, and marshes varying in size from 1.4 to 106 000 m 2 (see Fig. 1). Not all of the selected water bodies were available for each analysis because some were dry or did not contain tadpoles at the moment of sampling. monthly visits in which water depth >0 cm). In December 2010 and February 2011, respectively, at the beginning and middle of the rainy season, water bodies were subject to comprehensive habitat and physicochemical characterization, amphibian calling surveys, and freshwater community surveys. Timing of sampling surveys in this seasonal landscape was chosen based on the phenology of occurrence of water bodies (two-thirds of which only existed in the rainiest months) and their communities.

Physical and chemical characterization of water bodies
Prior to the community sampling surveys, we conducted a comprehensive characterization of water bodies. Habitat characterization comprised recording spatial positioning, land-use type, any ongoing activities of land management, water body dimensions, canopy cover over the pond basin (in %, with a spherical crown densiometer; Forestry Suppliers, Jackson, MS, USA), substrate type (rocks, pebbles, sand, clay, mud, wood, leaf-litter) and cover (each substrate type categorized as comprising 0-24%, 25-49%, 50-74%, 75-100% of the pond basin), aquatic vegetation type (i.e. terrestrial/paludicolous plants, floating macrophytes, submerged macrophytes) and cover (as above), and a description of the surrounding terrestrial vegetation. We then quantified in situ basic physicochemical water parameters including temperature, pH, conductivity (with an Orion 3-star pH meter and an Orion 3-star conductivity meter; Thermo Scientific, Waltham, MA, USA), and dissolved oxygen (henceforth DO; with a YSI Pro20 DO meter; Fig. 1 Effects of land use on selected environmental variables related to water body structure (left), water physicochemistry (middle) and productivity (right). See Table S1 for univariate statistical analyses. Data from both sampling surveys are lumped for each response variable. Because in each land-use water bodies were selected so as to cover the full spectrum of the hydroperiod gradient, the inclusion of area and hydroperiod in this figure is meant to demonstrate that this important gradient was controlled in our sampling surveys. Corrected a = 0.002. Yellow Springs, OH, USA), as well as turbidity (with an AQ 4500 Advanced Aquafast IV Turbidimeter; Thermo Scientific, Waltham, MA, USA) and phytoplankton standing crop (in triplicate, with an Aquafluor Handheld Fluorometer with in vivo chlorophyll-a and phycocyanin channels; Turner Designs, Sunnyvale, CA, USA). We also collected and preserved water samples for laboratory analyses of water hardness, total nitrogen (TN), and total phosphorus (TP). These parameters were analyzed with a HACH DR 2800 spectrophotometer and appropriate reagent kits (HACH, Loveland, CO, USA).

Calling activity of amphibians
Prior to the freshwater community sampling surveys, we recorded in each water body the calling activity of amphibian adults to provide a more comprehensive assessment of amphibian community composition, to assess the efficiency in larval sampling, and to aid in larval identification. Calling surveys were semiquantitative (i.e. in the following categories of abundance: 0, 1-2, 3-5, 6-10, 11-20, 21-50, or more than 50 calling individuals; Heyer et al., 1994) and conducted for a total of 2-4 nights in each water body.
Samples were taken by pipe sampling, dipnetting, and seining (e.g. Heyer et al., 1994;Werner et al., 2007). Pipe sampling provided quantitative per-unit-area information on species composition, abundances, and biomass. Pipes were 60 cm tall and 33 cm in diameter and therefore sampled 0.09 m 2 of water body habitat. The sample was taken by quickly thrusting the pipe through the water column and into the sediments to seal the sample area. We then removed by dipnetting all animals >5 mm from the sampled water column and the first centimeter of the sediments. In practice, we considered the pipe to be empty when at least 10 consecutive sweeps were made without capturing any more animals. Sampling effort was scaled to pond surface area: Five pipes were sampled in water bodies <10 m 2 , 15 pipes in water bodies >10 m 2 and <50 m 2 , 30 pipes in water bodies >50 and <200 m 2 , 60 pipes in water bodies >200 m 2 and <500 m 2 and 90 pipes in water bodies >500 m 2 . Pipes were distributed so as to representatively sample all microhabitats in a water body. In water bodies with deeper waters (i.e., deeper than sampling pipes; NAT01, NAT05, PAS09, SUG08, SUG09), we supplemented pipe sampling with three hauls of a 10-m seine in the deeper water. Tadpoles and fishes were immediately euthanized and preserved in 10% buffered formalin, invertebrates in 70% ethanol.
In the laboratory, all biological samples were sorted and identified to family (invertebrates) or species level (vertebrates), except for the amphibians Scinax and Rhinella, which were assigned to genus only because could not be confidently assigned to species (at our study site, there are two species in the genus Scinax, S. fuscovarius and S. similis, and two in the genus Rhinella, R. ornata and R. schneideri), counted and weighed. Fishes with detritivorous or herbivorous feeding habits were excluded from the analysis. Summarizing, in each water body we obtained as biological response variables community composition, species or family richness, densities, and per-unit area biomasses of tadpoles and predators.

Data analyses
Effects of land use on environmental and biological variables. We first conducted two-way analyses of variance followed by Tukey HSD post hoc tests for testing the effects of land use, survey (i.e. whether in the beginning or middle of the rainy season), and a land-use-by-survey interaction term on each individual variable. The effects of survey could not be assessed for hydroperiod as this is a year-round summary metric, and for conductivity and TN, as these variables were quantified in the first survey only due to technical issues. In all other cases, a land-use-by-survey interaction term was found nonsignificant and was dropped from the analysis. When the assumptions of normality and homogeneity of variances could not be met even after transformations, we employed Kruskal-Wallis nonparametric tests followed by pairwise comparisons. In all cases, we applied Bonferroni corrections to avoid Type I error inflation.
We then conducted two separate permutational multivariate analyses of variance (perMANOVA) to test for multivariate landuse differences in the environmental and community datasets. PerMANOVAs were run on the Euclidean distances of the logtransformed environmental data (Anderson et al., 2006) and on the Bray-Curtis distances of the Hellinger-transformed tadpole species density data (Legendre & Gallagher, 2001). Both transformations were implemented in the function decostand of the R package vegan (Oksanen et al., 2012). The homogeneity of variances among groups was checked previously to the perMA-NOVA using the function betadisper from the R package vegan (Oksanen et al., 2012). When significant differences among land-use types were detected, we ran multiple comparisons with Bonferroni-adjusted significance levels. R version 2.15.2 (R Core Team, 2012) was used in all statistical analyses.
Alpha, beta, and gamma diversity. We compared alpha diversity among land-use types by an analysis of variance on the total number of larval amphibian species, or predator families, observed in each water body. We also compared amphibian beta diversity among water bodies within each land-use type using the procedure proposed by Baselga (2010Baselga ( , 2012 and implemented in the R package betapart (Baselga & Orme, 2012). This procedure decomposes the overall beta diversity (b SOR ) into two components: (i) spatial turnover (b SIM ), which is the fraction due to replacement of some species by others along the sites, and (ii) nestedness (b NES ), which is the fraction due to difference in species richness, such as sites with lower richness being a subset of the richest site (Baselga, 2010(Baselga, , 2012. To test the difference between beta diversity components, we used a resampling procedure where we took 1000 random samples of five sites per area and calculated the average b-values. With the estimated distribution of b values, we assessed the between-area significance by calculating the probability of the area with smaller b-value being larger than the other just by chance (Baselga, 2010). We decided not to compare predator beta diversity within each land-use type because of the difficulty in interpreting this metric at such a crude taxonomic (i.e. family) level. We compared gamma diversity among land-use types by rarefaction curves extrapolating species richness in each land-use type to 11 samples, that is, the maximum number of water bodies sampled in a given land-use type (Colwell et al., 2012). Rarefaction was performed using the software ESTI-MATES (Colwell, 2013).

Predictors of amphibian species richness.
To assess what factors influence local amphibian species richness (i.e. alpha diversity), we performed a generalized linear model (GLM) with a model selection and model averaging approach considering land use and seven variables recognized in the literature as influencing tadpole richness. These were water body area, hydroperiod, canopy cover, the presence of fish, predator density, competitor density, and distance to forest edge (Wellborn et al., 1996;Werner et al., 2007). Quantitative variables ranged over several orders of magnitude and were log-transformed prior to the analysis (Anderson et al., 2006) using the function decostand of the R package vegan (Oksanen et al., 2012). We assessed the occurrence of multicollinearity among quantitative variables by the variance inflation factor (VIF). All VIFs were below the suggested threshold (VIF < 3; Zuur et al., 2009Zuur et al., , 2010. We fitted a global model to the data using a Poisson distribution (no overdispersion detected). The global model contained all eight explanatory variables and from that model we generated a set with all possible submodels using the function dredge from the R package MuMIn (Barton, 2012). We decided to run all possible models as we a priori selected only variables known to influence amphibian richness. We employed Akaike's information criterion corrected for small sample sizes (AICc) to rank the obtained models and selected models with ΔAICc < 2 (Burnham & Anderson, 2002) to compute the model-averaged parameters using the function model.avg in MuMIn package (Barton, 2012). We did the model averaging with shrinkage, where in models where a variable is absent its coefficient is set to zero (Burnham & Anderson, 2002). We also assessed the likelihood of the models using Akaike weights, which indicates how likely the model is the best among the ones considered.
Community structure, land use, and associated environmental variables. To determine the exclusive and shared effects of environmental variables and land use on community structure we conducted redundancy analyses (RDA) using the variation partitioning method with adjusted R 2 following Borcard et al. (2011). To objectively reduce dimensionality while maintaining interpretability, we first separated the 18 environmental variables in structural properties (area, depth, hydroperiod, canopy, vegetation cover, organic sediment cover), basic physicochemical water properties (temperature, pH, DO, conductivity, turbidity, hardness), properties related to ecosystem productivity (total nitrogen, total phosphorus, chlorophyll-a concentration, phycocyanin concentration), and predation pressure (density of fish, density of predatory invertebrates). We ran one PCA on each of these groups and retained the site scores of the first and second axes as explanatory environmental variables.
Tadpole density and environmental variables were, respectively, subject to the Hellinger (Legendre & Gallagher, 2001) and log-transformation (Anderson et al., 2006), both implemented in the function decostand of the R package vegan (Oksanen et al., 2012). We selected the environmental variables to be used in the variation partitioning by a forward selection method (Blanchet et al., 2008;Borcard et al., 2011). We tested by a permutation test with 999 permutations (Borcard et al., 2011) the significance of the testable fractions of the explanatory variables in the variation partitioning. All analyses were performed in the R package vegan (Oksanen et al., 2012).
Analyses using the first and the second sampling surveys, and a summary dataset comprising both surveys had similar results. Therefore, we present here only the analyses with the summary dataset. The summary environmental explanatory matrix was constructed using the average values of the two surveys for variables and water bodies sampled in both surveys; the only value available for variables or water bodies sampled in a single survey; and the maximum density of invertebrate and fish predators across surveys. The summary community matrix was constructed using the maximum density of each tadpole species per-pond across surveys.
Metacommunity structure, land use, and associated environmental variables. To analyze the structure of the amphibian metacommunity, we first summarized the distribution pattern of 21 amphibian species in 26 water bodies (three additional water bodies were excluded for not containing amphibians) by an incidence matrix representing the site-by-species presence/ absence data. For completeness, we included both calling survey and tadpole survey data in this incidence matrix; analyses of an incidence matrix comprising only tadpole survey data yielded qualitatively similar results.
We tested whether incidences differed from a null model of random placement of species in water bodies by ordering rows and columns of the incidence matrix by reciprocal averaging (i.e. correspondence analysis), such that water bodies with the most similar species lists and species with the most similar distributions were closest together. We then tested this ordinated matrix for its coherence, that is, whether embedded absences within species distributions were significantly less common than expected by chance, using the null model randomization method r0 (Presley et al., 2010). Significant coherence implies that there is a strong axis of variation in metacommunity structure, and justifies further characterization of elements of metacommunity structure in terms of species turnover or nested subsets, and if there is a significant degree of boundary clumping among species distributions (Leibold & Mikkelson, 2002). Altogether, these analyses permit testing whether our incidence matrix fits one of the several idealized patterns of metacommunity structure (Clemmentsian vs Gleasonian gradients, nested subsets, evenly spaced, checkerboards), which could in some cases provide insights into hypothesized underlying structuring mechanisms.
We performed Spearman rank correlations between site ordination scores and site richness and site environmental variables to interpret the axis of environmental variation associated with the latent gradient (as in Meynard et al., 2013) and tested by means of an ANOVA on site ordination scores whether land use significantly influences community arrangement along the ordination axes.

Land-use patterns of environmental properties of water bodies
Land use had a significant univariate effect on seven of 18 measured environmental variables. These were canopy cover, temperature, conductivity, turbidity, water hardness, and chlorophyll-a and phycocyanin concentrations ( Fig. 1, Table S1). Among the three forms of land use, water bodies in native habitats had the lowest turbidity, water hardness, and phycocyanin concentrations. In turn, water bodies in sugarcane plantations had the highest temperatures and chlorophyll-a concentrations. The only distinctive characteristic of water bodies in pastures were the lowest temperatures. Irrespective of land use, a significant effect of survey was detected for chlorophyll-a concentration (greater in December), and for TP and phycocyanin concentration (greater in February).

Land-use patterns of amphibians and their predators
Most water bodies in all three forms of land use contained amphibians and their predators (Fig. 3a). Tadpoles were found in 70% of the water bodies sampled in native habitats, but in 100% of the water bodies sampled in pastures and sugarcane fields. In turn, predators were found in 100% of the water bodies sampled in native habitats and pastures, but in 70% of the water bodies sampled in sugarcane fields. Interestingly, land use had no univariate effect on the density and per-unit area biomass of tadpoles or of their predators (Fig. 3b,c, Table S1).
Although the three water bodies with the highest amphibian species richness were all found in native habitats (12, 10 and 7 species), and although the richest water bodies in pastures and sugarcane fields had half as many species (6 species), no significant univariate effect of land use on amphibian species richness was detected either (Fig. 3d, Table S1). This occurred because of the high among-pond variability in amphibian species richness in all three land-use types (0-12 species in native habitats, 2-6 in pastures, and 1-6 in sugarcane fields). Likewise, the overall amphibian beta diversity (b SOR ) and its two components, spatial turnover (b SIM ) and nestedness (b NES ), were similar among land use In contrast, gamma diversity differed among land-use types. Total measured amphibian richness was 17 species in native habitats, 8 in pastures, and 10 in sugarcane fields (Fig. 3e). Correspondingly, extrapolation of rarefaction curves of species richness to 11 samples, that is, the maximum number of water bodies sampled per land-use type, indicated that species richness in native habitats (mean 17.62, 95% CI 13.45-21.79) was significantly higher than species richness in pastures (8.68, 6.24-11.12) and sugarcane fields (10, 8.71-11.29), which, in turn, did not differ from each other.
While the analysis of patterns of predator diversity should be exercised with caution due to the crude taxonomic resolution, no univariate effects of land use on the per-pond number of predator families were detected ( Fig. 3d; Table S1). Overall, the number of predator families declined from preserved to converted landscapes (Fig. 3e), but this occurred exclusively due to fish, which were more common in native habitats because of three water bodies seasonally linked to streams.

Predictors of amphibian species richness
From the set of models generated to infer the best combination of environmental variables predicting amphibian species richness, five models were best supported (five models had DAICc < 2; Table 1, Table S2). Hydroperiod and tadpole density were present in all five top models, and therefore, were assigned a relative importance of 1.00 after model averaging (Table 1). Tadpole species richness increased with tadpole density and hydroperiod, that is, water bodies that contained more tadpoles per-unit area and that held water for longer  Fig. 3 (a) Frequency of occurrence (b) Density (c) Per-unitarea biomass (d) Richness per pond and (e) Richness per landuse type of tadpoles (left) and predators (right) in water bodies distributed across a gradient in land-use intensity comprising native habitats, pastures, and sugarcane plantations. Here frequency of occurrence refers to the % of water bodies that contained tadpoles or predators in either the first (December) or second (February) community surveys. Density and per-unitarea biomass boxplots lump data from both surveys, and richness per pond and richness per land-use type consider the cumulative species list found in a given pond (or land use) across surveys. Richness for amphibians is species richness, for invertebrates (which include dragonflies, water bugs, beetles, and fish) is family richness. Because few fish species were captured by pipe sampling, estimates of predator richness include fish captured by seine nets as well. Land use had a significant effect in (e) but not in (b), (c), or (d). See Table S1 for univariate statistical analyses. periods had higher tadpole richness (Fig. S1). Because the 95% confidence intervals of all other variables included zero (Table 1), there is not enough evidence that distance to native habitats, fish presence, and canopy cover affect tadpole species richness in our study.
Land use, however, comprises several axes of environmental change that could, alone or in combination, influence amphibian community structure. Among the many descriptors of freshwater habitats, the forward selection procedure retained four composite environmental variables in our approach to reduce dimensionality. These were the first PCA axis of properties related to ecosystem productivity (prod1), the first PCA axis of water physicochemical properties (phys_chem1), and the first and second PCA axes of predation pressure (pred1, pred2) (see Table S3 and Fig. S2 for details). Axis prod1 explained 65.4% of the total among-site variation and was positively related to TN. Chlorophyll-a and phycocyanin concentrations, and TP, were weaker positive contributors to prod1. Axis phys_chem1 explained 69.5% of the total among-site variation and was positively related to turbidity. Conductivity and water hardness were weaker positive contributors to phys_chem1. Axes pred1 and pred2 explained all total among-site variation (57.3% and 42.7%, respectively), with pred1 being negatively related to both the maximum density of predatory invertebrates and fish, and pred2 being positively related to the maximum density of fish but negatively related to the maximum density of predatory invertebrates.
The above-mentioned selected composite environmental variables (axes of PCA) and land use together explained 48.6% of the total tadpole community variability (R 2 adj = 0.49, P < 0.01). From the total variability, 31% was explained exclusively by the selected environmental variables (R 2 adj = 0.31, P < 0.01), and over fourteen percent was explained by the shared component between environmental variables and land use (R 2 adj = 0.14, p not testable). Land-use alone explained less than three percent of the tadpoles community variability and was not significant (R 2 adj = 0.03, P = 0.11). The first RDA axis was negatively related to phys_chem1 and pred1, and less so to prod1 and pred2 (see Table S4 and Fig. S3 for RDA site and species scores, and biplots of individual selected composite environmental variables). Therefore, increasing score values along the first RDA axis, among other characteristics, described a gradient ranging from water bodies with high turbidity and TN, and less predatory invertebrates and fishes (water bodies mainly located in sugarcane fields), to water bodies with low turbidity and TN, and more predatory invertebrates and fishes (water bodies mainly located in native habitats and pasture). The second RDA axis was positively related to pred1 and negatively related to pred2 (Fig. S3). Therefore, increasing score values along the second axis described mainly a gradient from water bodies with more predatory fish and less predatory invertebrates to water bodies with less predatory fish and more predatory invertebrates. Most species were found in water bodies devoid of fish; the glaring exception was the bufonid Rhinella, typically found in water bodies with high densities of fish. Scinax and Physalaemus cuvieri were usually related to less productive water bodies, with lower TN and lower turbidity. By contrast, Eupemphix nattereri and Leptodactylus fuscus were usually associated to more productive water bodies with higher TN and higher turbidity.
Amphibian metacommunity structure, land use, and associated environmental properties The amphibian metacommunity exhibited strong coherence along the first axis of ordination (z = 3.21, P < 0.001; Fig. 4). There were 141 embedded absences in a 21 9 26 species-by-site incidence matrix (546 cells). The community exhibited a nonsignificantly positive turnover (z = À1.39, P = 0.164) and was significantly clumped (Morisita's index = 2.559, df = 18, P < 0.001). The first axis of ordination was positively correlated with land-use intensity, conductivity, turbidity, TN, TP, chlorophyll-a, and phycocyanin concentrations (Spearman rank correlation, P ≤ 0.027) and negatively correlated with area, hydroperiod, vegetation cover, and organic sediment cover (P ≤ 0.017). As hypothesized, land-use intensity had a significant effect on ordination scores (ANOVA F 2,26 = 3.51 P = 0.047): post hoc analyses indicated that ordination scores of native habitat communities were significantly different from those of sugarcane communities, whereas ordination scores of pasture communities were not significantly different from either those of native habitats or of sugarcane plantations. Species richness was significantly correlated with site scores of the first axis (Spearman rho = À0.47, P = 0.016).

Discussion
Land use and land management for the production of sugarcane had a clear signal on the physical, chemical and biological properties of freshwater systems. In a  Fig. 4 Amphibian species-by-site incidence matrix for the Estac ßão Ecol ogica do Jata ı and surrounding agricultural areas. The matrix is ordered by reciprocal averaging such that water bodies with the most similar species lists and species with the most similar distributions are closest together. Position of each water body across three important freshwater environmental gradients is represented in color coded rows at increasing shades of gray at the top of each matrix (intensity of land use: native habitats, pastures and sugarcane plantations; hydroperiod: 0-49%, 50-99%, and 100% days with water; canopy cover: 0-33%, 34-66%, and 67-100% canopy cover). Arrows display environmental factors significantly associated with site ordination scores according to a Spearman rank correlation. multivariate space defined by sixteen environmental variables, there was virtually no overlap among water bodies in preserved land and water bodies in sugarcane fields. Water bodies in pastures were significantly different from those in native habitats, but not from those in sugarcane fields. Likewise, land use had a clear signal on community and metacommunity organization. This signal was revealed both by guided (i.e. nMDS and associated perMANOVA) and by unguided ordination tools (i.e. RDA, reciprocal averaging in the analysis of elements of metacommunity structure) and emerged even after accounting for other gradients of known importance in structuring freshwater metacommunities such as hydroperiod, canopy cover, and water body surface area (Wellborn et al., 1996;Skelly et al., 1999;Werner et al., 2007). Amphibian communities in native habitats were significantly different from those in sugarcane fields; pastures were intermediate in community structure in that they were not significantly different neither from native habitats nor from sugarcane fields. Therefore, both abiotic and biotic data lend general support to our hypothesis of a gradient in land-use intensityand therefore environmental degradationranging from preserved cerrados and seasonal Atlantic Forests to pastures to sugarcane fields.

Land-use patterns of environmental properties of water bodies
The transition of native habitats to productive agricultural land was associated with a decrease in canopy cover and to an increase in turbidity, water hardness, and abundance of cyanobacteria in lentic water bodies. In turn, the transition of pastures to sugarcane plantations was associated with a further reduction in canopy cover, to an increase in phytoplankton standing crop, and to a sharp increase in water conductivity. Overall, the physical, chemical, and biological characterization of water bodies indicated that, in addition to an axis of habitat loss and fragmentation, the gradient in land-use intensity is clearly a gradient in ecosystem productivity. In fact, every multivariate analysis conducted in this study strongly pointed to a gradient in land-use intensity on indicators of productivity, especially phytoplankton standing crop. This enhanced productivity is likely a product of increased incident solar radiation due to canopy opening (Schiesari, 2006), and of the increase in nutrients from cattle grazing and defecation in pastures, and from fertilization in sugarcane fields. Indeed, sugarcane is an intense consumer of fertilizers, and overuse is the rule (Martinelli & Filoso, 2008). The likely fate of a considerable fraction of applied nutrients is water bodies, both by leaching (as is the case of nitrate, which is highly soluble) and by erosion (as is the case of phosphorus), which high rates are well documented for sugarcane and reinforced by our high measures of turbidity and water hardness. Although it may seem somewhat surprising that in this study TN or TP did not show a univariate increase with land-use intensity, 9 of the 10 highest values of TP and of TN were found in converted land, and landuse intensity was accompanied by a sharp increase in water conductivity.
The amphibian metacommunity exhibited strong coherence along the first axis of ordination, implying that, overall, species responded to the same environmental gradient. This environmental gradient ranged from larger water bodies with longer hydroperiods, richer in aquatic and semiaquatic vegetation and organic sediment, and containing water with lower conductivity, turbidity, nutrient concentrations, and phytoplankton standing crop, to water bodies with the opposite characteristics. This environmental gradient also ranged from water bodies embedded in preserved landscapes to those embedded in landscapes with high intensity of land use that, as hypothesized, was accompanied by decreasing species richness.
One might expect that intensification in land use would be associated with a cumulative pattern of species losses such that communities in intensively managed landscapes would be a subset of communities in moderately intensively managed landscapes, which, in turn, would be a subset of communities in preserved, native habitats. In operational terms, this would be the case if the metacommunity exhibited positive boundary clumping and negative turnover (Presley et al., 2010). This pattern was not strictly observed, as we found positive boundary clumping but positive turnoveralbeit nonsignificantly so. Such structure corresponds to a 'quasi-Clementsian' metacommunity structure. In a strictly Clementsian metacommunity structure (i.e. with significant positive turnover), multiple species show coincident range boundaries because of shared evolutionary history and interdependent ecological relationships. Therefore, groups of adjacent sites across the latent environmental gradient have shared species compositions termed 'compartments' (Lewinsohn et al., 2006). We identified two compartments. The first comprised 13-15 sites in all land use types and across the full range of hydroperiods that presented comparatively low conductivity, turbidity, and primary productivity, and up to 17 amphibian species. The second comprised eight water bodies in converted land, characterized by being nonpermanent, mostly devoid of vegetation, turbid, and productive. This compartment contained ten species including Physalaemus cuvieri, P. centralis, Eupemphix nattereri, Dendropsophus minutus, Scinax spp, Leptodactylus fuscus, L. mystacinus, and Dermatonotus muelleri. Five to seven species occurred in both compartments.
As is frequently observed in Clementsian metacommunity structures, many species ranges were very broad and had boundaries clumped at the ends of the gradient. This could be a result either from the sampled environmental gradient being too short, or from species presenting broad niches. Considering the remarkable breadth of lentic habitats sampled, we conclude that the niches of several species are indeed very broad, including their tolerances to land-use and land management practices. The range of Leptodactylus fuscus, Scinax spp (mostly S. fuscovarius), P. cuvieri and D. minutus comprised no less than 92% (counting embedded absences), 85%, 77%, and 73% of the latent gradient. This is in strong contrast with species such as Odontophrynus cultripes and Hypsiboas lundii, each found in a single site and therefore only 4% of the gradient.
In contrast to patterns in prey distributions, land use had no effect on predatory invertebrate faunasat least at a crude taxonomic resolution (i.e. family level). All three sampled bug families, all two beetle families and two of four dragonfly families were found across the gradient; the two remaining dragonfly families (Corduliidae, Gomphidae) were found in at least one converted land use type. The number of fish families in the cerrado (8) was much larger than that in pastures (4) and sugarcane fields (1), but this effect is likely due to some water bodies in native habitats being seasonally linked to lotic water bodies, and not to land use.

Potential drivers of biodiversity change
Considering that most changes in community composition occurred in the transition from cerrados and seasonal Atlantic Forests to pastures, and not in the transition from pastures to sugarcane fields, it is reasonable to hypothesize that habitat loss, habitat fragmentation and habitat splitthat is, the spatial segregation between adult and larval habitatare important contributors to species losses in productive land (Becker et al., 2007). Indeed, organisms with complex life cycles appear to be particularly vulnerable to land conversion not only because integrate environmental degradation in aquatic habitats as larvae, and in terrestrial habitats as juvenile and adults (each of which alone could contribute to population declines), but also because the spatial segregation of adult and larval habitat obliges adults and metamorphs to risky migrations.
Differences in land management between pastures and sugarcane fields also are likely contributors to changes in local species compositions and/or abundances. Sugarcane fields are subject to intense physical disturbance every year (the harvesting) and every 5-6 years (the field reform), plus the occasional pulses of chemical disturbance following the application of pesticides and fertilizers. Small ponds, puddles, and ditches in sugarcane fields were typically devoid of dragonfly larvae; we observed the same pattern in soybean fields, where dragonfly mortality was tightly synchronized with pesticide application (D. Negri and L. Schiesari, pers. obs.). Some species of tadpoles were frequent and apparently endured pesticide applications in both sugarcane and soybean fields, but we witnessed tadpole dieoffs in ponds adjacent to sugarcane fields which could plausibly be a result of pesticide application. Indeed, Moutinho (2013) found that glyphosate, ametryn, and acetochlor, three of the most common herbicides applied in sugarcane fields, caused mortality of frog larvae when manipulated at application doses recommended by the manufacturer. Likewise, ammonium concentrations measured in ditches embedded in sugarcane plantations reached concentrations sufficient to cause larval amphibian mortality (Ilha & Schiesari, 2014), especially at earlier developmental stages (Bellezi et al., 2015).

Implications for the conservation of aquatic and semiaquatic biodiversity in managed land
Increasing environmental concerns of further conversion of native habitats to agricultural land, and of competition between biofuel and food crops for land, led independent researchers (e.g., Lapola et al., 2010), NGOs (WWF, 2009) and governmental agencies (Manzatto et al., 2009) to argue for an expansion of sugarcane over pastures. While we agree that this is the best conservation option face to an expansion of the agroindustrial production of biofuel feedstocks, it is nonetheless important to acknowledge that large-scale conversion of pastures into sugarcane fields following prevailing land management practices is expected to lead to further terrestrial habitat homogenization, and increased siltation and eutrophication of freshwater systems. Expansion of sugarcane fields over pastures is also expected to lead to increased environmental contamination with pesticides (Schiesari & Grillitsch, 2011). Assuming conservation of wetlands and surrounding buffer strips, our data suggest no substantial effects of further conversion of pastures into sugarcane fields on amphibian frequency of occurrence, density, per-unit-area biomass, and alpha diversity. Likewise, conversion of pastures into sugarcane fields appears to have no effect on predator abundance or per-unit-area biomass. However, a more careful evaluation of the effects of land use on predator alpha diversity awaits for a species-level analysis, which is underway.
The observation that 100% and 73% of the water bodies in sugarcane plantations contained tadpoles and predatory invertebrates, respectively, and that tadpole and invertebrate predator densities and per-unit-area biomass were not affected by land use imply that (i) lentic freshwater habitat is limiting for organisms with complex life cycles in this markedly seasonal landscape but also that (ii) sugarcane fields are permeable to at least part of the freshwater community. It is important to realize, however, that distribution of biodiversity in sugarcane fields was highly uneven: three large (0.3-10.7 ha) marshes contained 11 of the 14 species of amphibians found in sugarcane fields, six of which found exclusively there, and 82% of all individuals sampled. Correspondingly, they contained all eight families of invertebrate predators found in sugarcane fields, two of which found exclusively there, and 94% of all individuals sampled. Contrary to small ponds, farm ponds, puddles and ditches, preservation of buffer zones around streams and wetlands greater than 1 ha are mandated by the Forest Code, federal legislation regulating land use in every private property in the country. Therefore, these water bodies were at least 30 m from sugarcane fields and, not surprisingly, then, had less extreme abiotic characteristics when compared with other water bodies embedded in sugarcane fields, including lower turbidity, conductivity, and chlorophyll-a and phycocyanin concentrations (Fig. S4). These larger wetlands very likely act as regional sources of biodiversity. As an example, multiplying organismal densities by surface area, in December 2010 a 10.7-ha marsh contained no less than~34 million tad-poles [31 million Scinax spp, 2.3 million Dendropsophus minutus, and 0.5 million Physalaemus cuvieri] and~3 million larval insect predators [1.8 million dragonflies, 0.6 million bugs, 0.5 million beetles]. Part of these individuals eventually reached metamorphosis and likely dispersed across the terrestrial environment colonizing other water bodies in the landscape.
Another evident argument for the importance of larger wetlands as sources of biodiversity in sugarcane fields is the observation that large-scale burnings were until recently the standard preharvesting practice in sugarcane fields throughout Brazil. Indeed, 43% of the area in Luis Antonio, including part of the land surrounding sampled water bodies, was burned in the year preceding our surveys (data provided by CANASAT; Rudorff et al., 2010;Aguiar et al., 2011;INPE, 2014). Under such conditions, buffered streams, streamside marshes, and permanent and semipermanent wetlands are more than likely important sources of biodiversity in sugarcane landscapes, especially to amphibians, which are at the same time sensitive to desiccation and dispersal limited in comparison with dragonflies and other predatory insects.
Whereas our observations reinforce the importance of the Forest Code for biodiversity conservation, farmers compliance with these the so-called Permanent Preservation Areas is well below 30% in the most important sugarcane-producing regions in the country (Sparovek et al., 2012). Therefore, in a country with relatively strong environmental legislation but subject to deficient enforcement, Bonsucro (Better Sugar Cane Initiative) and any other certification scheme that proposes to foster the sustainability of the sugarcane sector are helpful. Bonsucro Standards demand compliance with all relevant applicable laws, including the Forest Code and those regulating the storage, delivery, and application rates of vinasse and pesticides. Bonsucro Standards also prevent greenfield expansion over High Conservation Value Areas and require development and implementation of an Environmental Management Plan taking into account endangered species, habitats, and ecosystems in the farm. Finally, Bonsucro Standards include application of fertilizers according to soil or leaf analysis and set a moderately restrictive maximum dosage of inorganic fertilizers, pesticides, water consumption and of tilled land area.
Considering that most native habitats in Central-South Brazileven if fragmented and in early phases of regenerationare contained within private properties, enforcement of current legislation and recognition of better land management practices by means of certification systems are clearly important for the conservation of biodiversity and, more broadly, the environmental sustainability of biofuels.

Supporting Information
Additional Supporting Information may be found in the online version of this article: Table S1. Univariate analyses for the effects of land use and survey on 16 environmental variables and 8 community variables measured in sampled water bodies. Table S2. Top generalized linear models predicting tadpole species richness on 27 water bodies distributed across a gradient in land use intensity comprising native habitats, pastures, and sugarcane plantations. Table S3. Scores from the Principal Component Analyses used to reduce dimensionality in the 18 environmental variables measured in each water body. Table S4. Scores from the Redundancy Analysis (RDA) of 25 larval amphibian communities in water bodies distributed across a gradient of environmental degradation comprising native habitats, pastures and sugarcane plantations. Figure S1. Positive relationship between tadpole richness and hydroperiod, and tadpole richness and tadpole density on 27 water bodies distributed across a gradient in land use intensity (native habitats, pastures, and sugarcane plantations). Figure S2. Biplots of the Principal Component Analyses used to reduce dimensionality in the 18 environmental variables measured in each water body. Figure S3. RDA triplot of 25 larval amphibian communities in water bodies distributed across a gradient of environmental degradation comprising native habitats, pastures and sugarcane plantations. Figure S4. A comparison of representative physico-chemical properties of small and large water bodies in sugarcane fields.