An intermediate level of disturbance with customary agricultural practices increases species diversity in Maya community forests in Belize

The study addresses how customary Indigenous swidden (slash-and-burn) agricultural practices affect biodiversity at landscape scales, where complex interactions among forest succession, landscape variables, and socio-economic factors create multiple pathways of secondary forest development. Policymakers often question the sustainability of swidden due to perceived negative impacts, yet ethnographic and ecological evidence suggests Maya agroecology can enhance biodiversity and soil fertility. Extending Bliege Bird’s model of emergent anthropogenic disturbances, the authors test the hypothesis—grounded in the Intermediate Disturbance Hypothesis (IDH) and complex adaptive systems theory—that intermediate levels of swidden disturbance are associated with higher canopy tree diversity (inferred via spectral diversity) around two Q'eqchi' Maya villages in southern Belize. The research aims to clarify links between household-level land use decisions and large-scale ecological patterns to inform conservation policy and support Indigenous land rights and livelihoods.

The paper situates Maya swidden within debates about intentional conservation versus emergent, coevolved management. Prior Mesoamerican studies show Maya forest gardens, intercropping, and rotational strategies can enhance biodiversity and soil fertility, with legacy effects detectable in modern forests near archaeological sites. Evidence from other Indigenous systems indicates that beneficial ecological outcomes can emerge without explicit conservation intent, consistent with complex adaptive systems (e.g., Balinese irrigation scheduling, Martu burning in Australia). While swidden is often criticized for degrading forests, remote sensing in southern Belize indicates less permanent forest loss where swidden is concentrated. The IDH posits maximal diversity at intermediate disturbance levels; although debated, it remains plausible in human-modified tropical systems where household decisions on field placement, frequency, and size could generate mosaics that foster niche diversity and reduce competitive exclusion. The review highlights fragmentation-diversity relationships, functional redundancy loss with land-use intensification, and calls for landscape-scale studies spanning the full disturbance gradient.

Study area: ~18,000 acres of community forests around Crique Sarco (settled 1910) and Graham Creek (1998) in the Toledo District, Belize, encompassing swidden mosaics, undisturbed lowland tropical forest, secondary forests, riparian areas, and pastures. Household surveys (2018) indicated swidden as a primary livelihood in 75% of Crique Sarco and 93% of Graham Creek households. Data collection: During April 16–27, 2018 (dry season), long-range drones flew at ~1500 ft with a MicaSense RedEdge 3 multispectral sensor (5 bands: RGB, Red Edge, NIR) to scan ~18,000 acres. The dataset comprises 52,798 georeferenced 5-band images and 630 ground truth points. Local community mapping teams provided GPS points and land-use history (e.g., last clearing year), enabling training and validation. Image processing: Photogrammetry and radiometric correction were performed in Pix4D, using calibration panel and down-welling light sensor data, producing a radiometrically corrected 5-band raster mosaic. Land-use classification: Using community GPS training polygons and a semi-automated classification workflow, pixels were labeled into classes reflecting swidden history and land use: 2018 clearings, 1 year, 2–4 years, 5–11 years, 12–19 years, 20+ years (forest), and pasture. Confusion matrices guided temporal aggregation to optimize classification accuracy. Spectral species and diversity: Following the spectral variation hypothesis and optical types concept, principal component analysis reduced dimensionality of reflectance data, then k-means clustering assigned each pixel to a spectral species (biodivMapR). Candidate cluster counts k = 20, 35, 50 were evaluated against ground data; k = 20 yielded highest correlations with field diversity (including tallest 10% trees) and was used for Shannon entropy and effective species number calculations. Class-level alpha diversity was computed as Shannon diversity of spectral species per land-use class. Landscape metrics and sampling: Fragmentation statistics (e.g., edge density, patch metrics) were computed with landscapemetrics. To capture scale-dependent processes, a hexagonal grid sampling framework was applied, with sensitivity analysis identifying 70.15 ha as an asymptotic sampling unit capturing beta diversity while ensuring sufficient samples. For each hexagon, the team computed effective spectral diversity, edge density (m edge per ha), covariates including proportion pasture (PPAST), proportion of 20+ forest (PPF), and distance to nearest village or road. Validation of spectral diversity: In May 2019, five modified Gentry plots spanned four land-use classes (two plots in 20+ y, one each in 12–19 y, 5–11 y, 2–4 y). Within each plot, five 2×50 m transects recorded all vascular plants and trees ≥2.5 cm DBH, identified taxonomically/parataxonomically/Indigenously. A total of 1,847 individuals across 139 species were recorded. Spectral diversity from corresponding drone regions was correlated with field diversity (overall and tallest 10% of individuals). Modeling: A Bayesian multilevel model (brms) related effective spectral species diversity to edge density via a quadratic (nonlinear) term to test for an IDH-like convex relationship. The model included village as a grouping factor, an ICAR spatial term for areal spatial autocorrelation, and covariates (PPAST to account for cattle pasture effects differing from customary swidden; PPF). Informative priors constrained polynomial shapes to ecologically plausible forms. Model diagnostics, prior predictive checks, and sensitivity analyses are reported in Supplementary Information. Simulation: An intermediate disturbance model simulated landscapes across disturbance gradients, fitted with a second-order polynomial to provide reference expectations for interpreting posterior predictions.

- Spectral diversity by land-use class: Effective spectral species diversity (H, effective number) increased with fallow age and was highest in older forests: 1 y = 13.5; 2–4 y = 15.3; 5–11 y = 15.7; 12–19 y = 16.1; 20+ y forest = 16.5. Pasture was lower (12.8), and recently disturbed 2018 clearings lowest (7.35). Asymptotic accumulation curves indicate adequate sampling and distinct spectral species compositions per class. - Fragmentation patterns: Edge density generally decreased with distance from village/road, with variability near settlements due to mixed land uses, pastures, and natural disturbances. Pastures exhibited higher fragmentation than 2018 clearings, reflecting management (e.g., fencing, planted grasses) and retained shade trees/palms. - IDH test (main result): Multilevel Bayesian models showed a significant convex (nonlinear) relationship between edge density and effective spectral species diversity in both villages, consistent with the Intermediate Disturbance Hypothesis. The effect was stronger in Crique Sarco (more samples with higher edge density than Graham Creek). Spatial predictions indicated that hexagons with mean observed edge density had the highest mean spectral diversity and the lowest standard deviation compared to low or high edge density conditions. - Beta diversity: High spectral beta diversity associated with intermediate disturbance levels, indicating more heterogeneous patch combinations and potential niches for rare species. - Spectral–field validation: Positive correlations between spectral diversity and field-measured vascular plant diversity: Spearman’s rho = 0.77, p = 0.072 (all species); rho = 0.83, p < 0.05 (tallest 10% of individuals), robust across k choices for clustering. Overall, the highest canopy spectral diversity occurs at intermediate levels of swidden-driven landscape fragmentation, replicated across two Maya villages.

The findings support the hypothesis that household-level swidden decisions (field placement, size, and rotation timing) can produce intermediate levels of landscape fragmentation that maximize canopy spectral diversity, aligning with the Intermediate Disturbance Hypothesis and concepts from complex adaptive systems. Emergent properties of swidden mosaics—without centralized, top-down planning—can enhance ecosystem functions by creating diverse habitats and edge conditions that aid seed dispersal, maintain seed banks via embedded mature forest patches, and offer stepping-stone habitats, collectively promoting higher alpha and beta diversity. The results strengthen the case that Indigenous swidden practices, governed by customary norms and social institutions, can contribute to biodiversity and resilience at landscape scales. They also engage ongoing ecological theory debates by demonstrating IDH patterns in a human-dominated terrestrial system using remote sensing at appropriate spatial scales. The work implies that policies supporting Indigenous land rights and customary practices may align with biodiversity conservation and climate initiatives when mosaic structures remain within intermediate disturbance regimes.

The study demonstrates that in two Q'eqchi' Maya community forests in Belize, swidden practices can maintain landscape mosaics at intermediate disturbance levels that maximize canopy spectral diversity, an emergent outcome of decentralized household decisions interacting with tropical forest dynamics. This provides a scalable, remotely sensed framework to understand and potentially support sustainable swidden systems within broader conservation and climate policies. Future research could expand replication across regions and years to incorporate temporal dynamics, integrate more extensive ground inventories to link spectral with taxonomic/functional diversity, refine modeling to incorporate species-specific growth/productivity differences, and assess how policy or market shifts (e.g., pasture expansion) alter mosaic dynamics and biodiversity outcomes.

- Biodiversity inference relies on spectral diversity from a single-season, single-year drone survey; seasonality and plant physiological status may influence spectral signatures. - Spectral species primarily reflect canopy-visible flora; sub-canopy and faunal diversity are not captured. - Limited ground-truth biodiversity: five plots spanning four classes; while correlations are positive, sample size is small (p = 0.072 for overall diversity). - The approach does not distinguish species-specific growth rates/productivity, which can affect IDH detection and interpretation. - Field conditions required deviations from ideal acquisition timing (e.g., non-apogee sun angles), potentially adding radiometric variability. - Pasture dynamics differ from customary swidden and can confound fragmentation–diversity relationships; controlled via covariates but residual effects may remain. - Scale dependence: results hinge on chosen sampling scale (70.15 ha), though sensitivity analyses were conducted.