\documentclass[10pt]{article}

\usepackage{fullpage}
\usepackage{setspace}
\usepackage{parskip}
\usepackage{titlesec}
\usepackage[section]{placeins}
\usepackage{xcolor}
\usepackage{breakcites}
\usepackage{lineno}
\usepackage{hyphenat}





\PassOptionsToPackage{hyphens}{url}
\usepackage[colorlinks = true,
            linkcolor = blue,
            urlcolor  = blue,
            citecolor = blue,
            anchorcolor = blue]{hyperref}
\usepackage{etoolbox}
\makeatletter
\patchcmd\@combinedblfloats{\box\@outputbox}{\unvbox\@outputbox}{}{%
  \errmessage{\noexpand\@combinedblfloats could not be patched}%
}%
\makeatother


\usepackage{natbib}




\renewenvironment{abstract}
  {{\bfseries\noindent{\abstractname}\par\nobreak}\footnotesize}
  {\bigskip}

\titlespacing{\section}{0pt}{*3}{*1}
\titlespacing{\subsection}{0pt}{*2}{*0.5}
\titlespacing{\subsubsection}{0pt}{*1.5}{0pt}


\usepackage{authblk}


\usepackage{graphicx}
\usepackage[space]{grffile}
\usepackage{latexsym}
\usepackage{textcomp}
\usepackage{longtable}
\usepackage{tabulary}
\usepackage{booktabs,array,multirow}
\usepackage{amsfonts,amsmath,amssymb}
\providecommand\citet{\cite}
\providecommand\citep{\cite}
\providecommand\citealt{\cite}
% You can conditionalize code for latexml or normal latex using this.
\newif\iflatexml\latexmlfalse
\providecommand{\tightlist}{\setlength{\itemsep}{0pt}\setlength{\parskip}{0pt}}%

\AtBeginDocument{\DeclareGraphicsExtensions{.pdf,.PDF,.eps,.EPS,.png,.PNG,.tif,.TIF,.jpg,.JPG,.jpeg,.JPEG}}

\usepackage[utf8]{inputenc}
\usepackage[T2A]{fontenc}
\usepackage[ngerman,polish,greek,english]{babel}



\usepackage{float}








\begin{document}

\title{Climate variability and drought modulate the role of structural refuges
for arthropods: a global experiment}



\author[1]{Gustavo Romero}%
\author[2]{Thiago Gonçalves-Souza}%
\author[3]{Tomas Roslin}%
\author[4]{Robert Marquis}%
\author[5]{Nicholas Marino}%
\author[6]{Vojtech Novotny}%
\author[7]{Tatiana Cornelissen}%
\author[8]{Jerome Orivel}%
\author[9]{Shen Sui}%
\author[10]{Gustavo Aires}%
\author[11]{Reuber Antoniazzi}%
\author[12]{Wesley Dattilo}%
\author[2]{Crasso Breviglieri}%
\author[13]{Annika Busse}%
\author[14]{Heloise Gibb}%
\author[15]{Thiago Izzo}%
\author[16]{Tomas Kadlec}%
\author[17]{Victoria Kemp}%
\author[18]{Monica Kersch-Becker}%
\author[19]{Michal Knapp}%
\author[20]{Pavel Kratina}%
\author[21]{Rebecca Luke}%
\author[22]{Stefan  Majnarić }%
\author[23]{Robin  Maritz}%
\author[24]{Paulo Martins}%
\author[25]{Esayas Mendesil}%
\author[26]{Jaroslav Michalko}%
\author[27]{Anna Mrazova}%
\author[28]{Samuel Novais}%
\author[7]{Cássio Pereira}%
\author[22]{Mirela Perić}%
\author[29]{Jana Petermann}%
\author[30]{Servio Ribeiro}%
\author[31]{Katerina Sam}%
\author[32]{M. Kurtis Trzcinski}%
\author[2]{Camila Vieira}%
\author[33]{Natalie Westwood}%
\author[34]{Maria Bernaschini}%
\author[35]{Valentina Carvajal}%
\author[36]{Ezequiel González}%
\author[37]{Mariana  Jausoro}%
\author[9]{Stanis Kaensin}%
\author[35]{Fabiola Ospina}%
\author[38]{Jacob Cristóbal-Pérez}%
\author[39]{Mauricio Quesada}%
\author[32]{Pierre Rogy}%
\author[40]{Diane Srivastava}%
\author[41]{Scarlett  Szpryngiel}%
\author[42]{Ayco Tack}%
\author[43]{Tiit Teder}%
\author[34]{Martin Videla}%
\author[43]{Mari-Liis  Viljur}%
\author[44]{Julia Koricheva}%
\affil[1]{State University of Campinas (UNICAMP)}%
\affil[2]{UNICAMP}%
\affil[3]{Swedish University of Agricultural Sciences}%
\affil[4]{University of Missouri-St. Louis}%
\affil[5]{Universidade Federal do Rio de Janeiro}%
\affil[6]{Institute of Entomology, Czech Academy of Sciences}%
\affil[7]{UFMG}%
\affil[8]{CNRS, UMR Ecofog}%
\affil[9]{New Guinea Binatang Research Center}%
\affil[10]{Federal Rural University of Pernambuco}%
\affil[11]{Stephen F Austin State University}%
\affil[12]{Instituto de Ecología AC}%
\affil[13]{Bavarian Forest National Park}%
\affil[14]{La Trobe University}%
\affil[15]{UFMT}%
\affil[16]{Czech University of Life Sciences Prague Faculty of Environmental Sciences}%
\affil[17]{Queen Mary, University of London}%
\affil[18]{The University of Alabama System}%
\affil[19]{Czech University of Life Sciences Prague}%
\affil[20]{Queen Mary University of London}%
\affil[21]{Royal Holloway University of London}%
\affil[22]{University of Zagreb}%
\affil[23]{University of the Western Cape}%
\affil[24]{Universidade Federal Rural de Pernambuco}%
\affil[25]{Jimma University}%
\affil[26]{Slovak University of Agriculture in Nitra}%
\affil[27]{Biology Centre Czech Academy of Sciences}%
\affil[28]{Universidad Nacional Autónoma de México Escuela Nacional de Estudios Superiores Unidad Morelia}%
\affil[29]{Universitat Salzburg}%
\affil[30]{Universidade Federal de Ouro Preto}%
\affil[31]{University of South Bohemia and Institute of Entomology}%
\affil[32]{The University of British Columbia}%
\affil[33]{The University of British Columbia Faculty of Science}%
\affil[34]{Instituto Multidisciplinario de Biología Vegetal}%
\affil[35]{Universidad de Caldas}%
\affil[36]{University of Koblenz Landau Faculty 5 of Education}%
\affil[37]{Universidad Nacional de Chilecito}%
\affil[38]{Instituto de Investigaciones en Ecosistemas y Sustentabilidad, Universidad Nacional Autónoma de México.}%
\affil[39]{Universidad Nacional Autonoma de Mexico}%
\affil[40]{University of British Columbia}%
\affil[41]{Swedish Museum of Natural History}%
\affil[42]{Stockholm University}%
\affil[43]{University of Tartu}%
\affil[44]{University of London}%


\vspace{-1em}



  
  \date{May 08, 2024}


\begingroup
\let\center\flushleft
\let\endcenter\endflushleft
\maketitle
\endgroup





\selectlanguage{english}
\begin{abstract}
Current climate change is disrupting biotic interactions and eroding
biodiversity worldwide. However, species sensitive to drought, high
temperatures and climate variability might persist in microclimatic
refuges, such as leaf shelters built by arthropods. We conducted a
distributed experiment across an 11,790 km latitudinal gradient to
explore how the importance of leaf shelters for terrestrial arthropods
changes with latitude, elevation and underlying climate. Our analyses
revealed leaf shelters to be key facilitative elements for the diversity
of arthropods. Predator diversity and overall biomass within shelters
increased with local drought and temperature variability, regardless of
latitude and elevation. In contrast, shelter usage by herbivores
increased with abundance of predators on those same plants and in wetter
climates. Projected increase in climatic variability and drought in
certain geographic regions is therefore likely to enhance the importance
of biotic refuges, especially for predators, in mitigating the impact of
climate change on species persistence.%
\end{abstract}%



\sloppy


\textbf{INTRODUCTION}

Climate warming, extreme precipitation events and temperature
variability are all intensifying at unprecedented rates (IPCC Climate
Change, 2014; Donat \emph{et al.} , 2016; Bathiany \emph{et al} .,
2018). Such changes cause species--environment mismatches, and have
already led to shifts in the geographic ranges and elevational or
vertical distributions of species as they track their climatic niche
(V\selectlanguage{ngerman}ázquez\emph{et al} ., 2017; Essl \emph{et al} ., 2019; Trisos \emph{et
al} ., 2020). As other means to avoid adverse conditions, terrestrial
species may seek microclimatic refuges (Scheffers \emph{et al} ., 2014;
Suggitt\emph{et al} ., 2018) and even manipulate their own microclimate
(Pincebourde and Casas, 2019). Microclimatic refuges provided by a
variety of habitat structural components -- including plant
architectural complexity, tree bark, downed woody debris, leaf litter,
and leaf shelters -- reduce the exposure of their inhabitants to climate
variability and extremes (Scheffers \emph{et al} ., 2014;
Suggitt\emph{et al} ., 2018; Pincebourde and Casas, 2019; Pinsky
\emph{et al} ., 2019). Such refuges can be particularly important to
ectotherms (Pinsky\emph{et al} ., 2019), such as arthropods. Ectotherms
depend on external energy to thermoregulate and are amongst the taxa
most threatened by global change (García-Robledo \emph{et al} ., 2016;
Warren \emph{et al.,}2018; Wagner, 2020; Wagner \emph{et al} ., 2021).
The availability of these refuges, and their impact on insect use, may
be one factor mitigating the various negative climatic factors leading
to insect decline.

Protection against harsh environmental conditions is typically provided
by ecosystem engineers, which are organisms that modify their habitats
by creating physical structures that can be used as structural refuges
by other organisms (Jones et al. Hastings \emph{et al} ., 2007;
Romero\emph{et al} ., 2015). In terrestrial ecosystems, bark beetles and
stem borers make stem and trunk galleries, gallers manipulate plant
physiology to produce galls, and various arthropods (e.g., caterpillars,
aphids, mites and thrips) build leaf shelters, such as leaf rolls and
ties, thus providing refuge to many other plant-dwelling organisms
(Wang\emph{et al.,} 2012; Vieira and Romero, 2013; Cornelissen \emph{et
al} ., 2015; Priest \emph{et al} ., 2021). The role of leaf-rolling
engineers as biodiversity amplifiers can be even stronger in dry
seasons, and extend to the whole plant, thus influencing arthropod
assemblages at temporal and microspatial scales (Vieira and Romero,
2013). After such shelters are abandoned by their creators, they become
available for other arthropods for a long period of time (Vieira and
Romero, 2013). Similar to many other biotic interactions, which are
stronger at lower latitudes (Schemske \emph{et al.,} 2009), the
facilitative effects of vertebrate and invertebrate ecosystem engineers
seem particularly pronounced in the wet tropics, but also in arid
regions (Romero \emph{et al} ., 2015). Therefore, these latitudinal,
regional (Romero \emph{et al} ., 2015) and temporal patterns (Vieira and
Romero, 2013) suggest that climatic conditions may be a common driver
underlying such facilitative interactions.

Large-scale patterns in climatic signatures have been recently
investigated using standardized replicated experiments, which can
differentiate the direct effects of climate from the indirect effects of
latitude or elevation on biotic interactions (LaManna \emph{et al} .,
2017; Roslin \emph{et al} ., 2017; Romero \emph{et al} ., 2018). Most of
these macroecological studies have focused on antagonistic interactions,
such as predation and competition (LaManna \emph{et al} ., 2017;
Roslin\emph{et al} ., 2017; Romero \emph{et al} ., 2018). Although the
importance of facilitative interactions among organisms is expected to
increase under stressful conditions (He \emph{et al} ., 2013;
Romero\emph{et al} ., 2015), no empirical studies to date have examined
the impacts of current climate and climate variability, or predicted the
impacts of future climates, on facilitative interactions at global
scale.

Here, we test the effects of current and future climatic scenarios on
facilitation mechanisms provided by leaf-rolling ecosystem engineers. We
compare how current and projected climate affects the usage of leaf
rolls by terrestrial arthropods, and how shelter usage is affected by
body size and trophic position. We hypothesize that abundance, species
richness and biomass of arthropods within leaf rolls will increase with
increasing climate variability and drought. We also expect that species'
responses to climate will vary with their trophic position and body
size, both of which are known to influence the metabolic requirements of
ectotherms (Daufresne \emph{et al.,} 2009) and their capacity to
dissipate heat (Rubalcaba \emph{et al.,} 2019). Competition at higher
trophic levels (predators) selects for larger body size. Larger-bodied
organisms typically require more favourable climatic conditions than
their prey (Petchey \emph{et al.,} 1999; Voigt \emph{et al.,} 2003;
Brose \emph{et al.,} 2012). Therefore, we expect larger predators to be
dominant in refuges especially in arid and climatically variable regions
(Fig. 1). Finally, with continuing climate change, we expect current
climatic dependencies to reflect future changes in geographic patterns
of refuge use by arthropods.

To test the above predictions, we conducted a global, coordinated
experiment at 52 sites across an 11,790 km latitudinal gradient (from
45.8° S to 60.2° N, Fig. 1) and elevation spanning from 5 to 2900 m
above sea level. To measure the beneficial effects of leaf rolls, we
recorded abundance and diversity of arthropods colonizing manually
rolled and control (unrolled) leaves, and calculated the standardized
mean difference (Hedges' d) between arthropod abundance, species
richness, biomass and body size in rolled and control leaves for each
site. To test how climate influences shelter colonization, we
investigated the effects of different categories of moderators
(predictor variables) on the Hedges' d effect size, including both
geographic distance (absolute latitude, elevation) and direct measures
of climate (temperature and precipitation). Finally, we provide
geographical interpretations of climate change scenarios on the patterns
of refuge use and predict future changes in refuge use by projecting
effect sizes to future (2070) climatic conditions at the study
locations.

\textbf{MATERIAL AND METHODS}

\textbf{Study sites}

We performed a geographically-coordinated experiment simulating the
presence of ecosystem-engineering structures built by leaf-rolling
caterpillars on plants. The experiment was replicated in 2017-2019 in 52
sites across an 11,790 km latitudinal gradient (45.8° S to 60.2° N) and
elevation spanning from 5 to 2900 m above sea level, comprising diverse
types of biomes, including tropical and temperate forests, xeric
shrubland (caatingas), savannas etc. We considered our sites spatially
independent, as 94\% of them were at least 500 km apart from each other.
A detailed description of the sites, including environment and plant
features, is provided in Table S1.

\textbf{Global experiment and arthropod samples}

Experiments at all 52 sites were conducted following a standardized
protocol. Each experiment followed a randomized block design; for each
site, we randomly chose 10 to 20 paired trees, shrubs or herbs (blocks).
Paired plants were at least 2 m apart from each other. Each block was at
least 6 m apart from each other. In most sites, the experiment was
conducted using a single native plant species, typically the most common
species locally. However, in some tropical forests, due to the high
species diversity and low relative density, we used different plant
species among blocks, but the plants within each pair always belonged to
the same species. We only used broadleaf native plant species that did
not bear any apparent type of indirect defence (e.g., domatia,
extra-floral nectaries, glandular trichomes) (Romero and Koricheva,
2011).

One plant in each pair was randomly chosen as a control whereas the
second plant was used for the leaf-rolling (shelter addition) treatment.
On each plant, we selected 5-10 fully expanded leaves without apparent
damage. Prior to the experiment, any arthropods present were removed by
hand from leaves. Then, the leaves on the treatment plant were manually
rolled from the adaxial to the abaxial surface transversally to the leaf
axis to form a cylinder of roughly 0.6 cm in diameter (Fig. S1) to mimic
shelters built by caterpillars from at least 17 moth and butterfly
families, including Hesperiidae, Nymphalidae, Gelechiidae, Oecophoridae,
Lasiocampidae, Pyralidae, Gracillariidae, Tortricidae, Geometridae,
Erebidae (Fitzgerald \emph{et al} ., 1991; Fukui, 2001; Lill and
Marquis, 2007). The rolls were secured with a metal hairpin (see Fig. S1
and Vieira and Romero, 2013). The control plant of the pair had 5-10
unrolled leaves marked with a metal hairpin. Rolled and control leaves
were exposed for 10 days in the field. This was deemed sufficient as
previous bioassays showed that leaf shelters can be colonized very
quickly (within 24 hours) (Vieira and Romero, 2013). Maximum width of
both control and rolled leaves was measured as an estimate of leaf size.

After 10 days of the experiment, we collected rolled and control leaves,
and stored them grouped per replicate and treatment in small zip-lock
plastic bags. The leaves were either frozen for later sorting or
immediately sorted to collect the invertebrates. We collected all the
invertebrates visible to the naked eye (except mites) and stored them in
ethanol. We identified the invertebrates to the lowest taxonomic level
possible and classified them into morphospecies and feeding guild (i.e.,
predator, parasitoid, herbivore, detritivore, omnivore). Individual body
size (dry body mass) was estimated from the dry mass (dried at
70\textsuperscript{o}C for 24h) or by measuring total body length and
then calculating the dry mass through published taxon-specific
allometric equations (Hódar, 1996). Four dependent variables were used
for analyses: arthropod abundance, species richness, biomass and mean
body size. Arthropod biomass represented the sum of all individual body
masses. All variables were weighted by the number of sampled leaves per
plant.

\textbf{Climate, topography and productivity data}

We extracted different sets of current climate data for two purposes:
(i) to test the prediction that climate variability and drought would
influence the use of shelters by arthropods at the experiment scale, and
(ii) to predict the effects of future climate variability on refuge
usage (Fig. 1). Climate data and topographic data were extracted from
WorldClim (Fick and Hijmans, 2017) version 2 (http://www.worldclim.org/)
and ENVIREM (Title and Bemmels, 2017) (https://envirem.github.io/). For
each site we chose four variables for temperature (bio1, bio2, bio4,
bio7), four variables for precipitation (bio12, bio14, bio15, Aridity),
two variables related to topography (TRI and topoWET) and one variable
denoting site productivity (AnnualPET). A detailed description of the
variables is presented in the Table S2. WorldClim variables were
extracted at 30 arc-second, 2.5 arc-minute and 10 arc-minute
resolutions. Since the above variables were very strongly correlated
among these different resolutions (Pearson correlation, r \selectlanguage{english}[?] 0.97), we
decided to use only data on a 10 arc-minute resolution. The variables
TRI and topoWET from ENVIREM were available only for 30 arc-second
resolutions. The variables Aridity and AnnualPET were extracted at 30
arc-second and 10 arc-minute resolutions; since these variables were
strongly correlated between the two resolutions (Pearson correlation, r
\selectlanguage{english}[?] 0.98), we used only data on 30 arc-second resolutions.

Future climate database at local scale was extracted from WorldClim
version 2 using MIROC5 (RCP8.5) and CCSM4 (RCP8.5) as representative
concentration pathways of CO\textsubscript{2} emission projected for
2070 (Romero \emph{et al} ., 2018). Future bioclimatic variables are not
available in ENVIREM databases. Since the bioclimatic variables were
very strongly correlated between these two predictive climatic models
(Pearson correlation, r \selectlanguage{english}[?] 0.98), we focused our analyses only on MIROC5
(RCP8.5) database.

We measured average local temperature during the 10-day experiment. For
sites missing such data, we extracted the mean near-surface daily
temperature from the RNCEP database (Kemp \emph{et al.,} 2012). This was
done for all the 10 days of the experiment, from which data we then
calculated the average over the experiment. The RNCEP database has a
spatial resolution of 2.5 x 2.5\textsuperscript{o} and a temporal
resolution of 6 h (Kemp \emph{et al.,} 2012).

\textbf{Statistical analyses}

\emph{Effect size calculations}

Even though we have used a standardized protocol, the sampling was
conducted by several researchers in different biomes. Thus, we used a
meta-analysis (meta-regression) approach to control for biases in
sampling (Gurevitch, 2013). The magnitude of the leaf shelter effect on
arthropod community structure (measured as abundance, species richness,
biomass and body size) was calculated at the site level (n = 52 sites)
using Hedges' d (standardized mean difference; Rosenberg \emph{et al} .,
2013) following the formula:

\(Hedges^{\prime}\ d\ =\ \frac{Xe-Xc}{\text{SDpooled}}\ \),

where \emph{Xe} and \emph{Xc} are the mean community measures for plants
that had their leaves rolled and control plants, respectively,
and\emph{SDpooled} is the pooled standard deviation of the two groups.
Positive effects indicate that leaf rolls had higher arthropod species
abundance, richness, biomass or mean body size relative to the control
leaves. Values around 0.2, 0.5 and 0.8 are considered small, medium and
large effect sizes, respectively (Rosenberg \emph{et al.,} 2013). Larger
effect size values mean larger differences in occupation of rolled
compared to control leaves. Although we initially classified the
surveyed arthropods into different feeding guilds (i.e., predators,
herbivores, detritivores, omnivores and parasitoids), the sample size
was only sufficiently large for separate analyses of predators and
herbivores.

\emph{Moderators (predictors)}

We investigated different categories of moderators, including absolute
latitude and elevation, as well as climate, productivity (annualPET),
topography (TRI, topoWET), leaf size and average temperature during the
experiments. Additional local moderators included available potential
prey (abundance of herbivores) for predators and abundance of predators
on experimental leaves for herbivores.

To create the climate gradients, we projected separately precipitation
and temperature variables in multivariate space using principal
component analysis (PCA) (Romero \emph{et al} ., 2018). Scores of the
first axis of PCA for temperature variables
(PC1\textsubscript{temperature}, which explained 78\% of the total
variance) represented a gradient spanning from higher temperature
variability to warming and lower temperature variability. Conversely the
first axis of PCA for precipitation variables
(PC1\textsubscript{precipitation}: 64\% of the total variance)
represented a rainfall (dry to wet) gradient, at global (Fig. 1) and
local scales (Figs. S2-S3). PC1\textsubscript{precipitation} included
mean annual precipitation (bio12), precipitation of driest month
(bio14), index of the degree of water deficit below water need (aridity)
and precipitation seasonality (coefficient of variation) along the year
(bio15) (Table S1). The only variable denoting rainfall variability
(bio15) did not contribute much to the PCA1 variance (Fig. 1, Fig. S2).
Therefore, the PC1\textsubscript{precipitation} represented a rainfall
gradient, which varied from dry (negative values) to wet conditions
(positive values) (Fig. 1, Fig. S2). Sites defined as being climatically
more suitable were those characterized by higher
PC1\textsubscript{precipitation} scores.
PC1\textsubscript{temperature}included mean annual temperature (bio1),
mean diurnal range (bio2), temperature seasonality (bio4) and
temperature annual range (bio7). The signs of original
PC1\textsubscript{temperature} scores were multiplied by -1 to change
their graphic presentation for the ease of interpretation, matching with
precipitation trends, thus producing a gradient from temperature
variability (negative values) to warming and lower temperature
variability (positive values). Sites defined as being thermally more
variable were those characterized by higher intra-annual (e.g., bio4,
bio7) or diurnal (e.g., bio2) variation in temperature (Fig. 1, Fig.
S3).

Some moderators were collinear (Fig. S4), and were thus removed from the
models based on their variance inflation factor (VIF) (Zuur \emph{et
al.,} 2010). The removed moderators were Latitude, topoWET and
AnnualPET. All remaining moderators had VIF \textless{} 2. The
best-fitting linear models were then obtained through backward
selection; the best models retained were those to present the smallest
Akaike Information Criteria corrected for small sample sizes, AICc
(\selectlanguage{english}[?]AICc\textgreater{}2). Null models (intercept-only effects, i.e., no
moderators included) were contrasted with models containing moderators
plus the intercept.

\emph{Meta-analysis and structural equation modelling}

We tested the effects of the moderators on Hedges' d effect sizes using
meta-regression approach (linear models), conducted using
\emph{metafor}package (version 2.1-0) (Viechtbauer, 2010) implemented in
R environment (R Development Core Team, 2019). Prior to the analyses, we
inspected outliers through diagnostic plots (Cook's distance). We also
visually inspected normality of the residuals using normal
quantile--quantile (Q--Q) plots; our models presented appropriate
behaviour in terms of residuals and outliers. Overall effect size was
obtained through random-effects models, and the moderators were examined
through mixed-effects models (Viechtbauer, 2010). Such mixed-effect
models were obtained using a restricted maximum likelihood (REML)
estimator, assuming a random variation among replicates within a
treatment (control or experimental), and fixed variation among
treatments. The moderator effects were estimated using Cochran's Q-test
with model fitting through\emph{mods} argument to ignore the intercept
and \emph{rma} function (Viechtbauer, 2010). Effect sizes were
considered statistically significant if their 95\% confidence intervals
did not overlap with 0.

Piecewise structural equation models (\emph{pSEM} ) were used to test
direct and indirect effects of latitude and underlying climatic
components (PC1\textsubscript{temperature} and
PCI\textsubscript{precipitation}) (Romero \emph{et al} ., 2018) on the
Hedges' d effect size for predators and herbivores; for this, we used
the \emph{piecewiseSEM} package in R (Lefcheck, 2016), and evaluated
model fits using Shipley's test of d-separation through Fisher's C
statistics (Lefcheck, 2016).

\emph{Raw data analysis and visualisation}

Analyses of raw data were performed to graphically visualize the
patterns of refuge use by predators and herbivores. For this, we built
linear mixed effects models (LMEs) and evaluated the effects of
treatments (control vs rolled leaves) in raw data using
\emph{nlme}package in R (Pinheiro \emph{et al} ., 2020). Treatment was
the fixed effect, and block was the random effect. We then used these
LME models to generate the predictive regression trends graphically
using\emph{visreg} (Breheny and Burchett, 2017) and \emph{ggplot2}
(Wickham, 2016) packages in R environment.

\emph{Predicting spatial and temporal variations: maps}

We mapped the effect of future climate change on the magnitude of the
effect of ecosystem engineers on arthropods using the effect sizes
(Hedges' d) obtained at the experiment scale.

The magnitudes of the engineering effects on arthropod communities
observed at each site and their association with climatic variables (see
Table S1 for variable definitions) were tested using linear models
(Gaussian error distribution). We used the function \emph{predict.lm} to
obtain the predicted relationship between independent variables and the
Hedges' d effect size (abundance, richness and biomass) (n = 52).
Temperature (bio1, bio2, bio4, bio7) and precipitation (bio12, bio14,
bio15) were included as independent continuous variables in separate
models (see Table S1 for variable definitions). Separate climatic
variables rather than the principal component scores identified in our
meta-analyses were used because the correlations structure among
variables may shift with climate change, and thus parameter estimation
may be biased in different models (Hadi and Ling, 1998). In addition,
using the climatic variables in the predictive modelling allowed us to
inspect the reliability of the \emph{predict.lm} function. However, for
both temperature and precipitation climatic models most variables were
collinear (see Figs. S2 and S3). We therefore first eliminated the
influence of collinearity using the variance inflation factor (VIF)
(Zuur \emph{et al.,} 2010). We then reduced the model through backward
selection, retaining only the significant variables and those that
improved the model fit (i.e., higher R\textsuperscript{2} values). We
modelled the current climate scenario and those projected to 2070 using
rcp8.5 (MIROC5) as a general circulation model. We checked normality,
homoscedasticity and outliers through graphical inspections (e.g.,
qq-plots, Cook's d, and influence), and residuals of the response
variables were log transformed when necessary. Previous studies have
used this approach to modelling distribution maps to provide a
geographical interpretation of global ecological patterns
(Tallavaara\emph{et al} ., 2018, Gusmao \emph{et al} ., 2020).

Model performance was tested using a \emph{k} -fold cross-validation
approach (see Appendix S1).

\textbf{RESULTS}

\textbf{Current interactions}

Overall, predators, herbivores and the pooled arthropod assemblages
including all feeding guilds (predators, herbivores, detritivores,
parasitoids, omnivores) were strongly influenced by the presence of leaf
rolls on plants (Fig. 2). Arthropod abundance, species richness, biomass
and mean body size were significantly higher in rolled compared to
control leaves (i.e., positive Hedges' d values, with 95\% confidence
intervals not spanning 0).

The effects of leaf rolls on pooled arthropod assemblages were not
significantly influenced by any moderator (Table 1). The effects of leaf
rolls on predators were mediated solely by climate, but the effects of
precipitation and temperature varied depending on predator metrics
(abundance, richness, biomass). The magnitude of the effect of leaf
rolls on predator abundance and species richness increased with drought
(PC1\textsubscript{precipitation}, Table 1, Fig. 3). Predator abundance
on plants was strongly positively correlated with predator species
richness (r = 0.86, t=12.2, P\textless{}0.001), and thus the effects of
leaf shelters on both community metrics were similarly influenced by
precipitation gradients (Table 1, Fig. 3). In contrast, the positive
effects of leaf rolls on predator biomass increased with greater
temperature variability, but decreased in warmer regions (Table 1, Fig.
3). Latitude did not directly influence refuge usage by arthropods, but
had a positive indirect effect via temperature gradient on the increase
in predator biomass in leaf shelters (PC1\textsubscript{temperature}, \selectlanguage{greek}β
\selectlanguage{english}= 0.28, Fig. S5c). Predator biomass within leaf shelters and on control
leaves was positively correlated with predator body size (shelter: r =
0.82, t=9.6, P\textless{}0.001; control: r = 0.94, t =16.6,
P\textless{}0.001), but not with predator abundance (shelter: r = 0.23,
t = 1.6, P = 0.11; control: r = 0.14, t=0.98, P=0.33). Thus, shelter
usage, mostly by larger predators, increased toward higher latitudes,
where temperature variability is higher.

In contrast to the effect of leaf rolls on predators, positive effects
of leaf shelters on the magnitude of the effect on herbivore abundance
and biomass were stronger in wet sites, but also on plants with a higher
abundance of predators (Table 1, Fig. 3). Herbivore biomass within leaf
shelters and on control leaves was strongly positively correlated with
the mean body size of the individuals (shelter: r = 0.88, t=12.4,
P\textless{}0.001; control: r = 0.74, t=6.7, P\textless{}0.001), but not
with their abundance (shelter: r = 0.05, t=0.37, P=0.71; control: r =
0.13, t=0.88, P=0.38). Average herbivore body size was 41\% smaller than
average predator body size (Fig. S6).

Likewise, data visualization performed using linear mixed effects
modelling with raw data revealed contrasting patterns of leaf shelter
usage by herbivores and predators (Fig. S7). Predator abundance, species
richness and biomass increased within refuges in a drier climate and
with increasingly variable temperatures. In contrast, herbivores were
more likely to use shelter under wetter climates (Fig. S7). Other
moderators such as elevation, ecosystem productivity, topography and
local average temperature did not directly influence the effects of
shelters on arthropods (Table 1).

\textbf{Future predicted interactions}

The magnitude of the effect of refuge on species richness and abundance
of predators is expected to increase by up to 33\% and 38\%,
respectively, in areas where precipitation is predicted to decline
(e.g., French Guiana and Mexico). In contrast, with increasing
precipitation predicted for sites in Europe, USA, and Borneo, the
effects of refuge on predator species richness and abundance are
expected to decrease by up to 25\% and 34\%, respectively (Fig. 4a, Fig.
S8).

Future scenarios of the influence of temperature on refuge usage by
larger predators depended on the temperature component evaluated (i.e.,
temperature variability {[}bio7{]} versus mean annual temperature
{[}bio1{]}). Maps of effects for predator biomass revealed weaker
changes in response to temperature variability over time, compared to
responses to precipitation (Fig. 4b). Indeed, key bioclimatic variables
reflecting precipitation (e.g., bio14) are predicted to be more variable
over time than key temperature variables (e.g., bio7) (Appendix S2).
This weak effect is manifested as smaller differences in future versus
current effect sizes (ranging from a decrease of up to 2.3\% and an
increase of up to 6.7\%), likely stemming from the fact that temperature
variability itself will change little in the future. However, under a
scenario of increasing mean annual temperature, a general warming is
predicted to decrease (6-17\%) the effect of shelters on larger
predators, especially in Europe (Fig. S9).

\textbf{DISCUSSION}

Our results indicate that structural refuges provided by leaf rolls
increase the local abundance, diversity, biomass and mean body size of
arthropods from different trophic levels on plants, and they do so on
all four continents on which the experiment was conducted. However, the
pattern of refuge use by herbivores and predators in response to climate
differed: while predators used refuges more with increasing drought and
temperature variability, use of refuges by herbivores showed opposite
trends. This suggests that herbivores might avoid the use of refuges in
sites with increasing drought and temperature variability, where risk of
predation is higher (Wetzel \emph{et al.,} 2016). Such predation risk
might drive herbivores to suboptimal microhabitats in harsh conditions.
Herbivorous insects are often smaller than predatory arthropods (see
Fig. S6), and thus likely are less vulnerable to climatic adversities
(Petchey \emph{et al.,} 1999; Voigt \emph{et al.,} 2003; Brose \emph{et
al.,} 2012). In addition, in contrast to many predators, herbivores can
build their own refuges using plants (Pincebourde and Casas, 2019) and
access a more predictable source of water from food (green leaves).
Predators rely on more stochastic water sources (e.g., prey haemolymph),
which are dependent on prey capture rates. Conversely, larger herbivores
might seek refuges against vertebrate predators under more stable and
favourable climatic conditions (Tvardikova and Novotny, 2012), where
predation pressure on larger herbivores is often high (Romero \emph{et
al} ., 2018). It is therefore likely that herbivores use leaf shelters
more as a refuge against natural enemies (Baer and Marquis, 2020) than
for protection against harsh climatic conditions.

Climate variability and drought appear to be the key determinants of
refuge use, especially by predatory arthropods. However, different
metrics of predator communities (abundance, species richness, biomass)
were differentially influenced by components of precipitation and
temperature. It is likely that drought may shape species composition
within leaf shelters, as the refuges tend to be occupied by species
sensitive to desiccation. Indeed, predatory arthropod taxa have varying
susceptibility to desiccation (Edney, 2012); leaf shelters could provide
them with humidity via leaf transpiration, and offer stable
microclimatic refuges against long droughts. However, we found refuge
usage to be higher for larger predators. As larger predators show a
smaller surface-to-volume ratio than small ones (K\selectlanguage{ngerman}ühsel\emph{et al} .,
2017), they are likely less prone to desiccation. Despite being less
vulnerable to desiccation, larger predators can be competitively
dominant over smaller ones in using available refuges.

Our results also reveal that larger predators were most influenced by
temperature variability. This effect is likely mediated by large
predators being subject to the stronger physiological constraints
imposed by temperature and/or having the strongest competitive ability
under such climatic conditions. In fact, temperature variability imposes
greater risks to arthropods than warming (Vasseur \emph{et al} ., 2014;
Colinet \emph{et al.,} 2015), especially for species from higher trophic
levels (Stireman \emph{et al} ., 2005). Refuges might ameliorate
temperature oscillations and heat waves (Caillon \emph{et al.,} 2014),
which are harmful especially to larger organisms, such as predators
(Voigt \emph{et al.,} 2003; Brose \emph{et al.,} 2012; Horne \emph{et
al} ., 2017). It is known that heat transfer (or dissipation) is reduced
in larger organisms because of their smaller surface area-to-volume
ratios (Horne \emph{et al} ., 2017), causing a threat in more thermally
variable climates. Due to their higher competitive ability, large
predators in particular may thus rely on refuges against extreme
temperature oscillations, despite being less vulnerable to desiccation
(Kühsel \emph{et al.,} 2017).

The climatic variables that we considered consist of both annual means
and variances, thus characterizing the overall climatic conditions of
each site. We found that increasing variability in climate, and not
changes in the mean of the climate variable, was positively correlated
with predator response to shelter availability. If communities of
predators experience higher climatic variability during at least some
part of the year over time, then shelter use in general may be favoured
by selection. This hypothesis is supported by the fact that the climatic
variable spanning the duration of the experiment itself (average
temperature over 10 days) did not moderate the influence of leaf
shelters on predators. Predators (e.g., spiders) are well known to
actively select sheltered microhabitats in harsh conditions and
fire-prone environments (de Omena and Romero, 2010). Therefore, we
suggest that the behaviour of using leaf shelters might be an
evolutionary response to drought and climatic variability.

Arthropods comprise over two-thirds of all terrestrial species, and are
key elements in intricate food webs. They provide valuable ecosystem
services, such as biological control, nutrient cycling and pollination.
Because arthropods are ectothermic and closely dependent on external
energy, they are also among the most threatened organisms on Earth due
to ongoing global changes (van Klink \emph{et al} ., 2020; Wagner, 2020,
Wagner et al. 2021). With increase in climate variability and in
frequency of extreme droughts, terrestrial arthropod communities should
be more dependent on suitable microclimate (Pinsky \emph{et al} .,
2019), with body size and trophic position having a major influence on
microhabitat use. Understanding how arthropods deal with climate
oscillations and extremes is thus pivotal for improved conservation and
mitigation strategies. Efficient strategies should aim at identifying
and preserving local conditions, such as habitat structural components,
to buffer climatic adversities. Here, we present clear evidence that
refuges provided by ecosystem engineers may mitigate predicted increases
in climate variability and drought. As a direct projection from the
patterns observed, we predict that future changes in climatic conditions
could change the global pattern of refuge use by arthropods. In the
future, ecosystem engineers are expected to increase in importance,
especially for predators and in sites where rainfall is expected to
decrease. Increased future usage of shelters by predators might then
strongly alter predation risk on herbivores, and even strengthen
intraguild predation pressure upon plants.

Our results point to larger predators as the dominant organisms using
refuges. This dominance likely results from higher competitive dominance
and/or higher sensitivity to global changes. The refuges can mitigate
climatic variation, allowing predators to adjust their use of niche
space according to their physiological tolerances. The ability to forage
from such refuges may still have a strong effect on their overall
performance, and on the overall imprint of climate on arthropod
community change. Increasing climatic variability and drought may
decrease the ability of such predators to feed and thus to perform
important ecosystem services (Barton and Schmitz, 2009; Rosenblatt and
Schmitz, 2016) such as biological control. Likewise, increasing refuge
usage by predators with climate change could displace herbivores to
suboptimal microclimates (Barton and Schmitz, 2009). Asymmetric reliance
on refuge use by predator and prey, and mismatches in predator-prey
encounter rates, could cause asymmetric extinctions of larger herbivores
and predators (Clark \emph{et al.,} 2020). While the realism of these
scenarios remains to be seen, current patterns of climatic impacts on
refuge use and their extensions into the future flag a potential
ecosystem transition. Importantly, they identify a minuscule scale in
microclimatic research as an important avenue for future research on
ecosystem functioning and community ecology.

\textbf{ACKNOWLEDGMENTS}

We acknowledge financial support for research provided by the Sa\selectlanguage{english}o Paulo
Research Foundation (FAPESP: grants 2017/09052-4, 2018/12225-0, and BPE
Proc. 2016/01209-9) and CNPq-Brazil through research grants to GQR. GQR
and PK gratefully acknowledge funding from the Royal Society, Newton
Advanced Fellowship (grant no. NAF/R2/180791). TR was supported by the
Academy of Finland (grant 322266) and a senior career support grant from
the Vice Chancellor of the Swedish University of Agricultural Sciences.
RJM thanks the Missouri Department of Natural Resources for access to
Cuivre River State Park, Zhaid Khaja for help in the field, and the
Students and Teachers as Research Scientists (STARS) program for
financial support. NACM acknowledges CAPES for a post-doctoral research
grant (PNPD/CAPES \#2013/0877), VN acknowledges the CSF 20-10205S and
ERC 669609 grants. TC acknowledges CNPq (307210/2016-2) and CAPES. JO
acknowledges access to the Paracou research station in French Guiana,
which is supported by CIRAD (https://paracou.cirad.fr/) and financial
support provided by an Investissement d'Avenir grant of the Agence
Nationale de la Recherche (CEBA, ANR-10-LABX-25-01) and by the PO-FEDER
2014-2020, R\selectlanguage{ngerman}égion Guyane (BiNG, GY0007194). WD and RA are grateful to
Rosamond Coates and all the staff of the Estación de Biología Tropical
Los Tuxtlas (Universidad Nacional Autónoma de México), besides the staff
of the Santuario del Bosque de Niebla and Centro de Investigaciones
costeras La Mancha (Instituto de Ecología A.C.) who have made the
fieldwork possible. HG was supported by an Australian Research Council
Future Fellowship (FT130100821). TJI is supported by CNPQ (grants
309552/2018-4). TK was supported by the Czech University of Life
Sciences Prague (grant no. 42900/1312/3166). VK acknowledges funding
from grant NE/K016148/1/Natural Environment Research Council. MK was
supported by the Czech University of Life Sciences Prague (grant no.
42900/1312/3166). RM acknowledges Tswalu Kalahari for permission to
conduct fieldwork and Bryan Maritz and Graham Alexander for help in the
field. PMM acknowledges the Dunedin Botanic Garden's administration for
allowing data collection, and thank Penelope Jacks and Esther Dale for
their valuable support during different stages of the fieldwork. EM
acknowledges Tamiru Shimales. JM was supported by the Slovak national
grant VEGA2/0100/18. AM acknowledges the GAJU n. 04-048/2019/P grant.
CCP acknowledges CAPES (financial code 001). MSP and SM thank the Nature
Park Medvednica and Ministry of Environment and Energy of the Republic
of Croatia for allowing access to the Croatian study

sites. SPR thank CNPq (306572/2019-2). MKT acknowledges Danielle
Courcelles for assistance in field work. CV was supported by Coordenação
de Aperfeiçoamento de Pessoal de Nível Superior (CAPES) Proc.
88882.314749/2019-01. EG was supported by the Czech University of Life
Sciences Prague (grant no. 42900/1312/3166). MQ was supported by grants
from Universidad Nacional Autónoma de México (PAPIIT IV200418),
SADER-CONACYT 291333, CONACYT-UNAM to LANASE (2019-LN299033); MQ thank
Gumersindo Sanchez for help in the field. TT was supported by the
Estonian Research Council (grant no. PRG741), and by the Internal Grant
Agency of the Faculty of Environmental Sciences, Czech University of
Life Sciences Prague (grant no. 42900/1312/3166). M-LV was supported by
the Estonian Research Council (grant no. PRG741).

\textbf{REFERENCES}

Baer, C.S. \& Marquis, R.J. (2020). Between predators and parasitoids:
complex interactions among shelter traits, predation and parasitism in a
shelter\selectlanguage{english}-building caterpillar community. \emph{Functional Ecology} , 34,
2186-2198.

Barton, B.T \& Schmitz, O.J. (2009). Experimental warming transforms
multiple predator effects in a grassland food web. \emph{Ecology
Letters} , 12, 1317-1325.

Bathiany, S., Dakos, V., Scheffer, M. \&Lenton, T.M. (2018). Climate
models predict increasing temperature variability in poor
countries.\emph{Science Advances} , 4, eaar5809.

Breheny, P. \& and Burchett, W. (2017). Visualization of regression
models using visreg. \emph{The R Journal} , 9, 56-71.

Brose, U., Dunne, J.A., Montoya, J.M., Petchey, O.L., Schneider, F.D. \&
Jacob, U. (2012). Climate change in size-structured
ecosystems.\emph{Philosophical Transactions of the Royal Society B} ,
367, 2903---2912.

Caillon, R., Suppo, C., Casas, J., Woods, H.A. \& Pincebourde, S.
(2014). Warming decreases thermal heterogeneity of leaf surfaces:
implications for behavioural thermoregulation by
arthropods.\emph{Functional Ecology} , 28, 1449-1458.

Clark, J.S., Scher, C.L. \& Swift, M. (2020). The emergent interactions
that govern biodiversity change. \emph{PNAS} , 117, 17074-17083.

Colinet, H., Sinclair, B.J., Vernon, P. \& Renault, D. (2015). Insects
in fluctuating thermal environments. \emph{Annual Review of Entomology}
, 60, 123-140.

Cornelissen, T., Cintra, F. \& Santos, J.C. (2015). Shelter-building
insects and their role as ecosystem engineers. \emph{Neotrop. Entomol}
., 45, 1-12.

Daufresne, M., Lengfellner, K. \& Sommer, U. (2009). Global warming
benefits the small in aquatic ecosystems. \emph{Proc. Natl Acad. Sci.
USA} , 106, 12788--12793.

de Omena, P.M. \& Romero, G.Q. (2010). Using visual cues of microhabitat
traits to find home: the case study of a bromeliad-living jumping spider
(Salticidae). \emph{Behavioral Ecology} , 21, 690--695.

Donat, M. G., Lowry, A. L., Alexander, L. V., O'Gorman, P. A. \& Maher,
N. (2016). More extreme precipitation in the world's dry and wet
regions. \emph{Nat. Clim. Change} , 6, 508--513.

Edney, E.B. (2012). Water balance in land arthropods. Springer-Verlag,
Berlin.

Essl, F., Dullinger, S., Genovesi, P. \ldots{} Bacher, S. (2019). A
Conceptual Framework for range-expanding species that track
human-induced environmental change. \emph{BioScience} , 69, 908--919.

Fick, S. E. \& Hijmans, R. J. (2017). WorldClim 2: new 1-km spatial
resolution climate surfaces for global land areas. \emph{Int. J.
Climatol} ., 37, 4302--4315.

Fitzgerald, T.D., Clark, K.L., Vanderpool, R. \& Phillips, C. (1991).
Leaf shelter-building caterpillars harness forces generated by axial
retraction of stretched and wetted silk. \emph{Journal of Insect
Behavior} , 4, 21--32.

Fukui, A. (2001). Indirect interactions mediated by leaf shelters in
animal--plant communities. \emph{Population Ecology} , 43, 31-40.

Garcia-Robledo, C., Erin K. Kuprewicz, E.,K., Staines, C.L., Erwin, T.L.
\& Kress, W.J. (2016). Limited tolerance by insects to high temperatures
across tropical elevational gradients and the implications of global
warming for extinction. \emph{PNAS} , 113, 680-685.

Gurevitch, J. (2013). Meta-analysis of results from multisite studies.
In Koricheva, J., Gurevitch, J. \& Mengersen, K. (eds.). Handbook of
meta-analysis in ecology and evolution. Princeton University Press,
Princeton and Oxford.

Gusmao, R.A.F., Hernandes, F.A., Vancine, M.H., Naka, L.N., Dona, J.,
Goncalves\selectlanguage{english}-Souza, T. (2020). Host diversity outperforms climate as a
global driver of symbiont diversity in the bird\selectlanguage{english}-feather mite system.
Diversity and Distributions https://doi.org/10.1111/ddi.13201

Hadi, A.S. \& Ling, R.F. (1998). Some cautionary notes on the use of
Principal Components Regression. \emph{The American Statistician} , 52,
15-19.

Hastings, A., Byers, J. E., Crooks, J. A., Cuddington, K., Jones, C. G.,
Lambrinos, J. G., Talley, T. S. \& Wilson, W.G. (2007). Ecosystem
engineers in space and time. \emph{Ecology Letters} , 10, 153-164.

He, Q., Bertness, M.D, Altieri, A.H. (2013). Global shifts towards
positive species interactions with increasing environmental
stress.\emph{Ecology Letters} , 16, 695-706.

Hodar, J.A. (1996). The use of regression equations for estimation of
arthropod biomass in ecological studies. \emph{Acta Oecologica} , 17,
421-433.

Horne, C.R., Hirst, A.G. \& Atkinson, D. (2017). Seasonal body size
reductions with warming covary with major body size gradients in
arthropod species. \emph{Proc. R. Soc. B.} , 284, 20170238.

IPCC Climate Change 2014: Synthesis Report (eds Core Writing Team,
Pachauri, R. K. \& Meyer, L. A.).

Jones, C.G., Lawton, J.H., Shachak, M. (1997). Positive and negative
effects of organisms as physical ecosystem engineers. \emph{Ecology,}78,
1946-1957.

Kemp, M.U., van Loon, E.E., Shamoun\selectlanguage{english}-Baranes, J. \& Bouten, W. (2012).
RNCEP: global weather and climate data at your fingertips. \emph{Methods
in Ecology and Evolution} , 3, 65-70. (2012).

Kuhsel, S., Bruckner, A., Schmelzle, S., Heethoff, M. \& Bluthgen, N.
(2017). Surface area--volume ratios in insects. \emph{Insect Science} ,
24, 829-841.

LaManna, J. A. et al. (2017). Plant diversity increases with the
strength of negative density dependence at the global
scale.\emph{Science} , 356, 1389--1392.

Lefcheck, J. S. (2016). piecewiseSEM: piecewise structural equation
modelling in R for ecology, evolution, and systematics. \emph{Methods
Ecol. Evol} ., 7, 573--579.

Lill, J.T., Marquis, R.J. (2007). Microhabitat manipulation: ecosystem
engineering by shelter-building insects. In: K. M. D. Cuddington, J. E.
Byers, A. Hastings, and W. G. Wilson (eds.), Ecosystem engineers:
concepts, theory, and applications in ecology, pp. 107- 138, Elsevier
Press. San Diego, CA.

Petchey, O.L., McPhearson, P.T., Casey, T.M. \& Morin, P.J. (1999).
Environmental warming alters food-web structure and ecosystem
function.\emph{Nature} , 402, 69--72.

Pincebourde, S. \& Casas, J. (2019). Narrow safety margin in the
phyllosphere during thermal extremes. \emph{PNAS} , 116, 5588-5596.

Pinheiro J., Bates D., DebRoy S., Sarkar D., R Core Team. (2020). nlme:
Linear and Nonlinear Mixed Effects Models. R package version 3.1-148.

Pinsky, M.L., Eikeset, A.M., McCauley, D.J., Payne, J.L. \& Sunday, J.M.
(2019). Greater vulnerability to warming of marine versus terrestrial
ectotherms. \emph{Nature} , 569, 108--111.

Priest, G.V., Cameroata, F., Powell, S., Vasconcelos, H.L., Marquis,
R.J. (2021). Ecosystem engineering in the arboreal realm: Heterogeneity
of wood-boring beetle cavities and their use by cavity-nesting ants.
Oecologia (in press).

R Development Core Team (2019). R: A Language and Environment for
Statistical Computing (R Foundation for Statistical Computing).

Romero, G.Q. \& Koricheva, J. (2011). Contrasting cascade effects of
carnivores on plant fitness: a meta-analysis. \emph{J. Anim. Ecol} .,
80, 696-704.

Romero, G.Q., Goncalves-Souza, T., Kratina, P., Marino, N.A.C, Petry,
W.K., Sobral-Souza, T. \& Roslin, T. (2018). Global predation pressure
redistribution under future climate change. \emph{Nature Climate Change}
, 8, 1087--1091.

Romero, G.Q., Goncalves-Souza, T., Vieira, C., Koricheva, J. (2015).
Ecosystem engineering effects on species diversity across ecosystems: a
meta-analysis. \emph{Biological Reviews} , 90, 877-890.

Rosenberg, M.S., Rothstein, H.R. \& Gurevitch, J. (2013). Effect sizes:
conventional choices and calculations. In Koricheva, J., Gurevitch, J.
\& Mengersen, K. (eds.). Handbook of meta-analysis in ecology and
evolution. Princeton University Press, Princeton and Oxford.

Rosenblatt, A. E. \& Schmitz, O. J. (2016). Climate change, nutrition,
and bottom-up and top-down food web processes. \emph{Trends Ecol. Evol}
., 31, 965--975.

Roslin, T. et al. (2017). Higher predation risk for insect prey at low
latitudes and elevations. \emph{Science} , 356, 742--744.

Rubalcaba, J.G., Gouveia, S.F., Olalla\selectlanguage{english}-Tarraga, M.A. (2019). A
mechanistic model to scale up biophysical processes into geographical
size gradients in ectotherms. \emph{Global Ecology and Biogeography} ,
28, 793-803.

Scheffers, B.R. Edwards, D.P., Diesmos, A., Williams, A.E \& Evans, T.A.
(2014). Microhabitats reduce animal's exposure to climate
extremes.\emph{Global Change Biology} , 20, 495-503.

Schemske, D. W. et al. (2009). Is there a latitudinal gradient in the
importance of biotic interactions? \emph{Annu. Rev. Ecol. Evol. Syst} .,
40, 245--269.

Stireman III, J.O., L.A. Dyer, D.H. Janzen, M.S. Singer, J.T. Lill, R.J.
Marquis, R.E. Ricklefs, G.L. Gentry, W. Hallwachs, P.D. Coley, J.A.
Barone, H.F. Greeney, H. Connahs, P. Barbosa, H.C. Morais, and I.R.
Diniz. (2005). Climatic unpredictability and parasitism of caterpillars:
Implications of global warming. \emph{PNAS} , 102, 17384-17387.

Suggitt, A.J., Wilson, R.J., Isaac, N.J.B. et al. (2018). Extinction
risk from climate change is reduced by microclimatic
buffering.\emph{Nature Clim Change} , 8, 713--717.

Tallavaara, M., Eronen, J.T. \& Luoto, M. (2018). Productivity,
biodiversity, and pathogens influence the global hunter-gatherer
population density. \emph{Proc Natl Acad Sci USA} , 115, 1232-1237.

Title, P. O. \& Bemmels, J. B. (2017). ENVIREM: an expanded set of
bioclimatic and topographic variables increases flexibility and improves
performance of ecological niche modeling. \emph{Ecography} , 41,
291--307.

Trisos, C.H., Merow, C. \& Pigot, A.L. (2020). The projected timing of
abrupt ecological disruption from climate change. \emph{Nature} , 580,
496-501.

Tvardikova, K. \& Novotny, V. (2012). Predation on exposed and
leaf-rolling artificial caterpillars in tropical forests of Papua New
Guinea. \emph{Journal of Tropical Ecology} , 28, 331-341.

van Klink, R., Bowler, D.E., Gongalsky, K.B., \ldots{} Chase, J.M.
(2020). Meta-analysis reveals declines in terrestrial but increases in
freshwater insect abundances. \emph{Science} , 368, 417-420.

Vasseur, D. A. et al. (2014). Increased temperature variation poses a
greater risk to species than climate warming. \emph{Proc. R. Soc. B} ,
281, 20132612.

Vazquez, D.P., Gianoli, E., Morris, W.F., \& Bozinovic, F. (2017).
Ecological and evolutionary impacts of changing climatic
variability.\emph{Biological Reviews} , 92, 22-42.

Viechtbauer, W. (2010). Conducting meta-analyses in R with the metafor
package. \emph{Journal of Statistical Software} , 36, 1-48.

Vieira, C. \& Romero, G.Q. (2013). Ecosystem engineers on plants:
indirect facilitation of arthropod communities by leaf\selectlanguage{english}-rollers at
different scales. \emph{Ecology} , 94, 1510-1518.

Voigt W., Perner J., Davis A.J., Eggers T., Schumacher, J., Bahrmann R.,
Fabian B., Heinrich W., Kohler G., Lichter D., Marstaller R. \& Sander
F.W. (2003). Trophic levels are differentially sensitive to
climate.\emph{Ecology} , 84, 2444-2453.

Wagner, D.L. (2020). Insect declines in the Anthropocene. \emph{Annual
Review of Entomology} , 65, 457-480.

Wagner, D.L., Grames, E.M., Forister, M.L., Berenbaum, M.R., Stopak, D.
(2021). Insect decline in the Anthropocene: Death by a thousand
cuts.\emph{PNAS} 118 e2023989118.

Wang, H.G. Marquis, R.J., Baer, C.S. (2012). Both host plant and
ecosystem engineer identity influence leaf\selectlanguage{english}-tie impacts on the arthropod
community of Quercus. \emph{Ecology} , 93, 2186-2197.

Warren, R., Price, J., Graham, E., Forstenhaeusler, N. \& VanDerWal, J.
(2018). The projected effect on insects, vertebrates, and plants of
limiting global warming to 1.5\selectlanguage{english} degC rather than 2\selectlanguage{english} degC.~\emph{Science}
,~360, 791--795.

Wetzel, W.C., Screen, R.M., Li, I. \ldots{} \& Yang, L.H. (2016).
Ecosystem engineering by a gall\selectlanguage{english}-forming wasp indirectly suppresses
diversity and density of herbivores on oak trees. \emph{Ecology} , 97,
427-438.

Wickham, H. (2016). ggplot2: elegant graphics for data analysis.
Springer-Verlag New York.

Zuur, A.F., Ieno, E.N., Elphick, C.S. (2010). A protocol for data
exploration to avoid common statistical problems. \emph{Methods in
Ecology and Evolution} , 1, 3--14.

\textbf{SUPPORTING INFORMATION}

Additional supporting information may be found online in the Supporting
Information section at the end of the article.

\selectlanguage{english}
\FloatBarrier
\end{document}

