Latitudinally structured variation in the temperature dependence of damselfly growth rates

Ecology Letters, (2013) 16: 64–71
LETTER
Viktor Nilsson-Örtman,1* Robby
Stoks,2 Marjan De Block,2 Helena
Johansson1,3 and Frank
Johansson1,4
doi: 10.1111/ele.12013
Latitudinally structured variation in the temperature
dependence of damselfly growth rates
Abstract
The Metabolic Theory of Ecology predicts that the slope of the rate–temperature relationship, E, remains
consistent across traits and organisms, acting as a major determinant of large-scale ecological patterns.
Although E has recently been shown to vary systematically, we have a poor understanding of its ecological
significance. To address this question, we conducted a common-garden experiment involving six damselfly
species differing in distribution, estimating E at the level of full-sib families. Each species was sampled
throughout its latitudinal range, allowing us to characterise variation in E along a latitudinal gradient spanning 3600 km. We show that E differs among populations and increases with latitude. E was right-skewness across species, but this was largely an artefact of the latitudinal trend. Increased seasonality towards
higher latitude may contribute to the latitudinal trend in E. We conclude that E should be seen as a trait
involved in local adaptation.
Keywords
Environmental variation, growth rate, metabolic theory of ecology, thermal dependence, universal temperature dependence.
Ecology Letters (2013) 16: 64–71
INTRODUCTION
Temperature is one of the main forces acting on all organisms,
influencing the rates of most biological processes, from enzyme
kinetics to species interactions and population growth rates
(Angilletta 2009; Daniel & Danson 2010). Biological rate processes
generally increase with temperature, reaching a single thermal
optimum, Topt, and decrease rapidly at temperatures above the maximum. The shape of the initial, rising phase was recently used by
Gillooly et al. (2001) as the basis for the Metabolic Theory of
Ecology (MTE). They proposed that the rate at which a wide range
of traits increases with temperature during this phase follows simple
Boltzmann–Arrhenius kinetics:
R ¼ b0 M b e E=kT ;
ð1Þ
where R is the rate, b0 is a an organism-dependent scaling factor
describing the metabolic efficiency per unit mass, M is body mass, b
is a universal allometric constant (assumed to be 3/4), E is the activation energy, k is Boltzmann’s constant and T is the temperature
in Kelvin (Brown et al. 2004). This model is assumed to hold over
the physiologically relevant temperature range (PTR) below Topt
(Savage et al. 2004). The temperature-response term E thus represents a measure of phenotypic plasticity, defined within the PTR.
Using metabolic rate data from a wide range of species, Gillooly
et al. (2001) argued that the E remains largely stable across species,
initially suggesting that it ranges between 0.2 and 1.2 eV, but later
refining their estimate, suggesting that E is tightly constrained
between 0.6 and 0.7 eV, with an average value of 0.65 eV (Allen &
Gillooly 2006; Allen et al. 2006). This final prediction is known as
1
Department of Ecology and Environmental Science, Umeå University,
the universal temperature dependence (UTD; Gillooly et al. 2001).
The MTE relies critically on the validity of this prediction, as it postulates that ecological interactions across all levels of organisation,
from predator–prey interactions to patterns of biodiversity at the
global scale, are a consequence of the UTD (Allen & Gillooly 2006;
Allen et al. 2006).
Several recent studies have challenged the MTE’s reliance on the
UTD on both theoretical and empirical grounds (Clarke 2004;
Glazier 2005). From the perspective of the MTE, the temperature
dependence of biological traits is a static, non-evolving relationship
with most of the variation in absolute rates being explained by variation in body size and the normalisation constant, b0 (Gillooly et al.
2001; Allen & Gillooly 2007). The MTE shares this focus on how
evolution cannot easily overcome biochemical constraints with the
‘hotter-is-better’ hypothesis (Frazier et al. 2006; Angilletta Jr et al.
2010). An emerging, contrasting view is to consider a trait’s temperature response as a dynamic, evolving relationship, with E as one
of several potentially useful ways of describing it. A number of
recent reviews of published temperature responses have revealed
systematic structure in the variation in E that can be explained by
differences in trait motivation [e.g. traits that promote survival vs.
traits involved in food acquisition (Dell et al. 2011; Englund et al.
2011)], latitude (Irlich et al. 2009) and the level of organisation at
which traits are expressed (Dell et al. 2011), strongly suggesting that
the amount of temperature-induced plasticity that a trait display,
described by the parameter E, is of ecological importance to
individual organisms.
A major source of controversy stems from the fact that E can be
measured at multiple hierarchical levels. Firstly, E can be estimated
3
Centre of Excellence in Biological Interactions, University of Helsinki, PO Box
SE-90187, Umeå, Sweden
65, 00014, Helsinki, Finland
2
4
Laboratory of Aquatic Ecology and Evolutionary Biology, University of
Leuven, Ch. Deberiotstraat 32, BE-3000, Leuven, Belgium
Department of Ecology and Genetics, Uppsala University, SE-75236, Uppsala,
Sweden
*Correspondence: E-mail: [email protected]
© 2012 Blackwell Publishing Ltd/CNRS
Letter
among species: a trait is then measured at each species’ optimal temperature and data from many species are plotted together to reveal a
shared pattern (i.e. an interspecific approach). Secondly, E can be
estimated within species: a trait is then measured in many individuals
at different temperatures and used to estimate a species-level
response (i.e. an intraspecific approach). Interspecific studies have
typically found relatively little variation in E and tended to support a
UTD close to 0.65 eV (Brown et al. 2004; Allen & Gillooly 2007; but
see Terblanche et al. 2007), whereas intraspecific studies have repeatedly detected significantly more variation in E than predicted by the
MTE (Irlich et al. 2009; Dell et al. 2011; Englund et al. 2011). To
some extent, this variation may reflect limited precision when estimating E, for example due to the choice of temperature range (Knies
& Kingsolver 2010), how biological rates are measured (Bennett
1990) and differences in whether organisms can acclimate to experimental conditions (Fedorenko 1975; Huey & Berrigan 1996). Regardless, these studies suggest a need to study variation in E over
multiple hierarchical levels, both above and below the species level
(Chown 2001; Enquist et al. 2003; Terblanche et al. 2007).
To provide a direct test of these contrasting views of E and to
test the hypothesis that E is involved in local adaptation, we investigated within-species variation in E in a clade of European damselflies. Although earlier studies have considered repeated estimates of
E from the same species as pseudoreplicates (Frazier et al. 2006;
Dell et al. 2011), we explicitly focused on variation in E at the level
of genotypic groups. We will hereafter refer to variation in E
among latitudinal populations of the same species as interpopulation
variation in E. The main novelty of this study lies in the focus on
variation at the level of genotypes, the standardised sampling
scheme and large geographical scope: each species is sampled from
populations in the north, centre and south of each species’ latitudinal range. If E is strongly related to individuals’ fitness and evolves
in response to environmental variation, we expect to find evidence
of local adaptation with respect to E at this level.
METHODS
We here study six species of damselflies (Odonata: Coenagrionidae)
in the genus Coenagrion. The distribution patterns of the six species
clearly differ, occurring in northern, central and southern Europe
(Fig. 1). Together, the ranges of these species cover more than 30°
of latitude and single species cover up to 16°. Larvae are aquatic
predators on smaller invertebrates. The northern species require
2 years to finish larval development, whereas central species are univoltine (Norling 1984; Corbet 1999). The life cycles of southern
species are poorly known, but they are presumably univoltine in
most areas.
Common-garden experiment
Each studied species was sampled at sites close to its northern and
southern range margins, and at the centre of its range (hereafter
termed latitudes and denoted N, S and C; Fig. 1). Three populations were sampled at each latitude to minimise the influence of
maternal and environmental effects. The total number of populations per species was typically 9. From each population, we reared
the offspring from between 2 and 12 field-collected mated females
(hereafter called ‘families’) at three or four of the following temperatures: 16.3, 19.5, 21.5 or 24.0 °C. Although we cannot fully
Thermal dependence of damselfly growth rates 65
exclude the possibility that some of the offspring of the same
female represent paternal half-sibs, we will refer to these as fullsib families as damselfly males typically remove all of the sperm
from previous matings (Miller & Miller 1981; Waage 1986). The
temperatures used were chosen to reflect typical field water
temperatures in most areas (Nilsson-Örtman et al. 2012). The fieldwork and rearing was divided across two laboratories and 3 years:
C. armatum Charpentier (all populations) and C. puella Linnaeus
(C, N) were collected and reared during 2008 at 16.3, 19.5 and
21.5 °C; C. johanssoni Wallengren (all) was collected and reared during 2008 at 19.5, 21.5 and 24.0 °C; C. mercuriale Charpentier (all),
C. puella (S), C. pulchellum Vander Linden (C, N) and C. scitulum
Rambur (C, S) were collected and reared during 2009 at 16.3,
19.5, 21.5 and 24.0 °C; C. puella (additional N), C. pulchellum
(S, additional N) and C. scitulum (N) were collected and reared
during 2010 at 16.3, 19.5, 21.5 and 24.0 °C. We will refer to the
temperature range 16.3–21.5 °C as the ‘lower’ temperature range,
19.5–24.0 °C as the ‘upper’ temperature range and 16.3–24.0 °C as
the ‘full’ temperature range. Biases introduced by estimating E
based on different temperature ranges are treated in detail below.
Only two central C. johanssoni populations were sampled due to
rainy and windy conditions during fieldwork. A list of sampled
populations and the number of families per population is presented in the Supporting information, Table S1.
The rearing conditions and growth rate calculations have been
described in detail elsewhere (Nilsson-Örtman et al. 2012) and only
the most important details are summarised here. Mated females
or their eggs were brought to the laboratory in Umeå, Sweden
(C. armatum, C. pulchellum, C. puella and C. scitulum) or Leuven,
Belgium (C. johanssoni and C. mercuriale). Immediately following
hatching, 20 full-sib offspring from each female were transferred
to individual 100-mL rearing containers, and five individuals were
placed at random in each of the three or four climate chambers.
The photoperiod was kept at a fixed 14:10 h light:dark (L:D) ratio.
On initiation, the experiment consisted of 5060 individuals from
304 families. Larvae were fed 6 days a week on laboratory-reared
brine shrimp. Each individual was photographed at 0, 42, 84 and
126 days after hatching (hereafter called ‘measurement events’). At
the same time, water was changed and the rearing containers were
cleaned. The maximum distance between the distal parts of the
eyes of each individual was measured from these photographs as
an approximation of overall size, as this trait display less allometric
variation during ontogeny than other measures of size (Corbet
1999). Relative growth rates (RGRs) were calculated by modelling
individual growth trajectories as third-degree polynomial functions
of head width, hw, and time, t, of the form hw(t) = at + bt2 + ct3,
with a, b and c being polynomial coefficients estimated from the
data. Size-corrected RGRs were then calculated as the slope of the
curve at the point in time, when 1/5 of the adult size of each
species (based on field-collected adults from central populations)
had been reached. This method accounts for differences in both
initial and final size (Rose et al. 2009).
Modelling temperature responses
To estimate the activation energy E of growth rates for each fullsib family, we fitted the least-squares regression of ln(RGR) and
1/kT (the inverse rearing temperature in Kelvin) to the data from
each family group. The slope of this regression provides an estimate
© 2012 Blackwell Publishing Ltd/CNRS
66 V. Nilsson-Örtman et al.
Letter
(a)
(b)
(c)
(d)
(e)
(f)
Figure 1 European distributions of the study species and interpopulation variation in the activation energy term E (eV). Shown are the locations of the sampled
populations (red squares) and family level temperature responses, E. Asterisks mark populations that differ significantly compared with other populations of the same
species. Note that the activation energy of growth rates increases with latitude for all species (the trend is significant in Coenagrion johanssoni, C. puella, C. mercuriale and
C. scitulum). Open circles are outliers. Maps reproduced from Askew (2004) with permission from Apollo Books.
of the coefficient E in the Arrhenius equation (eqn 1). The inverse
of the temperature term was used to produce positive estimates of
E.
To ensure that subsequent analyses only included families where
we could obtain a meaningful estimate of E, we used the following
criteria for which families to include: (1) a family must have at least
two RGR estimates from each of at least three temperatures and (2)
the small-sample AIC of the model containing the linear Arrhenius
coefficient must be lower than that of an intercept-only model.
Families that did not fulfil these criteria were excluded.
The Arrhenius model and the UTD prediction are based on the
assumption that the slope of the rate–temperature relationship
remains linear over the entire PTR. In reality, temperature responses
are often curved, with the slope typically decreasing towards Topt
(Knies & Kingsolver 2010). This is important here for two reasons:
the first is that if some species have lower Topt, estimates of E will
© 2012 Blackwell Publishing Ltd/CNRS
be biased downwards for these. The second is that this can introduce a bias as not all our estimates of E were based on measurements from the exact same temperatures (due both to mortality and
logistics). Estimates of E based on the upper range of temperatures
will likely be biased downwards relative to those estimated based on
the lower range of temperatures. These issues were dealt with in
two ways, described next.
First, we assessed whether each temperature response was linear
within the experimental temperature range by comparing the AIC
scores of the linear Arrhenius model against a model including the
second order polynomial of the relationship between ln(RGR) and
1/kT. Note that our power to detect significant curvature is
limited with only three or four temperatures. The direction of
curvature was calculated from the second order derivative of the
polynomial model (Knies & Kingsolver 2010). This analysis
indicated that the slope decreased towards higher temperatures in
Letter
many families. Six families displayed lower growth rates at 24.0
than 21.5 °C, indicating that their thermal optima were located
within the experimental temperature range. For these families, we
re-estimated E on the lower range of temperature (Supporting
Information, Fig S1).
Secondly, for those species where E was estimated based on different temperature ranges, we quantified the strength and direction
of this bias and corrected for it directly. We fitted species-specific
ANCOVAs with E as the response variable, latitude as a continuous
fixed effect and temperature range (lower, upper, all) as a fixed,
categorical effect. The fitted coefficients describing the effects of
temperature range on E were then used to adjust the family level
estimates of E. This analysis confirmed that the inclusion of the
highest temperature caused a downward bias in E, but this was corrected for through this approach. No correction was applied to data
from C. armatum, C. johanssoni and C. scitulum as the same temperature range was used to estimate E for all families of these species.
Thus, the distribution of E within these species was not affected by
this bias, but the mean value may have been affected. As C. johanssoni was reared at the upper range of temperatures, estimates of
E are likely conservative. Conversely, C. armatum was reared at the
lower range of temperatures, which may have resulted in an overestimation of E. Including species as a random effect when testing
for trends in E (see below) did not alter the outcome of subsequent
analyses.
The final data set (families with at least two RGR estimates at
three or more temperatures) consisted of 2924 individuals from 207
families. On average, there were approximately 3.9 individuals per
combination of family and temperature.
Climatic predictors of E
To test the local adaptation hypothesis, we related the observed
variation in E to measures of environmental conditions in each
sampling area. We extracted environmental data from the WorldClim database (http://www.worldclim.org) for the 2.5 arc minute
grid (c. 5 9 5 km) that contained the sampled sites. The environmental variables used were mean annual temperature (MAT), mean
annual precipitation (MAP), annual seasonality of temperature
(AST), annual seasonality of precipitation (ASP), mean temperature
of the warmest quarter (TWQ), mean precipitation of the warmest
quarter (PWQ) and mean diurnal range (MDR).
Statistical analyses
The average value of E for each species and species 9 latitude was
tested against the UTD prediction using two-sided one-sample
t-tests. Differences in E among latitudinal populations of the same
species were tested using two-sided two-sample t-tests. To test for
skewness in the distribution of E values, we used Lilliefors normality test, as implemented in the R package nortest (Gross 2006).
The relationship between E and the explanatory variables were
modelled using a linear mixed model approach (LME). We first
tested whether E displayed a significant latitudinal trend across
species, while also testing for the effects of possible confounding
factors. This model included latitude, latitude2 (to test for a nonlinear latitude trend), altitude and sample size (number of individuals of each family when fitting eqn 1) as continuous fixed factors.
Rearing year was included as a random effect to account for possi-
Thermal dependence of damselfly growth rates 67
ble differences in environmental and experimental conditions. After
fitting the full models, fixed model terms were removed sequentially
on the basis of AIC scores. The presented AIC scores were calculated from the final, minimal models. Following model simplification, we also tested for a significant interaction between latitude and
species. To test for latitudinal trends within each species, we used
one-way ANOVAs.
Next, as latitude is merely a composite variable of underlying
environmental variables, we ran a second set of models to test for
the independent effects of the environmental variables. As sample
size and altitude were non-significant in the latitude model, they
were not included in these models. Rearing year was included as a
random effect. This model initially included temperature and precipitation averaged over the whole year (MAT and MAP), the seasonality of these variables (AST and ASP), average temperature and
precipitation during the warmest quarter (TWQ and PWQ) and
mean diurnal range (MDR). All predictor variables were centred and
normalised to 1 SD prior to the analyses (Schielzeth 2010). Model
selection was carried out as for the latitude model.
All statistical analyses were performed in R 2.14.0 (R core team.
2008). The LME models were implemented using the function lmer
in the R package lme4 (Bates et al. 2011). Confidence intervals and
significance levels of fixed effects were calculated with a Markov
chain Monte Carlo (MCMC) sampling approach, using the function
pvals.fnc in the package LanguageR (Baayen 2008). The number of
MCMC iterations was set to 10 000. Significance levels of terms
that were not included in the reduced models were calculated by reintroducing these terms singly into the reduced models. The final
models and their parameter estimates were very similar when
excluding the random effect, including species as a random effect
or using a forward stepwise model selection approach. Significance
of the Species 9 Latitude interaction term was calculated using the
ANOVA function in the package car, but note that there is no general
agreement about how to calculate degrees of freedoms in mixedeffects models (Pinheiro & Bates 2000). All analyses were performed with E as the response variable and this value was calculated at the family level, constituting this experiment’s lowest
replication level.
RESULTS
Our data on interpopulation variation in E offers several novel
insights, which we summarise here and justify below. Firstly, E was
on average significantly higher than the UTD prediction of 0.65 eV
(t-test: d.f. = 193, t = 8.68, P < 0.001). Secondly, E displayed considerable interpopulation variation (Fig. 1). Specifically, E increased
with latitude both within and across species (Fig. 2; coefficient:
0.022, t = 9.94, P < 0.001). On the basis of these observations, we
reject the hypothesis that E is evolutionarily stable.
Family level activation energies of growth rates
Of the 207 family level responses that had at least two measurements of growth rates at a minimum of three temperatures, the linear Arrhenius model provided a meaningful description in 194 cases
(94%), as judged by AIC scores (Supporting Information, Table S3).
Those 13 families that were excluded generally had small sample
sizes [per-family sample size of excluded families was 10.4 ± 2.0
(95% CI) compared with 14.4 ± 0.7 for included families]. The
© 2012 Blackwell Publishing Ltd/CNRS
Table 1 Average values and skewness of activation energies (E) of growth rates.
Measurements are averaged over 8–10 populations from across the latitudinal
range of each species. The sample size is given in parentheses. (a) Mean activation energy ± 95% CI. (b) g1 skewness of E based on the original data. (c) g1
skewness of E after correcting for the latitudinal bias by subtracting each population’s mean value of E from the corresponding family level responses. P-values
are calculated from Lilliefors non-normality tests. Note that the data have been
corrected for the effect of using different temperature ranges to estimate E
Arm
Joh
Pue
Pul
Mer
Sci
1.0
1.5
Letter
0.5
Activation energy (eV)
2.0
68 V. Nilsson-Örtman et al.
40
45
50
55
60
65
Latitude (°N)
Figure 2 Latitudinal trend in the activation energy, E (eV), of growth rates.
Activation energy plotted against latitude for 194 family level temperature
responses. Solid line show the best-fitting linear regression (coefficient: 0.022,
t = 9.94, P < 0.001). Values on the y-axis have been jittered for clarity.
Arrhenius model predicts temperature responses below Topt to be
linear when plotted at the scale of 1/kT. Of the 194 responses
fulfilling the above criteria, 26 displayed significant downward curvature (concave downward) and 16 displayed significant upward
curvature (concave upward; Supporting information, Fig S1). Six of
the former also displayed reduced growth rates at the highest temperature (Supporting information, Fig S1) and the lower range of
temperature was instead used to estimate E for these. As five of
these families belonged to the southern species C. mercuriale, it
seems likely that the thermal optimum of this species is lower than
for the other species. Because of this, estimates of E may be biased
downwards for this species. Except for these families, however, the
nonlinear coefficients were small and only the results from the linear Arrhenius model (adjusted for the effect of the temperature
range used to estimate E) are analysed in detail below.
The mean value of E for growth rates across all species was
0.84 ± 0.04 eV (we report mean activation energy values ± 95% CI
throughout this article), significantly higher than the UTD
prediction of 0.65 eV. In four of the six species, activation energy
(averaged over their latitudinal range) was significantly higher than
0.65 eV, with the two southern species, C. mercuriale and C. scitulum
being the exceptions (Table 1).
Investigation of interpopulation variation in E revealed striking
patterns. The most noticeable is that E increases with latitude
(Fig. 2). The LME model revealed this trend to be significant and
linear (Table 2a) and that E increase with 0.022 eV for every degree
of latitude (95% CI: 0.017–0.027). There was also a significant
Species 9 Latitude interaction (d.f. = 5, v2 = 34.15, P < 0.001),
driven by a steeper latitudinal trend in C. johanssoni. When looking
at each species separately, all species displayed a positive correlation
between E and latitude (Fig. 1). Species-specific ANOVAs revealed
this relationship to be significant (P < 0.05) for four species
(C. johanssoni, C. puella, C. scitulum and C. mercuriale), but not significant for C. armatum (t = 1.595, P = 0.12) and C. pulchellum
(t = 1.077, P = 0.287).
© 2012 Blackwell Publishing Ltd/CNRS
Species
(a) Mean
E (eV)
All species (194)
Coenagrion armatum (29)
C. johanssoni (20)
C. puella (35)
C. pulchellum (50)
C. mercuriale (25)
C. scitulum (35)
0.84
1.04
1.30
0.73
0.88
0.65
0.62
±
±
±
±
±
±
±
0.04
0.10
0.19
0.08
0.04
0.08
0.05
(b) Raw data
(c) Biascorrected
Skew
P
Skew
P
1.10
0.36
0.15
0.63
0.17
0.29
0.08
<0.001
0.382
0.632
0.466
0.722
0.591
0.830
0.48
0.36
0.56
0.42
0.26
0.65
0.40
0.145
0.720
0.192
0.753
0.434
0.514
0.048
Table 2 Results from linear mixed model approach models on the effect of latitude and environmental variables on the slope of temperature responses, E. AIC
values were calculated after removing non-significant terms from the models
Predictor
(a)
Latitude
Latitude2
Altitude
N larvae
(b)
AST
TWQ
MAT
ASP
PWQ
MAP
MDR
Estimate
SE
t
P
AIC
0.2008
0.0000
0.0155
0.0297
0.0203
0.0003
0.0214
0.0222
9.8700
0.0280
0.7260
1.3400
0.0001
0.9936
0.4622
0.2762
4.74
0.7498
0.8912
1.3229
0.0018
0.0123
0.0449
0.0258
0.1881
0.2339
0.3796
0.0344
0.0216
0.0211
0.0274
3.9870
3.8100
3.4850
0.0510
0.5680
2.1260
0.9430
0.0001
0.0001
0.0004
0.9490
0.5870
0.0328
0.3686
7.26
AIC, Akaike Information Criterion value; AST, annual seasonality of temperature; TWQ, average temperature of the warmest quarter; MAT, mean annual
temperature; ASP, annual seasonality of precipitation; PWQ, average precipitation
of the warmest quarter; MAP, mean annual precipitation; MDR, mean diurnal
range. All predictor variables were centred and normalised to 1 SD.
One aspect of the variation in E that has been previously
reported (Dell et al. 2011) is a strong right-skewness in species-level
data. We observed the same pattern in our family level responses
when pooling data from all species (Fig. 3a). Analysed this way, the
distribution of E was significantly right-skewed (Lilliefors normality
test, D = 0.09, P < 0.001; g1 skewness = 1.10) and the overall median value of E was 0.81, 0.03 eV lower than the overall mean of
0.84 eV (Fig. 3A). When breaking down variation in E by species,
no species displayed significant skewness (Table 1b) and individual
values ranged between 0.17 and 0.63. As the majority of the
responses analysed originated from the southern half of the latitudinal range (the median latitudinal position was 50.5° N, whereas the
mean was 52.5° N), the observed latitudinal trend in E may have
contributed to the overall right-skewness of E. To test this hypothesis, we corrected for the latitudinal bias, normalising each family’s
value of E by subtracting the mean value at the latitude from where
Letter
Thermal dependence of damselfly growth rates 69
Environmental predictors of E
(a)
In the environmental model, only the three variables relating to
temperature, AST, TWQ and MAT, were retained in the final
model (Table 2b). These variables were strongly correlated (AST vs.
MAT: Pearson’s q = 0.89, P < 0.001; AST vs. TWQ: Pearson’s
q = 0.65, P < 0.001) so we also tested them singly. AST then
retained the positive correlation with E (coefficient = 0.14,
t = 6.19, P < 0.001; AIC = 19.54) and TWQ retained the negative
correlation
(coefficient = 0.13,
t = 5.55,
P < 0.001;
AIC = 24.63). MAT, however, then displayed a negative correlation
(coefficient = 0.16, t = 6.93, P < 0.001; AIC = 13.10). Neither
of the single-parameter models nor the full model with AST, TWQ
and MAT had lower AIC or BIC scores than the model containing
only latitude (Table 2).
DISCUSSION
(b)
Figure 3 Histogram showing the frequency distribution of family level
temperature responses before (a) and after (b) correcting for the latitudinal
sampling bias. The grey columns contains all responses and the semi-transparent,
coloured columns show the same data classified by latitudinal origin: red
represent families originating from below the mean latitude of 52.5° N, blue
represent families originating to the north of 52.5° N. Solid and dashed vertical
lines denote the mean and median value respectively. E was calculated by fitting
eqn 1 to data from full-sib families. Note that the data have been corrected for
the effect of using different temperature ranges to estimate E.
it originated. This resulted in a reduction of the skewness of E.
Although a tendency of right-skewness remained, this was no longer
significant (Lilliefors test, D = 0.07, P = 0.145; g1 skewness = 0.48;
Fig. 3b). C. scitulum displayed significant right-skewness after correcting for the latitudinal bias (Table 1c).
The largest significant differences in E between adjacent latitudinal populations were found for C. puella, where E ranged from
0.57 ± 0.08 to 0.77 ± 0.13 eV from southern to central populations
(separated by c. 1000 km) and C. scitulum, where E ranged from
0.56 ± 0.05 eV to 0.75 ± 0.08 from central to northern populations,
over a similar distance (Fig. 1). Note that the ANCOVAs used to correct for the effect of estimating E based on different temperature
ranges did not cause the observed differences among latitudinal
populations or skewness of E as it only affects the distribution of
E values within each latitudinal population.
Our study shows that E, the temperature-response term, can differ
strongly among populations from different parts of species’ ranges,
suggestive of local adaptation. Our results challenge the UTD prediction that underlies the MTE and lend support to the interpretation of E as a trait with ecological importance for individual
organisms.
We find that E is on average significantly higher across all species
(0.84 ± 0.04 eV) than the value of 0.65 eV predicted by the UTD,
and that populations of a single species separated by fewer than 9°
of latitude may differ by as much as 0.20 eV in average activation
energy. The most striking aspect of our results is that E increases
with latitude, both across and within species. We show that this latitudinal trend in E may bias inferences regarding E, contributing
both to the right-skewness that has been observed, here and in
other studies, for many traits (Dell et al. 2011), and to the average
value of E, forming the basis of the UTD prediction of 0.65 eV
(Allen et al. 2006).
The overall distribution of E in our data is significantly rightskewed (Fig. 3a), but this skewness largely arises as an artefact of
the latitudinal trend, as the majority of our responses originate from
the southern part of the latitudinal gradient. By correcting for this
bias, the distribution of values becomes approximately normal
(Fig. 3b). A striking observation, however, is that across- and
within-species values of skewness tends to converge towards the
average value of 0.48 after correcting for the latitudinal bias, with
the exception of C. pulchellum (Table 1). Although our results are
inconclusive, this could indicate that there exist some mechanisms
generating and maintaining a slight right-skewness also at the genotype level. As low E requires biological rates to be uniformly high
over a broader range of temperatures, generalist–specialist trade-offs
relating to thermal breadth (Angilletta et al. 2003) are one such factor that could constrain the evolution of very low values of E.
Those species and populations in our study that fell closest to the
value of 0.65 eV were sampled in the southern part of the latitudinal
gradient, between 35° N and 50° N latitude. This matches the latitudinal band (30° N–45° N) from where the majority of species-level
thermal responses reviewed by Irlich et al. (2009) originated. That
review, concerning activation energies of metabolic and development
rates in 129 insect species, also detected a weak but significant latitudinal trend, but did not find average levels of E to be significantly
different from 0.65 eV. The inclusion of more high- and mid-latitude
© 2012 Blackwell Publishing Ltd/CNRS
70 V. Nilsson-Örtman et al.
populations in our study most likely contributed to our finding that
E is significantly higher than predicted by the UTD.
Our interpretation of this is that the suggested UTD of 0.65 eV
may reflect not only an average over different traits and hierarchical
levels but also a geographical bias in the choice of species and populations studied. Although Glazier (2005) discussed this problem,
geographical biases have received little attention in the literature
with respect to temperature responses. This implies that if activation
energies in organisms distributed at random across the earth (or
reflecting biodiversity patterns) were measured, the predictions from
the UTD would likely be dramatically different. Considering the
uneven distribution of biodiversity (Gaston 2000), the consequences
of these geographical biases are unlikely to be trivial, and we caution against using average values of E for predicting global diversity
patterns (Allen et al. 2002) or responses to climate change (Dillon
et al. 2010) without carefully considering this potential bias.
Enquist et al. (2003) showed that the vertical elevation of temperature responses (parameter b0 in eqn 1) change across latitudes, possibly reflecting selection on increased growth potential (through
changes in metabolic efficiency or intensity) in areas with a short
growth season (Enquist et al. 2007). Our results strongly suggest
that also the slope of the rate–temperature relationship, E, can
evolve across broad environmental gradients. But is the observed
variation in E adaptive? The LME model with environmental variables suggested that variation in E likely reflects some aspect of the
thermal environment (Table 2). However, that E is higher in northern, thermally variable environments conflicts with the expectation
from the physiological literature, that organisms inhabiting thermally
variable environments should evolve enzymes with wider thermal
tolerances, equivalent to lower E (Somero & Low 1977). Empirical
support for this idea mainly comes from studies comparing tropical
and temperate taxa (e.g. Cunningham & Read 2002). That we find
the exact opposite latitudinal pattern here may be related to the fact
that even our southernmost populations are located at higher latitudes than the ‘high-latitude’ populations used in many other studies. To further explore this paradoxical result, we performed a
hypothetical ‘transplantation experiment’. Using simulated water
temperature data from high, low and mid-latitudes (Supporting
information, Figure S1), we calculated what the consequences would
be, in terms of realised annual growth rate, if we would introduce
an organism with a southern or central physiology into a high-latitude environment, and vice versa. This clearly showed that it would
rather be beneficial to have low E in cold and variable high-latitude
environments (Supporting information, Table S2). Instead, the
latitudinal trend within the temperate zone could possibly reflect a
phenological effect, that is, that it may be increasingly important for
damselfly larvae to exploit short periods of beneficial conditions in
more seasonal environments (Forrest & Miller-Rushing 2010). At
present, however, the selective forces responsible for this latitudinal
trend remain elusive.
The differences we observe in the northern and southern species
pairs are intriguing, where the studied species differ in the slope of
the latitudinal trends in E, despite having similar latitudinal distributions (Fig. 1). Differences in habitat preference could be important
in the case of the southern species C. mercuriale. In contrast to the
other species’ preference for shallow ponds, this habitat specialist
exploits small, slow-flowing streams that may provide a more constant temperature regime across latitudes. We are less confident,
however, that differences in the thermal characteristics of the
© 2012 Blackwell Publishing Ltd/CNRS
Letter
preferred habitat can explain differences between the two northern
species. Of these, C. johanssoni is restricted to more acidic, speciespoor, low-predation sites than C. armatum, raising the possibility that
differences in ecological factors such as resource availability, life history variation, phenology or seasonal variation in mortality, predation
and cannibalism (Norling 1984; Johansson 1993; Suhling & Lepkojus
2001; Enquist et al. 2003; Stoks et al. 2005) could be involved in
producing the steep relationship between E and latitude shown by
the former species. Although we cannot draw strong conclusions
regarding the causative mechanisms behinds the observed variation,
comparative data sets such as ours can guide future theoretical and
empirical explorations of how variation in temperature responses
depend upon the interaction between ecological and environmental
factors.
Our study is the first to characterise geographical variation in the
activation energy of a physiological trait below the species level and
for full-sibs. This approach revealed previously unrecognised patterns of variation – in particular a latitudinal trend in E across and
within species. These results challenge the idea that evolution cannot easily overcome thermodynamic constraints. We also highlight a
latitudinal bias in the literature that previous studies have not
addressed, which could contribute to the UTD prediction of
0.65 eV, as well as the right-skewness observed in species-level data.
Our findings suggest that temperature responses must be studied
over many different hierarchical and spatial scales if we are to gain
a deeper understanding of how organisms are affected by, and
evolve in response to, changes in temperature.
AUTHOR CONTRIBUTIONS
VNÖ, RS, MDB, HJ and FJ conceived the study and performed
the laboratory experiment; VNÖ developed the computer models,
analysed data and wrote the manuscript; VNÖ, RS, MDB, HJ and
FJ contributed to writing. MDB is a postdoctoral fellow of the
Fund for Scientific Research-Flanders (FWO). Financial support to
RS came from FWO and the KULeuven research fund, to VNÖ,
HJ and FJ from the Swedish Research Council Formas and to HJ
from Umeå University and the Academy of Finland.
REFERENCES
Allen, A.P. & Gillooly, J.F. (2006). Assessing latitudinal gradients in speciation
rates and biodiversity at the global scale. Ecol. Lett., 9, 947–954.
Allen, A.P. & Gillooly, J.F. (2007). The mechanistic basis of the metabolic
theory of ecology. Oikos, 116, 1073–1077.
Allen, A.P., Brown, J.H. & Gillooly, J.F. (2002). Global biodiversity, biochemical
kinetics, and the energetic-equivalence rule. Science, 297, 1545–1548.
Allen, A.P., Gillooly, J.F., Savage, V.M. & Brown, J.H. (2006). Kinetic effects of
temperature on rates of genetic divergence and speciation. Proc. Natl Acad. Sci.
USA, 103, 9130–9135.
Angilletta, M.J. (2009). Thermal Adaptation: A Theoretical and Empirical Synthesis.
Oxford University Press, New York.
Angilletta, M.J., Wilson, R.S., Navas, C.A. & James, R.S. (2003). Tradeoffs and
the evolution of thermal reaction norms. Trends Ecol. Evol., 18, 234–240.
Angilletta, M.J. Jr, Huey, R.B. & Frazier, M.R. (2010). Thermodynamic effects
on organismal performance: is hotter better? Physiol. Biochem. Zool., 83, 197–
206.
Askew, R.R. (2004). The dragonflies of Europe. Revised edition. Harley Books,
Colchester.
Baayen, R.H. (2008). Analyzing Linguistic Data: A Practical Introduction to Statistics
Using R. Cambridge University Press, Cambridge.
Letter
Bates, D., Maechler, M. & Bolker, B. (2011). lme4: Linear mixed-effects models
using S4 classes.
Bennett, A.F. (1990). Thermal dependence of locomotor capacity. Am. J. Physiol.
Regul. Integr. Comp. Physiol., 259, R253–R258.
Brown, J.H., Gillooly, J.F., Allen, A.P., Savage, V.M. & West, G.B. (2004).
Toward a metabolic theory of ecology. Ecology, 85, 1771–1789.
Chown, S.L. (2001). Physiological variation in insects: hierarchical levels and
implications. J. Insect Physiol., 47, 649–660.
Clarke, A. (2004). Is there a universal temperature dependence of metabolism?
Funct. Ecol., 18, 252–256.
Corbet, P.S. (1999). Dragonflies: Behaviour and Ecology of Odonata. Cornell University
Press, Ithaca, NY.
Cunningham, S. & Read, J. (2002). Comparison of temperate and tropical
rainforest tree species: photosynthetic responses to growth temperature.
Oecologia, 133, 112–119.
Daniel, R.M. & Danson, M.J. (2010). A new understanding of how temperature
affects the catalytic activity of enzymes. Trends Biochem. Sci., 35, 584–591.
Dell, A.I., Pawar, S. & Savage, V.M. (2011). Systematic variation in the
temperature dependence of physiological and ecological traits. Proc. Natl Acad.
Sci. USA, 108, 10591–10596.
Dillon, M.E., Wang, G. & Huey, R.B. (2010). Global metabolic impacts of
recent climate warming. Nature, 467, 704–706.
Englund, G., Öhlund, G., Hein, C.L. & Diehl, S. (2011). Temperature
dependence of the functional response. Ecol. Lett., 14, 914–921.
Enquist, B.J., Economo, E.P., Huxman, T.E., Allen, A.P., Ignace, D.D., Gillooly,
J.F. et al. (2003). Scaling metabolism from organisms to ecosystems. Nature,
423, 639–642.
Enquist, B.J., Kerkhoff, A.J., Huxman, T.E. & Economo, E.P. (2007). Adaptive
differences in plant physiology and ecosystem paradoxes: insights from
metabolic scaling theory. Glob. Change Biol., 13, 591–609.
Fedorenko, A.Y. (1975). Feeding characteristics and predation impact of
Chaoborus (Diptera, Chaoboridae) larvae in a small lake. Limnol. Oceanogr., 20,
250–258.
Forrest, J. & Miller-Rushing, A.J. (2010). Toward a synthetic understanding of
the role of phenology in ecology and evolution. Philos. Trans. R. Soc. Lond. B
Biol. Sci., 365, 3101–3112.
Frazier, M.R., Huey, R.B. & Berrigan, D. (2006). Thermodynamics constrains the
evolution of insect population growth rates: “warmer is better”. Am. Nat.,
168, 512–520.
Gaston, K.J. (2000). Global patterns in biodiversity. Nature, 405, 220–227.
Gillooly, J.F., Brown, J.H., West, G.B., Savage, V.M. & Charnov, E.L. (2001).
Effects of size and temperature on metabolic rate. Science, 293, 2248–2251.
Glazier, D.S. (2005). Beyond the “3/4-power law”: variation in the intra-and
interspecific scaling of metabolic rate in animals. Biol. Rev., 80, 611–662.
Gross, J. (2006). nortest: Tests for Normality. R package version 1.0.
Huey, R.B. & Berrigan, D. (1996). Testing evolutionary hypotheses of acclimation.
In: Animals and Temperature: Phenotypic and Evolutionary Adaptation (ed Johnston, I.
A. & Bennett, A.F.). Cambridge University Press, Cambridge, pp. 205–237.
Irlich, U.M., Terblanche, J.S., Blackburn, T.M. & Chown, S.L. (2009). Insect
rate-temperature relationships: environmental variation and the metabolic
theory of ecology. Am. Nat., 174, 819–835.
Johansson, F. (1993). Intraguild predation and cannibalism in odonate
larvae: effects of foraging behaviour and zooplankton availability. Oikos, 66,
80–87.
Knies, J.L. & Kingsolver, J.G. (2010). Erroneous Arrhenius: modified Arrhenius
Model best explains the temperature dependence of ectotherm fitness. Am.
Nat., 176, 227–233.
Thermal dependence of damselfly growth rates 71
Miller, P.L. & Miller, C.A. (1981). Field observations on copulatory behaviour in
Zygoptera, with an examination of the structure and activity of male genitalia.
Odonatologica, 10, 201–218.
Nilsson-Örtman, V., Stoks, R., De Block, M. & Johansson, F. (2012). Generalists
and specialists along a latitudinal transect: patterns of thermal adaptation in six
species of damselflies. Ecology, 93, 1340–1352.
Norling, U. (1984). Life history patterns in the northern expansion of
dragonflies. Adv. Odonatol., 2, 127–156.
Pinheiro, J.C. & Bates, D.M. (2000). Mixed-Effects Models in S and S-PLUS.
Springer Verlag, New York.
R core team. (2008). R: A Language and Environment for Statistical Computing. R
Foundation for Statistical Computing, Vienna, Austria.
Rose, K.E., Atkinson, R.L., Turnbull, L.A. & Rees, M. (2009). The costs and
benefits of fast living. Ecol. Lett., 12, 1379–1384.
Savage, V.M., Gillooly, J.F., Brown, J.H., West, G.B. & Charnov, E.L. (2004).
Effects of body size and temperature on population growth. Am. Nat., 163,
429–441.
Schielzeth, H. (2010). Simple means to improve the interpretability of regression
coefficients. Methods Ecol. Evol., 1, 103–113.
Somero, G.N. & Low, P.S. (1977). Eurytolerant proteins: mechanisms for
extending the environmental tolerance range of enzyme-ligand interactions.
Am. Nat., 111, 527–538.
Stoks, R., Block, M.D., Van De Meutter, F. & Johansson, F. (2005). Predation
cost of rapid growth: behavioural coupling and physiological decoupling.
J. Anim. Ecol., 74, 708–715.
Suhling, F. & Lepkojus, S. (2001). Differences in growth and behaviour influence
asymmetric predation among early-instar dragonfly larvae. Can. J. Zool., 79,
854–860.
Terblanche, J., Janion, C. & Chown, S. (2007). Variation in scorpion metabolic
rate and rate–temperature relationships: implications for the fundamental
equation of the metabolic theory of ecology. J. Evol. Biol., 20, 1602–1612.
Waage, J.K. (1986). Evidence for widespread sperm displacement ability among
Zygoptera (Odonata) and the means for predicting its presence. Biol. J. Linn.
Soc., 28, 285–300.
SUPPORTING INFORMATION
Additional Supporting Information may be downloaded via the online
version of this article at Wiley Online Library (www.ecologyletters.com).
As a service to our authors and readers, this journal provides supporting information supplied by the authors. Such materials are
peer-reviewed and may be re-organised for online delivery, but are
not copy-edited or typeset. Technical support issues arising from
supporting information (other than missing files) should be
addressed to the authors.
Editor, Brian Enquist
Manuscript received 5 June 2012
First decision made 14 July 2012
Manuscript accepted 5 September 2012
© 2012 Blackwell Publishing Ltd/CNRS