Introduction

As climate change and anthropogenic activity alter ecosystems at unprecedented rates, it has become critical to understand the consequences of biodiversity loss on ecosystem processes and the maintenance of ecosystem processes through species interactions. A complex mix of anthropogenic forces are eroding multiple dimensions of global biological diversity, including plant intraspecific, interspecific, and functional diversity (1,2). Plant diversity affects herbivore abundance and diversity, thereby influencing biomass allocation and energy fluxes between trophic levels (3). Because losses of plant diversity can destabilize the flow of resources to higher trophic levels (3,4), understanding the connection between biodiversity and trophic interactions is necessary to predict the consequences of the loss of primary producer biodiversity on ecosystem traits such as resilience to environmental perturbation (5,6). While climate change is associated with increased absolute precipitation in some regions and decreased precipitation in others, IPCC models predict mid-century increases in the frequency of extreme precipitation events in Central and South America (79), a phenomenon already observable in many ecosystems (10). Measuring how precipitation change will affect the relationship between biodiversity and ecosystem function is therefore increasingly important, particularly in diverse tropical systems.

Very few multi-site, manipulative diversity experiments have been reported from tropical areas compared to temperate environments (11), limiting our knowledge of the role of biodiversity in ecosystem function in the most species-rich regions of the planet (12). In the context of several established hypotheses (Table 1) we investigate how multiple dimensions of plant diversity affect ecosystem processes at five neotropical sites and explore how diversity modulates how ecosystems respond to changes in water availability at three of those sites. We focus on ecosystem responses that represent changes in energy fluxes between trophic levels as measured by herbivory, herbivore diversity and plant mortality. By altering plant uptake of nutrients and plant defense production, abnormal levels of precipitation can alter herbivore pressure, affecting the movement of resources into higher trophic levels (13). As such, extreme dry and wet periods of climate are expected to strongly perturb plant-insect interactions and thereby alter ecosystem function (1315). The insurance hypothesis suggests that greater biodiversity can act to stabilize ecosystems and improve their resilience to environmental change (4,16). While greater interspecific plant richness is expected to lead to increased diversity in higher trophic levels due to the accumulation of specialist herbivores, field studies have demonstrated both positive and neutral effects of interspecific plant richness on ecosystem resilience (17,18). Despite traditional views that interspecific richness has a greater impact on ecosystem processes than intraspecific diversity (19), recent research suggests that plant intra- and interspecific richness can have similar effects on ecosystem productivity and consumer abundance (20,21). As changes to intraspecific richness can alter the diversity of resources available to herbivores (22), it is necessary to consider both inter and intraspecific diversity when investigating the effects of biodiversity loss on ecosystem processes.

Path models and explanatory hypotheses.

Black arrows indicate causal paths between variables

We conducted common garden experiments at five sites across Central and South America to test the insurance hypothesis by quantifying 1) the relative strength of intra- and interspecific plant richness in driving ecosystem function and 2), the effects of increased water availability on ecosystem function. Using 33 species in the genus Piper (Piperaceae) as a model system (23), we manipulated intra- and interspecific plant richness in Costa Rica, Ecuador, Peru, and two sites in Brazil. We additionally manipulated water availability in Costa Rica, Ecuador, and Peru (Figs. 1 & S1 A). We predicted that reduced Piper diversity would lead to reduced diversity of higher trophic levels, that water addition would lead to altered herbivore pressure, and that lower Piper diversity would be associated with more extreme changes in herbivory and plant mortality in response to water addition. Finally, we predicted that changes in intraspecific and interspecific plant richness would affect ecosystem processes, including herbivory, herbivore diversity, and plant mortality, with similar magnitudes.

Treatments of intraspecific richness, interspecific richness, and water addition used in each of the five study sites.

White tiles represent treatment combinations which were not tested in a given site

Methods

Study sites and focal plant genus

We conducted a large-scale transplant experiment replicated across five sites spanning 42° latitude in the Neotropics (Fig. 1 & Table S1) encompassing a range in annual precipitation from 1271 mm to 4495 mm (Table S1). At each site, we studied herbivory on planted individuals in the genus Piper (Piperaceae) in response to experimental treatments. Study sites included lowland equatorial humid forest at La Selva Biological Station, Costa Rica; high elevation equatorial humid forest in Yanayacu Biological Station, Ecuador; high elevation equatorial humid forest at El Fundo Génova, Peru; lowland seasonally dry gallery forest in the cerrado within the phytogeographic domain of the Atlantic Forest in Mogi-Guaçu Biological Reserve, Brazil; and lowland seasonally semideciduous forest in the transition between the Atlantic Forest and the cerrado phytogeographic domains in Uaimii State Forest, Brazil. Climate classifications follow the Köppen-Geiger climate model (Table S1; 22).

Multiple species of Piper are found at all sites, ranging from 11 species in Mogi-Guaçu to 50 species in Costa Rica (25). Piper is an ideal genus for large-scale comparative studies, as it is found across the Neotropics and subtropics, from ∼10° N to about ∼32° S. In addition to being widespread, Piper is abundant and diverse across its range, encompassing ∼1000 species in the Neotropics (26). Piper has been the subject of detailed studies of herbivory, and its herbivore fauna has been surveyed across its range (2738). Members of the genus host both specialist caterpillars and beetles as well as generalist caterpillars and Orthopterans (23,2730,32). Leaf damage patterns produced by the different classes of herbivores are well documented and allow for the determination of unique taxa of herbivores (32).

Experimental design

Executing experiments across the Americas presents challenges, including nuanced variations in methodologies at each site. Nonetheless, the advantages of this expansive and consistent approach provides a greater understanding in the role of biodiversity in ecosystem function than examining isolated single-site studies. Here we describe the experimental design applied across sites, see Figures 1 & S1, and Tables S1 & S2 for site specific details.

At each study site, a factorial experiment was implemented to test the effects of plant interspecific and intraspecific richness on herbivory, variation in herbivory, and insect richness. Experimental plants were propagated from cuttings of naturally occurring Piper plants (typically with three nodes and zero to one leaves). At each site a subset of naturally occurring Piper species was selected to act as a species pool for each experimental plot (Table S2). These species were selected to constitute a breadth of genetic and functional diversity representative of the Piper community present at each site. In all sites except Uaimii, experimental plots measured 4 m in diameter and contained 12 Piper individuals, planted either in pots with locally derived soil (for experiments manipulating water availability) or directly into the ground (for experiments without water additions in Mogi-Guaçu and Uaimii; Table S1) and cultivated without fertilizer or irrigation. Interspecific richness treatment levels consisted of single species monoculture plots, two species plots, and high richness plots with the maximum number of species available to produce cuttings (Table S1). The minimum number of Piper species in a high richness plot was three in Uaimii, and the maximum was 12 in Costa Rica and Ecuador. In plots with more than one species, species were randomly sampled from the species pool for that site.

An intraspecific richness treatment was crossed with the interspecific richness treatment. Intraspecific richness was manipulated by clones taken from a single mother plant (low richness) or individual cuttings all taken from unique mother plants (high richness). In the high intraspecific richness treatment, cuttings were from different Piper individuals growing at least 10 m apart, to eliminate the likelihood of a shared root system thus representing genetically unique individuals. In Costa Rica and Ecuador, the high interspecific richness plots only included individuals from unique mother plants because all cuttings were taken from different species.

A water addition treatment was crossed with the inter- and intraspecific richness treatments in Costa Rica, Ecuador, and Peru in order to examine how IPCC predicted regional mid-century increases in precipitation in these regions (9) may affect herbivory and trophic interactions. Plants were cultivated in five to six liter pots with drainage holes. Each potted plant under the water addition treatment was watered with 2 L of water twice per month in Costa Rica and Ecuador. In Peru, 2 L of water were added to plants every three months. Water was rapidly applied as a flooding event to completely saturate the soil. Plant cultivation periods lasted in excess of one year at all sites, and as such water was added in both the wet season and dry season. Experimental plots were randomly located in the study sites; in Mogi-Guaçu and Uaimii they were organized in three replicate blocks. There was a minimum of 20 m between plots in blocks and 100 m between blocks. Randomly located plots were separated by a minimum of 50 m. Experimental periods lasted between 1.4 and 2.8 years, depending on the site (Table S1), and all measurements of herbivory were conducted on leaves that were initiated during these periods. Cuttings were replaced if they died during the first three months of the experiment. One of the authors was present at all of the sites for multiple visits to ensure as much standardization of treatment applications as possible.

Piper mortality resulted in a reduction in species richness in many plots and the loss of some treatment combinations, notably the loss of all Piper in unwatered, low intraspecific richness plots in Peru (Fig. S1B). Plots originally planted as high interspecific richness treatments had a final richness of between 1 and 12 species in Costa Rica, 9 and 11 species in Ecuador, 1 and 3 species in Peru, 2 and 4 species in Mogi, and only 2 species in Uaimii. To determine the effects of site variation in natural levels of precipitation on the outcomes of the water addition treatment, the absolute level of precipitation and precipitation anomalies relative to climate normals were collected for each month of the experimental periods in Costa Rica, Ecuador, and Peru (Fig. S2). Climate data for Costa Rica were obtained from La Selva Biological Station. Data for Ecuador and Peru were obtained using TerraClimate (39), and an interpolation error in the precipitation for Peru in February 2016, was corrected using data from the National Service of Meteorology and Hydrology of Peru (SENAMHI).

Measures of herbivory and insect richness

Plots were open to naturally occurring herbivores, and herbivory was recorded by taking photographs of all the experimental leaves at the end of the experimental period. Additional photos were taken every three months in Uaimii and Mogi-Guaçu, and in the first five months of the experiment in Costa Rica. These data were used to measure herbivory and to determine the types of herbivores feeding on leaves based on patterns of damage. When herbivores were observed on plants, they were photographed but were left to continue feeding so as to not interfere with the experiment. Herbivory was quantified by eye for each type of herbivore on each leaf by a single parataxonomist with extensive experience measuring herbivory on Piper following established protocols (32). The amount of leaf area consumed was measured in relation to the total leaf area by visually dividing the leaf into equal sized segments to determine the percent area missing. This was measured as a continuous value to the greatest possible accuracy, typically 1-5% of the total leaf area. Insect herbivores were identified to the lowest taxonomic level possible based on their damage patterns (genus for specialist Lepidoptera, family for generalist Lepidoptera, family for Coleoptera, order for Orthoptera). Direct observations of herbivores were rare, so only damage patterns were used in analyses. As insect damage patterns on leaves are tightly correlated with insect richness in tropical forests (40), the different types of damage recorded were used as proxy for the richness of above ground insect herbivores on Piper. Hereafter, the term “insect richness” refers to the richness of insect herbivore damage patterns on plants.

Data analyses

The percentage of leaf tissue consumed by each insect taxon (based on damage patterns) on individual leaves was summed to determine the total percentage of herbivory on each leaf, and the variance in herbivory was calculated as the variance in leaf damage within each plant in a plot. Due to high mortality rates across sites, Piper interspecific richness in each plot was the final number of Piper species present in each plot at the end of the experiment (rather than the number of species planted) and was analyzed as a continuous covariate rather than a categorical treatment. Because most plant deaths occurred early in the experiments, final Piper interspecific richness more accurately reflects the local plant richness experienced by herbivores. Interspecific diversity for each plot was analyzed as the proportion of the interspecific diversity present in the species pools at each site. This enabled easier comparisons between intra and interspecific richness, as intraspecific richness was only quantified as low and high, where high treatments represent the maximum intraspecific richness available at each site. The effect of intra and interspecific diversity on each response variable is reported as the change in response as diversity changes from the lowest possible value at each site to the highest.

The effects of intraspecific richness, interspecific richness, and water addition on percent leaf area consumed, the percentage of damaged leaves, variance in herbivory, and insect richness were analyzed across and within sites using hierarchical Bayesian models (HBM). This framework acts as the Bayesian equivalent of a random-effects model where site is a random effect, allowing for generalized parameter estimates across sites. Analyses testing Piper intraspecific richness and interspecific richness were conducted using data from all sites; the effects of water addition and its interactions with richness were conducted using data from Ecuador, Costa Rica, and Peru. Due to high mortality in Peru, interactions between the water addition and intraspecific richness treatments could not be modeled for measures of herbivory or insect richness, and were only tested for survival.

We acknowledge that treatment level combinations were not the same across the different sites, but experimental designs of this nature are encompassed within the framework of random-effects models, where different levels of random factors, such as site or year, consist of treatment levels that are unique to that site or year. This type of experimental design goes back to the origins of mixed models (41,42), and the lack of interactions between fixed and random effects increases generality in these models (43). Even at a single site, manipulated variables in ecological experiments do not even consist of the same manipulations across all of the units of replication, as they suffer from problems such as multiple versions of treatments, interference, and noncompliance (44)

Bayesian structural equation models (BSEMs) were constructed for each site using all treatments as exogenous variables and insect richness and herbivory as endogenous variables. Three path models were developed under the assumptions that 1) intraspecific and interspecific richness may influence herbivory both directly and indirectly by modulating insect richness 2) water addition may influence both herbivory and insect richness and 3) insect richness may influence herbivory directly (Table 1). Models I-III were tested in Ecuador, Costa Rica and Peru. while a model without the water addition variable was used to analyze data from Mogi-Guaçu and Uaimii. Additional models incorporating interactions between intraspecific richness, interspecific richness, and water availability were tested in Costa Rica, Ecuador, and Peru.

Models were run at the leaf level for herbivory and insect richness, and at the plant level for the percentage of leaves with damage, variance of herbivory, and Piper survival. For all HBMs and BSEMs, model convergence was estimated visually using traceplots and an RLJ discriminatory threshold of 1.1 (45). Model fit was determined via posterior predictive checks (PPC) using the sum of squares of the residuals as the discriminatory function (46). A PPC near 0.5 indicates a high model fit, while values near one or zero indicate poor fit. BSEM models were further compared using the Deviance Information Criterion (DIC).

All Bayesian models were written in JAGS via the jagsUI package in R (47) using Markov chain Monte Carlo (MCMC) sampling with weakly-informative priors. Residuals were modeled as normally distributed based on PPC comparisons between models. Models using gamma distributions for herbivory and binomial distributions for damage presence and mortality were found to consistently underestimate the magnitude of variance in the data based on PPC. For the majority of models, MCMC runs were conducted for 10,000 iterations using the first 1,000 iterations as a burn-in phase to generate posterior distributions of parameter estimates for each response variable. HBMs modeling interactions required 20,000 iterations with the first 5,000 as burn-ins for all models to consistently converge. Mean parameter estimates and 95% credible intervals (CIs) were calculated for all responses. 95% CIs which do not cross the y-axis are typically associated with less than a 2.5% type S error rate (48). A posterior probability of direction (PD) was calculated based on the percentage of the posterior distribution responding in the same direction as the median response. A PD of 95%, for example, indicates that the same direction of response (e.g. a positive or negative response) was observed in 95% of iterations, regardless of the magnitude of the response (49). PDs less than 95%, indicate lower confidence that a relationship exists, but can still be interpreted as the probability that an effect exists.

To analyze the effects of the experimental treatments on plant survival, survivorship curves were constructed for all sites. Analyses of Piper survival were based on the initially planted interspecific richness treatment of each plot because mortality occurred early in the experiment. The effects of water addition, intra- and interspecific richness and species identity on Piper survival were analyzed using Cox proportional hazard models. All data were analyzed using R software v4.1.0 (50).

Results

Considerable variation in herbivory and plant mortality was observed among study sites. Percent herbivory was lowest at Mogi where only 9% of leaf tissue was consumed by herbivores compared to 22% of leaf tissue consumed in Uaimii. Piper mortality was highest in Peru (89%), likely due to El Niño related drought, and was lowest in Ecuador (27%; Fig. S3A). We identified herbivore damage from 13 insect taxa at the five study sites, as well as damage by leaf miners of unknown orders. Damage from a total of 10 taxa were observed in Costa Rica, Ecuador, and Peru, 9 taxa were observed in Mogi, and 8 taxa were observed in Uaimii (Fig. S3 B). The majority of taxa in Peru, Mogi, and Uaimii were generalist herbivores, while 60% of taxa in Costa Rica and 58% of taxa in Ecuador were Piper specialists. The proportion of leaf tissue consumed by specialists was only greater than generalist damage in Costa Rica.

Our experiments revealed pronounced heterogeneity in ecosystem responses to water availability and Piper diversity between sites (Figs. S4 & S5). Posterior predictive checks (PPCs) for all hierarchical Bayesian models (HBMs), and for models I, and III were within 0.03 of 0.5, indicating models fit well. Model III (Table 1) was selected as the most parsimonious causal model for BSEMs in Costa Rica (PPC = 0.499), Ecuador (PPC = 0.499), and Peru (PPC = 0.498). Fit was high for models in Mogi (PPC = 0.5) and Uaimii (PPC = 0.497). Across all sites where the water addition treatment was applied, percent herbivory was 4.2±3.6% (mean ± 95% CI) lower in plots that received additional water (PD = 98.7%; Fig. 2 A-B). Greater Piper interspecific richness was associated with a 15±18.6% increase in the richness of insect herbivores (PD = 95.0%; Fig. 2E) and an indirect increase in herbivory was mediated by insect richness (Fig. 3). Insect richness was associated with an 8.8±2.8% increase in herbivory per insect taxon present (PD = 100%) and a 6.7±6.9% increase in the percentage of leaves with any damage (PP = 97%) (Figs. 2A, B). Intra- and interspecific richness affected herbivory, although effects varied in strength and direction across sites (Fig. 3, Table S3). Intraspecific richness had similar or greater effects on plant mortality and insect richness than interspecific richness. However, intra- and interspecific richness often had contrasting directions of effect on insect richness and measures of herbivore pressure (Figs. 3, S4B-C & S5B-C). For example, in Costa Rica insect richness was 6.0 ± 5.9% (PD = 93.7%) lower in plots with high interspecific plant richness, while high intraspecific richness increased insect richness by 7.2 ± 6.9% (PD = 98.1%). In contrast, in Ecuador high interspecific richness was associated with a 43.6 ±6.4% (PD = 100%) increase in insect richness, and high intraspecific richness decreased insect richness by 13.6 ±5.2% (PD = 100%; Fig. S5B-C).

Hierarchical Bayesian model parameter estimates for the effects of water availability, as well as intraspecific and interspecific Piper richness on herbivory, herbivorous insect richness, and Piper survivorship across sites

Violins represent the cross-site posterior parameter distribution for each relationship in site level hierarchical Bayesian models. Black lines represent the median posterior estimate and white bars represent 95% credible intervals. Percentages above violins indicate the probability of an effect being positive or negative (as indicated by a negative probability) in response to an increase of the independent variable. Distributions for water addition compare watered and control plots; distributions for interspecific richness compare Piper species richness standardized as the proportion of the maximum richness used at a site; distributions for intraspecific richness compare low and high intraspecific richness treatments; distributions for insect richness compare responses per insect taxon present on an individual leaf

Direct and indirect effects of plant diversity and water availability on insect herbivores at five study sites.

Bayesian structural equation models comparing effects of different drivers of herbivorous insect richness and herbivory at five sites. Standardized path coefficients are means of the posterior distribution for the effects estimated at each causal path. Positive relationships are indicated in blue with triangular heads, and negative relationships are indicated in red with circular heads. Black arrows indicate path coefficients of zero magnitude. Dot plots summarize the standardized mean of the posterior distribution for each causal path with 95% and 80% credible intervals. Asterisks indicate causal paths where the probability of an effect being positive or negative is >95%

The effects of water addition were altered by Piper intra- and interspecific richness at all sites (Fig. 4). Water availability reduced herbivory in Costa Rica (4.7 ±2.5%, PD = 100%) and Peru (5.1 ±3.7%, PD = 100%), but this effect was only present in Ecuador when intraspecific richness was high (Fig. S6). Across sites, water addition had negligible or negative effects on insect richness at low interspecific richness, but this pattern was reduced or reversed when interspecific richness was high. Insect richness increased by 20.8 ±18.9% more in interspecifically diverse plots water compared to unwatered plots in Costa Rica (PD = 99%). Insect richness increased by 15.2 ±16.8% more in rich, watered plots in Ecuador (PD = 96%), and 29.3 ±31.5% more in Peru (PD = 98%). Additionally, water addition had a negligible effect on insect richness in Ecuador when intraspecific richness was low, but increased insect richness by 12.1% when intraspecific richness was high (PD = 99%; Fig. 4B).

Interactions between intraspecific or interspecific richness and water availability on herbivory, herbivorous insect richness, and Piper survival

Bars indicate mean response and standard error of the mean. Percentages above each sub-figure indicate the probability that the two slopes are different as calculated using site level hierarchical Bayesian models. Due to high mortality in Peru, interactions between water and intraspecific richness could not be compared for any responses except mortality. Single and two species richness plots have been combined for visualization purposes only

Water addition had a negative effect on Piper survival in Costa Rica when intraspecific richness was low, but improved survival by 12.1 ±8.5% in high intraspecific richness plots (PD = 99%; Fig. 4C). Interactions with water in Peru may have been influenced by an El Niño related drought which resulted in high Piper mortality, while the typically wetter sites in Costa Rica and Ecuador experienced greater precipitation than average (Fig. S2).

Analysis of plant survival

Water addition did not affect Piper survival in Costa Rica (z = -1.2, P > 0.2), or in Ecuador (z = 0.2, P > 0.8), but survival was reduced by 48% in plots without water in Peru (z = 3.2, P = 0.001). Intraspecific richness increased survival by 38% in Costa Rica (z = -3.3, P < 0.001) and by 32% in Ecuador (z = -4.9, P < 0.001), but had no effect on survival in Peru (z = - 1.55, P > 0.1), Mogi (z = -0.22, P > 0.8), or Uaimii (z = 0.18, P > 0.8). Interspecific richness had no effect on survival in Costa Rica (z = -1.5, P > 0.1), Ecuador (z = -0.56, P > 0.5), Peru (z = - 0.49, P > 0.6), or Uaimii (z = 0.21, P > 0.8). In Mogi, survival was reduced by 13% in plots with higher interspecific richness (z = -2.16, P = 0.031). There was an interaction between intraspecific richness and water addition in Costa Rica and Peru. Survival in Costa Rica increased in response to water in high intraspecific richness plots, and decreased in response to water in low intraspecific richness plots (z = 3.8, P < 0.001), while the opposite pattern was observed in Peru (z = -2.4, P = 0.02; Fig. S7). Although plant die-offs cause a loss in richness in some plots, plant species identity was not related to survival with the exception of P. peltatum, which had the lowest survival rate of any species planted in Costa Rica. Further statistical results are available in the supplement.

Discussion

Our results demonstrate two key patterns. First, the strength of effects of intraspecific richness on higher trophic levels is comparable to that of interspecific richness, supporting our predictions and corroborating recent studies demonstrating the importance of intraspecific richness (20,21,51). However, the direction of the effect of intraspecific richness on herbivory, insect richness, and plant survival can vary dramatically from the direction of effect of interspecific richness, in contrast to our predictions. Second, we found that perturbations in water availability can have complex effects on herbivores and plant survival, and that these effects can be modulated by plant diversity in a context sensitive manner. While our prediction that water availability would influence herbivory across sites was supported, our results suggest that biodiversity loss and climatic perturbations may have dramatically different effects on ecosystem function at local scales, which may diminish our ability to predict how local communities will change as anthropogenic stressors increase.

As we did not directly measure plant stress or nutrition, it is difficult to determine the exact mechanism through which water addition reduced herbivory. The presence of an El Niño weather pattern during the course of the Peru experiment may have led the water addition treatment to relieve plants from drought stress, while water addition in Costa Rica and Ecuador may have added to water stress in treated plants as these sites received an above average level of rainfall during the course of the experiment. Despite this, water addition consistently suppressed herbivory in Costa Rica, Ecuador, and Peru under natural levels of Piper diversity, suggesting that predicted increases in precipitation in the next century (52) will dramatically alter the flow of resources from primary producers. Although we were only able to record the richness of insect herbivory patterns, this measure is indicative of the functional diversity of insect herbivores on Piper and changes to this value represent changes in interactions between Piper and higher trophic levels (32,40). As such, the additional reduction in effects of water addition on insect richness when Piper richness was low suggest that biodiversity loss in tropical systems will alter the ability of higher trophic levels to respond to environmental perturbations.

While our prediction that increased Piper interspecific richness would lead to increased insect diversity was met in the majority of sites, interspecific richness was associated with decreased insect diversity in Costa Rica and Peru. As herbivore taxa can be differently affected by manipulations of diversity (53), variation in the direction of the effect of intra- and interspecific richness may be due in part to variation in the composition of insect communities and herbivore pressure measured across sites. Changes in neighborhood effects when small numbers of plant species dominate a community can lead either to the reduction or increase of herbivore pressure, dependent both on the nature of plant species lost and herbivore species present in a given community (54). As such, local variation in community composition has the potential to greatly alter the effects of both climate change and biodiversity loss on ecosystem function.

Although experimental methods varied somewhat between study sites, this cannot fully explain the level of heterogeneity observed in ecosystem response. For example, the methods employed in Costa Rica and Ecuador were nearly identical, and yet the directions of effect of intra- and interspecific diversity on insect richness were reversed in these sites (Fig. S5 B-C). There was considerable variation in both biotic (Fig. S3 B) and abiotic factors (Table S1) across sites, which may have contributed to the heterogeneity observed in ecosystem response to Piper diversity. Regardless of how biodiversity loss effects ecosystem function at large scales, variation in abiotic and biotic factors at locals scales can alter these effects, reducing our ability to predict how anthropogenic activity will alter ecosystem function. This is especially relevant in tropical systems, which have been the subject of far fewer studies of ecosystem function than temperate ecosystems. As such, our knowledge of local effects on the relationship between biodiversity and ecosystem function remains limited in these systems.

A long-standing question in ecology has been the extent to which ecosystem function increases with biodiversity and if this relationship plateaus at a level past which ecological redundancy predominates. Recent results from less complex temperate grassland systems suggest that these ecosystems can be described by a mostly linear relationship between richness and function, where even rare species make unique contributions to ecosystem function (55). In these systems, high contingency can be expected, where ecosystem-level effects will depend on most of the interacting species. In contrast, we might expect that diverse, tropical communities could be characterized by greater ecological redundancy and thus be subject to less contingency (56,57). Despite these expectations, our results demonstrate heterogeneity in ecosystem response to changes in both intraspecific and interspecific richness in five tropical sites, suggesting that complexity in these systems may not reduce the contingency effects of biodiversity loss. Understanding the impact of biodiversity loss in tropical forests is fundamental to our ability to conserve those systems, and our findings highlight the importance of approaching the study of changes in ecosystem function as context sensitive responses in complex ecosystems.

Acknowledgements

The authors would like to thank Hacienda San Isidro, Ecuador for access to their forest reserve for experiments, the Meteorological Service from Peru (SENAMHI) for access to their weather data, and IEF, Minas Gerais, for access to Uaimii Forest Reserve during experiments. The collecting permit in Peru was issued by Peruvian Ministry of Agriculture (RDG N 288-2015-SERFOR-DGGSPFFS).

Funding sources

Funding was supplied through NSF grants DEB-1442103 to LAD, LAR, TLP, AMS, and DEB-1442075 to EJT

Author contributions

L.A.R., L.A.D., A.M.S., T.L.P., E.J.T., J.O.S., M.J.K., and T.J.M. developed the original idea for the research and secured funding. L.A.S., T.J.M, L.A.D, E.J.T., and H.G.L. collected Peru data. H.G.L. and L.A.D. collected Costa Rica data. W.R.S., T.R.W., and L.A.D. collected Ecuador data. T.J.M and A.R.N collected Mogi data. Y.A., L.B., T.J.M, and A.R.N. collected Uaimii data. W.R.S measured herbivory for all sites. A.J.G. analyzed the data and wrote the first draft of the manuscript. All authors contributed to analysis ideas and to writing revisions.

Competing interests

The authors declare no competing interests