Sampling design may obscure species–area relationships in landscape-scale field studies Research

We investigated 1) the role of area per se in explaining anuran species richness on reservoir forest islands, after controlling for several confounding factors. We also assessed 2) how sampling design affects the inferential power of island species–area relationships (ISARs) aiming to 3) provide guidelines to yield reliable estimates of area-induced species losses in patchy systems. We surveyed anurans with autonomous recording units at 151 plots located on 74 islands and four continuous forest sites at the Balbina Hydroelectric Reservoir landscape, central Brazilian Amazonia. We applied semi-log ISAR models to assess the effect of sampling design on the fit and slope of species–area curves. To do so, we subsampled our surveyed islands following both a 1) stratified and 2) non-stratified random selection of 5, 10, 15, 20 and 25 islands covering 1) the full range in island size (0.45–1699 ha) and 2) only islands smaller than 100 ha, respectively. We also compiled 25 datasets from the literature to assess the generality of our findings. Island size explained ca half of the variation in species richness. The fit and slope of species–area curves were affected mainly by the range in island size con-sidered, and to a very small extent by the number of islands surveyed. In our literature review, all datasets covering a range of patch sizes larger than 300 ha yielded a positive ISAR, whereas the number of patches alone did not affect the detection of ISARs. We conclude that 1) area per se plays a major role in explaining anuran species richness on forest islands within an Amazonian anthropogenic archipelago; 2) the inferential power of island species–area relationships is severely degraded by sub-optimal sampling designs; 3) at least 10 habitat patches spanning three orders of magnitude in size should be surveyed to yield reliable species–area estimates in patchy systems.


Introduction
The species-area relationship (SAR) is the earliest and best-documented pattern in spatial ecology (Rosenzweig 1995, Tjørve et al. 2018). The 'obvious fact that the larger the area taken the greater the number of species' (Arrhenius 1921) has led to a number of 108 refinements to understand such a pattern. Many mathematical models (Tjørve 2003, Triantis et al. 2012 in tandem with several types of sampling schemes (Scheiner 2003) and concurrent underlying mechanisms have been invoked to explain SARs (Hill et al. 1994).
SAR models were firstly developed in the 1920s for contiguous areas (i.e. mainland SAR; Arrhenius 1921, Gleason 1922. Accordingly, larger areas are more species-rich because they have more individuals and contain a wider spectrum of habitats (Rosenzweig 1995). Thus, given a random abundance distribution, the larger number of individuals recorded over larger areas should imply more species (i.e. sampling effect; Hill et al. 1994). Meanwhile, the greater variety of habitats encompassed by larger areas supports species restricted to specific habitats and those requiring a combination thereof (i.e. habitat diversity effect; Connor and McCoy 2001).
Subsequently, Wilson (1961) showed that the rate of species loss as a function of area reduction is higher for archipelagos (i.e. island SAR) than for contiguous areas. This occurs because the number of species on islands is also mediated by the dynamic of extinction and colonisation as postulated in the island biogeography theory (MacArthur and Wilson 1963Wilson ,1967, a paradigm also attributed to E. G. Munroe in 1948 (Tjørve et al. 2018). Accordingly, larger islands have larger population sizes resulting in lower extinction rates (i.e. area effect), and islands closer to a mainland source of species experience higher immigration rates (i.e. distance effect). Less isolated larger islands are therefore more species-rich than more isolated smaller islands (Fig. 5 in MacArthur and Wilson 1963).
Since any reduction in island area depresses species richness more than a similar reduction in contiguous areas, mainland and island SAR can be seen as extremes of a continuum that is determined by the 'islandness' (i.e. functional connectivity) of surveyed areas (Rosenzweig 1995). Such a property is arguably mediated by the dispersal ability of any given species group and the hostility of the intervening matrix in patchy systems (Bueno and Peres 2019). For example, in a forest archipelago induced by a hydroelectric dam in French Guiana, raptors were more prone to move between islands by traversing the water matrix than small mammals (Cosson et al. 1999), so the same archipelago is more functionally connected for raptors than for small mammals. Meanwhile, forest fragments surrounded by cattle pastures (i.e. less hostile matrix; Lees and Peres 2008) experience lower rates of bird species loss as a function of area reduction than forest islands within a hydroelectric reservoir (i.e. more hostile matrix; Bueno et al. 2018). Collectively, this means that islands and mainland areas will converge in their SARs at lower levels of landscape 'islandness'. Moreover, matrix type (Kennedy et al. 2010), history of disturbance (e.g. clear-cut or burned forests; Stouffer and Bierregaard 1995), time since habitat patch isolation (Jones et al. 2016), and direct human disturbance (e.g. hunting pressure; Canale et al. 2012) all mediate the number of species in habitat remnants embedded in human-modified landscapes.
Even though positive SARs appear to be ubiquitous (Connor and McCoy 2001), some studies have found a non-significant or even negative relationship (Lövei et al. 2006), with smaller patches harbouring more species than larger ones. Such unexpected results can emerge for several reasons. First, surveyed patches often span a modest size range (Watling andDonnelly 2008, Lion et al. 2014) and are, therefore, exposed to the 'small island effect' (Wang et al. 2018), where a modest variation in patch size does not affect species richness (Lomolino and Weiser 2001). Second, few patches are surveyed, thereby reducing SAR model fits (Triantis et al. 2012). Finally, the species assemblage under consideration includes both habitat (e.g. forest dwellers) and non-habitat affiliated species (e.g. matrix dwellers), resulting in compensatory dynamics whereby any loss of the former is either compensated for (Russildi et al. 2016) or exceeded (Lövei et al. 2006) by any gain of the latter.
Here, we investigated 1) the role of area per se in explaining anuran species richness on Amazonian forest islands induced by river damming, after controlling for several confounding factors (Table 1). We also assessed 2) how sampling design affects the inferential power of island species-area relationships (type IV curve sensu Scheiner 2003) aiming to 3) provide guidelines to yield reliable estimates of areainduced species losses in patchy systems. We took advantage of passive acoustic monitoring (Deichmann et al. 2018) to survey anurans at a large number of forest islands (n = 74) covering a wide size range (0.45-1699 ha) within the Balbina Hydroelectric Reservoir in central Brazilian Amazonia. We used anurans as a model group because each species has a distinct, simple and relatively stereotyped vocalisation, thereby permitting reliable species identification even in megadiverse regions (Marques et al. 2013, Ribeiro et al. 2017. Moreover, anurans generally show pronounced site fidelity and limited dispersal ability (Smith and Green 2005), allowing us to largely control for distance effects (Lima et al. 2015, Palmeirim et al. 2017. As a vast 'real-world' experimental landscape, the Balbina forest archipelago is a unique setting to examine habitat area per se effects on species assemblages because 1) it provides over 3500 replicated forest islands varying widely in size (Benchimol and Peres 2015a); 2) all forest islands were created simultaneously 28 yr ago (Fearnside 2016), having therefore been subjected to an uniform and relatively long relaxation time; 3) the open-water matrix is equally hostile; 4) forest islands span similar elevations, are restricted to upland forest and lack perennial streams, ultimately reducing habitat diversity; 5) adjacent control sites in undisturbed continuous primary forest are widely available; and the 6) de facto protection from any human disturbance covering most of the archipelago.

Study area
This study was carried out within the vast Balbina Hydroelectric Reservoir (BHR) and adjacent areas of continuous forest, located in central Brazilian Amazonia (1°40′S, 59°40′W; Fig. 1). The BHR was formed in 1987 by the damming of the Uatumã River, a tributary of the Amazon River, and covers an area of ca 300 000 ha (Fearnside 2016). The aftermath of dam construction created over 3500 islands (Benchimol and Peres 2015a) derived from former hilltops of the once primary continuous forest. To offset the environmental impacts of the dam, 938 720 ha were set aside as the Uatumã Biological Reserve (IUCN category Ia; Fig. 1b), the largest biological reserve in Brazil. Moreover, the left bank of the former Uatumã River, including all islands, has also been effectively protected (Fig. 1b).
The vegetation is characterised by submontane dense ombrophilous (terra firme) forest, although seasonally flooded igapó forest formerly occurred along the margins of the Uatumã River before damming. Islands span a wide range in size (0.2-4878 ha; Benchimol and Peres 2015a), virtually all of which lack perennial streams because lowland areas were submerged following the rise of floodwaters. Forest structure in larger islands resembles that of the continuous forest with a higher dominance of large-seeded and canopy tree species, whereas smaller islands are dominated by pioneer tree species due to unavoidable edge-mediated forest disturbance (Benchimol and Peres 2015a). According to the Köppen classification, the climate is equatorial fully humid (Af ), with mean annual precipitation and temperature of 2464 mm and 26.5°C, respectively (Alvares et al. 2013). Table 1. Confounding factors affecting area per se effects on species richness and how they were controlled for in this study. Note that, as an observational study exploiting a natural field experiment, confounding factors could not be entirely removed but were minimised to a large extent.

Factor
How the factor was controlled for Sampling effect We used the rarefied number of species as our response variable, rather than the observed number of species Habitat diversity effect Forest islands resulted from the rise of the reservoir floodwaters. Because lowlands are invariably flooded, only upland areas lacking streams persist in any one isolate Distance effect We focused on anurans because they show high site fidelity and limited dispersal ability. Also, the hostility of the open-water matrix further hampers anurans from moving across forest islands Matrix, history, time since isolation, direct human disturbance All forest island islands are surrounded by a lentic-water matrix and were created at the same time (1987).
The study region has a low human population density and negligible direct anthropogenic impact, and most islands are within a large strictly-protected area Species assemblage Only forest species were recorded across surveyed sites -i.e. species with 'Forest -Subtropical/Tropical Moist Lowland' listed as a 'Suitable' habitat according to IUCN (2018)

Sampling design
We surveyed 151 plots at 78 sites, including 74 islands and four continuous forest sites (Fig. 1c). We attempted to survey a similar number of plots in riparian (i.e. along streams) and non-riparian habitats (i.e. away from streams) within continuous forest, and all available riparian habitats on islands, but only two islands had streams. Accordingly, we surveyed 13 riparian and 10 non-riparian plots in continuous forest, and 4 riparian and 124 non-riparian plots on islands (Fig. 1c). The number of plots per survey site was defined according to island size and presence of streams and varied from 4 to 10 in continuous forest sites and from 1 to 7 on islands ( Fig. 1c; Supplementary material Appendix 1 Table A1). The overall study meta-landscape encompassed 253 951 ha in which plots were spaced apart by an average distance of 32.63 km (SD = 18.83 km; range = 0.19-84.60 km; Fig. 1c).

Anuran surveys
We surveyed anurans from July to December 2015 using autonomous recording units (ARUs) developed by the Automated Remote Biodiversity Monitoring Network (ARBIMON, < https://www.sieve-analytics.com >). Each ARU consists of an LG smartphone enclosed within a waterproof case with an external connector linked to an omnidirectional microphone. At each of the 151 plots, we deployed one ARU attached to a tree trunk 1.5 m above ground with the microphone pointing downward. ARUs were left unattended at each plot for five consecutive days and programmed to record 1 min in every 5-min interval using the ARBIMON Touch application. All recordings are archived at the ARBIMON II web platform and are freely available at < https://arbimon.sieve-analytics.com/project/balbina >.
We selected a subset of 62 1-min recordings per plot (n = 151) to identify all anuran species occurring therein, totalling 9362 1-min recordings. These recordings were derived from the following schedule: the first 1-min recording segment every 10 min over a 5-h period (from 17:00 to 22:01) during sample days 2 and 4. Anuran species were identified by GSM who inspected all the recordings both aurally and visually using the ARBIMON II Visualizer tool. Species identifications were validated thereafter by ILK as a procedure to ensure accuracy. During this validation procedure, species records had to be either readily identified by hearing or inspecting the sonograms to be included in the analysis. Species records were discarded if they could not be readily identified and/or if clearly audible sonogram acquisitions were inadequate (e.g. faint vocalisations too far away from the microphone).

Caveat on anuran surveys
We acknowledge that some species may have been missed during surveys for several reasons. First, species may have not vocalised over the two sampling days at each survey plot. Second, complementary methods are required to yield complete species assessments (Condrati 2009). Third, the number of species detected at a given site may vary across seasons and years (Menin et al. 2008). Thus, all surveys should be ideally carried out using complementary methods at the same time over ample sampling periods, which is virtually unfeasible.
To minimise these issues, we relied on auditory encounters, which have been shown to detect more species than visual encounters at the Uatumã Biological Reserve (Condrati 2009), and surveyed all island size categories throughout the sampling period (11 July to 4 December 2015) virtually randomly (Supplementary material Appendix 1 Fig. A1) to avoid temporal bias. We stress that our sampling design focused on favouring the spatial extent of sampling over temporal repetition, a general tradeoff faced in landscape-scale field studies (Banks-Leite et al. 2014). Most importantly, our aim was not to yield complete but satisfactory species assessments while accounting for differences in sampling effort to produce species-area curves, something still rare in studies on species-area relationships.

Response variable
Before accounting for differences in sampling effort (i.e. number of 1-min recordings) across survey sites, we inspected the degree to which the observed number of species was correlated with sampling effort. We then calculated the rarefied number of anuran species using sampling-unitbased incidence data with 1000 bootstrap replicates using the inext package (Hsieh et al. 2016) within the r software (R Core Team). To accomplish this, we created a species-bysampling-unit matrix per survey site, in which each species (row) was assigned as present (1) or absent (0) and each sampling unit (column) corresponded to a 1-min recording. We standardised the sampling effort to the statistical mode, the most frequent number of 1-min recordings across survey sites (n = 62). We did so because inext calculates both the interpolated and extrapolated number of species. Accordingly, we used the interpolated, observed and extrapolated number of species for sites with a sampling effort above (n = 33), equal to (n = 43) and below (n = 2) the statistical mode, respectively (all of which are hereafter referred to as the rarefied number of species). We also used inext to calculate the sample coverage to assess whether survey sites were sufficiently inventoried on the basis of 62 1-min recordings.

Predictor variable
Island size corresponds to the total insular forest cover and was calculated in Qgis software (QGIS Development Team 2016) using a classified image (Collection 2, 2015, Amazon) derived from 30-m resolution Landsat imagery downloaded from the Brazilian Annual Land Use and Land Cover Mapping Project (available at < http://mapbiomas.org >). Forest cover was defined as 'dense forest' (pixel value 3), because other pixel values effectively represent either heavily degraded forests or non-forest land cover types. Accordingly, the size of our 74 surveyed islands ranged from 0.45 to 1699 ha (Supplementary material Appendix 1 Table A1).

Island species-area relationships
To depict island species-area relationships (ISARs), we used the exponential equation (semi-log model; Gleason 1922) to fit simple linear regression models as follows: where S = rarefied number of anuran species, z = regression slope, A = island size (ha), c = regression intercept. In this equation, z indicates the rate of species loss as a function of island size reduction, whereas c indicates the carrying capacity per unit area (Triantis et al. 2012). Despite the fact that ISARs can be fitted with dozens of alternative models (Triantis et al. 2012), the semi-log model was chosen because it is widely used (Tjørve 2003), easy to interpret and allows the inclusion of sites at which no species was recorded (S = 0; Supplementary material Appendix 1 Table A1).
To assess the degree to which shortening the range in island size changes the fit (r 2 ) and the slope (z) of the ISAR for anurans at the BHR landscape, we first classified the survey islands into five size categories: very small (< 4 ha, n = 23), small (4-20 ha, n = 20), medium sized (20-100 ha, n = 17), large (100-400 ha, n = 7), very large (> 400 ha, n = 7). We then fitted semi-log models to islands classified as 1) very small + small + medium sized + large + very large (n = 74); 2) very small + small + medium sized + large (n = 67); 3) very small + small + medium sized (n = 60); 4) very small + small (n = 43); and 5) very small (n = 23). Because no species was exclusively recorded on the four riparian plots located on islands, anuran records from riparian and non-riparian plots were combined to produce species-area curves.

Tradeoff between replication power and extent of the gradient
Ideally, biodiversity surveys should include many sites covering a wide variation along any given gradient. However, logistical, financial or landscape (e.g. few and small habitat patches remaining) constraints, or combinations thereof, prevent attempts to adopt an ideal sampling design. With this in mind, we investigated the role of island-scale replication and the range in island size in detecting a positive ISAR for anurans at the BHR landscape. To do so, we subsampled our surveyed islands (n = 74) following both a stratified and nonstratified random selection of 5, 10, 15, 20 and 25 islands. In the stratified random selection, an equal number of islands belonging to each size category (see above) was selected to cover the full range in island size (0.45-1699 ha). Accordingly, for 5, 10, 15, 20 and 25 islands selected, each island size category was represented by 1, 2, 3, 4 and 5 islands, respectively. In the non-stratified random selection, subsets of 5, 10, 15, 20 and 25 islands smaller than 100 ha (n = 60) were selected, thereby covering a short range in island size.
We also carried out a literature review (Supplementary material Appendix 1 Data sources A1) focused on both tropical and temperate anuran studies worldwide to assess 1) how prevalent positive ISARs are at a global scale and 2) the role of the number of patches and range in patch size (largest minus the smallest) in detecting ISARs. Herein, the term 'patch' is used to encompass both 'fragment' and 'island'. Since results may be affected by the analytical approach employed (Bueno et al. 2018), we reanalysed data from each study based on our literature review using the semi-log model as described above. Note that for the global analysis, the response variable is the observed number of anuran species.

Results
Considering all 151 plots at 78 survey sites, we recorded 37 anuran species representing 18 genera and nine families (Supplementary material Appendix 1 Table A2). The most ubiquitous species was Ameerega trivittata (n = 54 sites), whereas five species were only recorded at one site (Supplementary material Appendix 1 Table A2). At the four continuous forest sites (n = 23 plots), we recorded 27 species from 15 genera and eight families; the number of species per continuous forest site ranged from 13 to 20 (mean ± SD = 15.75 ± 3.10; Supplementary material Appendix 1 Table A1). On the 74 islands (n = 128 plots), we recorded 34 species from 18 genera and nine families, and the number of species per island ranged from 0 to 21 (6.12 ± 4.46; Supplementary material Appendix 1 Table A1).

Species richness and sampling effort
The observed number of species was strongly and positively correlated with sampling effort (i.e. number of 1-min recordings; r = 0.82; Fig. 2a). However, the observed and rarefied number of species standardised to 62 1-min recordings were also strongly and positively correlated (r = 0.98; Fig. 2b). Therefore, high levels of local species packing was not artificially inflated by higher sampling effort, and sample coverage was adequate, exceeding 90% at 75 survey sites (Fig. 2c), indicating that a sampling effort per site of 62 1-min recordings was overall sufficient.

Island species-area relationships for anurans at Balbina
Island size explained ca half of the variation (r 2 = 0.49) in the rarefied number of species considering all 74 islands (Fig. 3). Regression slopes were flattened, and model fits dramatically reduced as the range in island size was narrowed down, leading to a non-significant species-area relationship for islands smaller than 4 ha (p = 0.90; Fig. 3). Importantly, the similar rarefied number of species calculated for very large islands (> 400 ha) and continuous forest sites (Fig. 3) suggests that a further increase in island size would not imply more species, as long as the sampling effort is standardised.
For the stratified random selection of islands covering the full range in island size (0.45-1699 ha), positive ISARs always held true for 15, 20 or 25 islands selected, usually for 10 islands, but only sometimes for 5 islands (Fig. 4). Slope deviances, measured as the degree to which the angle of a regression line deviates from that derived from all 74 islands, approximated one (i.e. no deviance) on average (Fig. 4a), while model fits (r 2 ) were increased to about 62% on average (Fig. 4c), regardless of the number of islands (5 to 25) selected. In none of the cases, ISARs were significantly negative.
For the non-stratified random selection of islands covering the short range in island size (< 100 ha), positive ISARs failed to hold true in the vast majority of cases, regardless of replication power; i.e. the number of islands selected (5 to 25; Fig. 4). Slope deviances were reduced in 35% to 65% on average (Fig. 4b), while model fits (r 2 ) were reduced to about 8% on average (Fig. 4d), regardless of the number of islands (5 to 25) selected. In none of the cases, ISARs were significantly negative.

Prevalence of island species-area relationships for anurans worldwide
We compiled 25 datasets from 23 anuran studies in fragmented landscapes representing 12 countries worldwide (Fig. 5). Our reanalysis of the original data using the semilog model revealed a positive ISAR for 18 datasets (mean r 2 adj ± SD = 0.45 ± 0.32) and a non-significant ISAR for seven datasets (Fig. 6). In none of the cases, ISARs were significantly negative. Remarkably, all 17 datasets that spanned a range in patch size (largest minus the smallest) larger than 300 ha yielded a positive ISAR, but the number of patches alone, which was widely variable (5 to 24), did not affect the detection of ISARs (Fig. 6).
In four instances, the results and the analytical approach from our reanalysis were in contrast with the original studies (Supplementary material Appendix 1 Table A3). The authors used a model selection approach (Bickford et al. 2010) or multiple linear regression models (Pineda and Halffter 2004, Figure 3. Anuran species-area relationships on forest islands surveyed at the Balbina Hydroelectric Reservoir landscape across five sets of islands (indicated by grey dashed lines and colour circles). The regression line for islands smaller than 4 ha (n = 23) is represented in red; up to 20 ha (n = 43) in orange; up to 100 ha (n = 60) in purple; up to 400 ha (n = 67) in blue; and up to 1699 ha in green (n = 74). Continuous forest sites (CF, black circles) were not included in the regression fits. Note that model fits (r 2 ) and regression slopes tend to be reduced as the range in island size is narrowed down, thereby decreasing the estimated impact and the explanatory power of forest shrinkage on the rarefied number of species. Islands larger than 400 ha yielded similar values of species richness as continuous forest sites, indicating that, by controlling for sampling effort, further increases in the range of island size would not necessarily increase the number of species detected on islands. Herrera 2011, Ferrante et al. 2017), whereas we used simple linear regression models.

Discussion
Our results indicate that habitat area per se plays a major role in explaining anuran species richness on forest islands induced by a large hydroelectric dam in lowland Amazonia. However, the fit and the slope of the island species-area relationship (ISAR) derived from the semi-log model [S ~ log 10 (A)] were affected mainly by the range in island size considered, and to a very small extent by the number of islands surveyed. Hence, reported failures in detecting positive ISARs (Watling andDonnelly 2008, Lion et al. 2014) could be attributed to sub-optimal sampling designs brought about inevitable realworld landscape constraints, rather than an inherent absence of area effects on species richness. (c) (d) Figure 4. Results of anuran species-area relationships (SARs) derived from random sampling of 74 forest islands (Fig. 3) surveyed at the Balbina Hydroelectric Reservoir landscape (see 'Material and methods: Tradeoff between replication power and extent of the gradient' for further clarification). Each dot corresponds to a single SAR. Slope deviance was measured as the degree to which the angle of a regression line deviates from that derived from all 74 islands (green line in Fig. 3): 1.0 indicates no deviance (i.e. same slope), and values smaller and larger than 1.0 indicate lower and higher slopes, respectively. Red circles and red lines show means and standard deviations. Box-andwhisker plots show median (at notch), lower and upper quartiles and 1.5 × interquartile ranges. Note that plots on the same row (a and b; c and d) are on the same scale to allow direct comparisons.

Area represents a supra-variable
Habitat patch size or forest remnant area can explain virtually the entire variation (r 2 = 0.99) of anuran species richness across Amazonian terra firme forest fragments, but this does not necessarily imply that area per se is the only underlying mechanism driving ISARs (Zimmerman and Bierregaard 1986). Because larger areas usually accommodate more individuals and more habitats, more species are expected to be recorded therein due to both the 1) sampling effect and 2) habitat diversity effect (Rosenzweig 1995). Accordingly, once these two effects are controlled for, ISARs tend to become weaker (lower fit and slope), relegating area per se to a lesser role in explaining overall species richness. Since the sampling effect is purely governed by the laws of probability, it should be refuted before ecological processes can be examined (Hill et al. 1994). These authors reported a strong fit (r 2 = 0.80 in a second-order polynomial regression) of the ISAR for birds in forest fragments in Ghana, which was largely attributed (r 2 = 0.16-0.37) to a sampling effect. Likewise, the shallow slope of the ISAR (log-log model) for birds in an anthropogenic forest archipelago in China was attributed to the low habitat diversity among islands, which were dominated by Masson pine (Yu et al. 2012).

Controlling for the sampling effect
Since poorly standardised sampling effort in biodiversity inventories is likely to produce misleading results, one should compare species richness among sites using either individualbased (abundance data) or sample-based (presence-absence data in a given sample) rarefaction curves (Gotelli and Colwell 2001). Here, we used a sample-based rarefaction procedure because of the nature of passive acoustic monitoring data, in which individuals recorded in consecutive samples (1-min recordings) are not independent. Accordingly, a rarefaction procedure (Hsieh et al. 2016) based on species incidence in 1-min recordings (samples) was used to divorce the sampling effect from the area effect. Meanwhile, we assessed the degree to which a standardised sampling effort of 62 1-min recordings yielded sufficient sampling effort to quantify the number of anuran species across survey sites. The fact that sampling effort was standardised to calculate species richness and that sample coverage was over 90% in virtually all surveyed sites (Fig. 2) provide robust support for an area effect that is independent of sampling effect.

Controlling for habitat diversity
In Amazonian terra firme forests, habitat diversity in terms of vegetation structure is associated with hydrological features of the terrain including elevation (Castilho et al. 2006), belowground vertical distance to the water table (Schietti et al. 2014) and horizontal distance to perennial streams (Drucker et al. 2008). Both elevation and distance from streams -two variables that are typically correlated -have been shown to shape anuran assemblages in continuous forest adjacent to our study landscape (Condrati 2009). Accordingly, low-elevation sites near streams are more species-rich and harbour a distinct anuran species composition compared to high-elevation sites far away from streams (Condrati 2009). However, our forest islands are, by definition, upland remnants induced by the flooding of lowland areas in the once continuous forest, thereby lacking perennial streams. Therefore, within-island habitat diversity associated with area is greatly reduced, therefore providing evidence of an area effect that is independent of habitat diversity.
The mechanisms underlying ecological patterns can only be properly inferred from field experiments, despite their limitations in isolating co-varying mechanisms (McGarigal and Cushman 2002). As a mensurative experiment (sensu McGarigal and Cushman 2002), the relationship between island size and habitat diversity could not be entirely removed. For example, one key feature of how habitat structure affects anuran assemblages is the presence of breeding sites (Hillers et al. 2008, Bickford et al. 2010, so a higher species richness at larger habitat patches (i.e. fragments or islands) often results from higher diversity of breeding sites (Almeida-Gomes and Rocha 2015, Almeida-Gomes et al. 2016). In Amazonian terra firme forests, peccary wallows, which are used by some species as a breeding site (Zimmerman and Bierregaard 1986), are unlikely to be found on islands smaller than 100 ha where the occupancy probability of peccaries is much lower (Benchimol and Peres 2015b), thereby decreasing the number of breeding sites. Moreover, the number of anuran reproductive modes (sensu Haddad and Prado 2005) -a proxy for breeding sites -was positively related to log-transformed island size (r 2 = 0.45, Supplementary material Appendix 1 Fig. A2), but not to the same extent as reported for in Atlantic Forest fragments of southeastern Brazil (r 2 = 0.87; Almeida-Gomes and Rocha 2015). Therefore, the effect of area per se on species richness was significant but probably weaker than that estimated in general (r 2 = 0.49 considering all 74 islands; Fig. 3).

Insular species assemblages
In fragmented landscapes, species assemblages in habitat patches are composed of relict species (that were present before fragmentation), matrix-derived species and interpatch dispersers (Watson 2002). The relative contribution of these groups depends on the type of matrix surrounding habitat patches. On the one hand, species richness in forest fragments embedded in terrestrial matrices is the result of 1) loss of relict species (Stouffer and Bierregaard 1995), 2) influx of matrix-derived species (Lövei et al. 2006) Figure 6. Prevalence of island species-area relationships (ISARs) derived from 25 anuran datasets worldwide (Fig. 5), showing how the number of habitat patches surveyed and the range in patch size (largest minus the smallest) affect the significance and fit of ISAR models. Blue and grey circles indicate positive and non-significant relationships, respectively. Circle sizes are proportional to the magnitude of r 2 -values. 'Many' and 'few' correspond to studies with more or less than 15 survey patches, respectively. 'Broad' and 'narrow' correspond to studies with survey patches spanning more or less than 300 ha, respectively. All studies covering a sufficiently wide range in patch size (> 300 ha; dashed vertical grey line) resulted in positive ISAR semi-log models.
colonisation of inter-patch disperses following matrix regeneration (Stouffer et al. 2009). On the other hand, on landbridge forest islands induced by hydroelectric dams, there are no matrix-derived species and the colonisation of interpatch disperses is largely prohibitive in the case of anurans. Therefore, we surmise that anuran assemblages inhabiting our forest islands primarily consist of a subset of relict species whose fate is determined by the area effect and aggravated by the distance effect.

The form of island species-area relationships
Regardless of the underlying mechanisms, larger areas tend to harbour more species. Combining data from two extensive syntheses on island species-area relationships (Triantis et al. 2012, Matthews et al. 2016, there were 584 out of 808 cases of significant ISARs (log-log model) and highly variable fits (r 2 ) and slopes (z) among those relationships that were significant. Many factors could account for such a variation in significance and model parameters (fit and slope). For example, the patch:matrix contrast (Kennedy et al. 2010), disturbance severity (Stouffer and Bierregaard 1995), relaxation time (Robinson 1999) and accessibility for hunters (Canale et al. 2012) can all affect the number of species in habitat patches, thereby modulating ISARs. Furthermore, the influx of matrix-derived species into habitat patches may either attenuate (Matthews et al. 2014) or even reverse (Lövei et al. 2006) the estimated impact of patch size on species richness. However, at the Balbina forest archipelago, all of these confounding factors were controlled for (Table 1), allowing us to depict unbiased patterns of ISARs for anurans in an Amazonian fragmented landscape. The form of ISARs is also affected by issues of sampling design. At Balbina, both the model fit and the regression slope of the ISAR were reduced as the largest islands were progressively removed from the analysis (Fig. 3), thereby jeopardising inferences on area-driven anuran species losses. Thus, the truncated size range of habitat patches surveyed in any give study is likely the main reason for the lack of a significant ISAR for anurans in Neotropical forests in Brazil (n = 23 fragments, range = 1.71-27.41 ha; Lion et al. 2014) and Bolivia (n = 24 fragments, range = 0.6-8.5 ha; Watling and Donnelly 2008). Accordingly, species richness below a certain threshold can vary independently of area because smaller patches are more susceptible to environmental stochasticity (i.e. the small island effect; Lomolino and Weiser 2001). Although not necessarily ubiquitous (Wang et al. 2016), the small island effect was detected for anurans in an anthropogenic forest archipelago in China under a threshold of ca 40 ha (n = 23 islands, range = 0.59-1289 ha; Wang et al. 2018). Likewise, the ISAR for anurans at Balbina was either non-significant or yielded a very weak inferential power (r 2 ≤ 0.09) considering only islands smaller than 100 ha.
Two studies in the Brazilian Atlantic Forest illustrate the role of sampling design in detecting ISARs for anurans. In the first (Almeida-Gomes and Rocha 2014), the authors surveyed 12 patches ranging from 4.7 to 272 ha and failed to detect a significant ISAR. In the second (Almeida-Gomes et al. 2016), they surveyed 21 patches ranging from 1.9 to 619 ha and subsequently detected a significant ISAR. These authors then concluded that failures to detect ISARs should be attributed to an insufficient number of patches and range in patch size. However, our results suggest that the number of patches is not as critical as the range in patch size (Fig. 4). At Balbina, 10 islands covering the full range in island size (0.45-1699 ha) yielded similar ISARs (in terms of model fit and regression slope) compared to 15, 20 or 25 islands. Conversely, the vast majority of ISARs were not significant regardless of the number of islands (5 to 25) if only a short range in island size (< 100 ha) had been sampled. Such a pattern was corroborated in our global review, which revealed that all datasets yielding a significant ISAR spanned a meaningful range in patch size larger than 300 ha, whereas all but one dataset covering a shorter range (< 300 ha) yielded nonsignificant ISARs (Fig. 6). Given the realities of field studies, the sampling tradeoff between number of sample replicates and extent of the gradient covered should therefore favour the latter to derive more reliable inferential relationships (Eigenbrod et al. 2011, Kreyling et al. 2018. This is in fact good news for field investigators who often face severely limited human, time and/or financial resources and can only survey a small number of sites.

Conclusions
We conclude that 1) habitat area per se plays a major role in explaining anuran species richness on Amazonian forest islands within one of the largest anthropogenic archipelagos on Earth; 2) the inferential power of island speciesarea relationships is clearly weakened by sub-optimal sampling designs; and 3) at least 10 habitat patches spanning three orders of magnitude in size should be surveyed to yield reliable estimates of area-driven species losses in patchy systems.