Adaptation is thought critical for the survival of species under global change, but our understanding of human-induced evolution in the wild remains limited. Here we show that widespread deforestation has underpinned repeated color shifts in wild insects. Loss of forest cover has led to drastic color evolution across lineages that mimic the warning coloration of a toxic forest stonefly. Predation experiments indicate that the relative fitness of color phenotypes varies markedly between forested and deforested habitats. Genomic and coloration analyses of 1200 wild specimens reveal repeated selection at the ebony locus controlling color polymorphism across lineages. These findings represent a compelling case of rapid, replicated anthropogenic evolution linking both the agents and targets of natural selection, revealing the potential for wild populations to adapt in the wake of human-driven change.

Colour measurements

Digital photography was used to obtain color measurements as reflectance in the visible spectrum adequately captures variance in melanism (60), which underpins the mimetic resemblance between Zelandoperla and Austroperla. Dry-frozen stoneflies (866 Zelandoperla, 68 Austroperla) were photographed on an 18% neutral grey reflectance standard at a standardized distance and angle. Uniform lighting was attained using a fiber optic illuminator split through two guides diffused through a polystyrene foam cylinder 60mm in diameter. Photos were taken in RAW format using a Nikon D5200 with a Nikon AF-S DX 40mm micro lens under fixed ISO sensitivity (320), aperture (f/5) and shutter speed (1/8) using automatic focus. Raw image files were subsequently normalized to the neutral grey reflectance standard and linearized using the Multispectral Image Calibration and Analysis (MICA) toolbox  and DCRAW plugins in ImageJ, generating a multispectral image file.

Color measurements were extracted from the calibrated multispectral images using the RGB Measure function in ImageJ across three regions of interest (ROIs): (1) head, (2) pronotum, and (3) femora. The ROIs are representative of overall body melanism underpinning mimetic resemblance. The mean absolute normalized reflectance (0 – 100%) of the red, green, and blue (RGB) channels were sampled for each ROI using the circular selection tool in ImageJ, with a diameter of 30 pixels for the head and femora, and 40 for the pronotum. RGB values for each ROI were subsequently converted to achromatic greyscale values or “brightness”, ranging from 0% (black) to 100% (white) which was quantified as the arithmetic mean of the three channels. ROI brightness variables were standardized (centered on zero with a unit variance) and we performed a principal component analysis (PCA) using R version 4.0.5 (63) to obtain a univariate measure of overall body melanism. Overall, 73.94% of the total variance in brightness was attributed to the first principal component, with approximately equal loadings for each ROI (z₁ = 0.5813head + 0.5968pronotum + 0.5531femora). We also included the Leith Zelandoperla population from Foster et al. 2023. RAW image files (338 Zelandoperla, 22 Austroperla) were processed using the methodology outlined above (PCA: z₁ = 0.5899head + 0.5747pronotum + 0.5672femora; 82.44%). 

Predation experiments

We constructed 410 Zelandoperla color-models from odorless, non-toxic polymer modelling clay (Du-kit). The bodies of melanic models were constructed entirely from black clay, whereas non-melanic models were constructed from a mixture of two parts dark-brown to one part yellow clay. Each model (body length ~25 mm) comprised distinct ‘head’, ‘thoracic’ and ‘abdominal’ segments shaped and assembled by hand. Model legs were shaped by hand from black or brown solid core hookup cable (0.75 mm; Duratech, Electus Distribution; WH3032), and paired antennae/cerci from the plastic bristles of a Raven Hearth Brush. The distinctive white (wing-base) and yellow (tibiae) ‘warning’ color-features of Zelandoperla  were added to models with a fine paint brush, using white and yellow acrylic paint, respectively. In summary, black and brown models were constructed from identical materials , and identical proportions, differing only in the color of modelling clay used for the body and cable used for the legs.

Using not-toxic glue, models were attached to cleaned, flat river cobble stones (typical dimensions ~80 mm x ~50 mm) of either schist or greywacke, both of which are common in rivers of southeastern New Zealand. The black and brown color models for each site were attached to a random assortment of stone types, but with each model pair on matching stones (fig. S5). Glued specimens were left to dry in the laboratory for 48 hours to ensure that any glue odors were absent during subsequent field experiments. Unique codes were added to the underside of each stone using a permanent marker pen. Once attached to stones, models were photographed and transported to and from the field unstacked in the base layer of storage boxes.

Both Zelandoperla and Austroperla stoneflies have lengthy adult emergence seasons extending at least from early spring to late autumn. To ensure that predation experiments were ecologically relevant, we conducted them at four of our Zelandoperla sampling sites, during mid-Autumn 2024 (mid April to early May), at the tail end of the stonefly emergence season. We selected two forested (Leith: elevation 100-200 m; Whisky:elevation 200-220 m) and two deforested (Pigroot: elevation 500-520 m; Black Jacks: elevation 70-100 m) streams. At each site, approximately 50 pairs of Zelandoperla color models mounted on cobble stones were placed on riparian rocks/boulders at intervals of ~5-10 m, in total spanning several hundreds of metres of riparian habitat. For each model pair, black and brown models were placed 30-50 cm apart on a large stone/boulder. Models were left at their field locations for 72 hours prior to their recollection and returned to the laboratory for inspection and photography. Longer-running field trials were considered unsuitable given the potential risks of model damage associated with precipitation and fluctuating river levels linked to New Zealand’s oceanic climate.

Models were inspected for damage and photographed, with distinctive damage attributed to candidate predator groups (birds, rodents, or insects) based on the approach of Zvereva and Kozlov 2021 and Guerra et al 2024. We recorded attack data for each model, and tested for associations between attack rates and mimicry/coloration across all predator categories, sites, and habitat types, using Fisher’s exact test.

Statistical analysis and code

Phylogenetic logistic regression (phylogenetic generalized linear mixed model) was fit using the Bayesian software JAGS and executed from R via the R interface rjags, where parameters were estimated using Markov Chain Monte Carlo (MCMC) sampling with the Gibbs sampler. The number of individuals classified as mimics from population  was modelled as a binomial distribution with total number of individuals  and frequency of mimics.