Checking a thermal performance analysis

Three checks for a thermal performance curve in R: the extrapolated critical maximum, curve-family dependence, and measurement error in body temperature.

A thermal performance curve delivers three headline numbers, the optimum, the breadth, and the critical thermal maximum, and they carry very different weights of evidence. The optimum and breadth sit inside the data. CTmax usually does not. This post is three checks that a performance-curve result should survive before those numbers are used to rank species or predict range limits.

Check 1: how much data supports CTmax?

CTmax is where the fitted curve reaches zero, and that point is almost always past the hottest temperature anyone measured. Refit the curve after dropping the warmest data and watch the estimate move.

fitTm <- function(T, yy) {
  f <- nls(yy ~ aa*T*(T-t0)*sqrt(pmax(tm-T, 0)),
           start = list(aa = 6e-4, t0 = 6, tm = 44), control = nls.control(warnOnly = TRUE))
  s <- summary(f)$coefficients
  c(est.tm = s["tm","Estimate"], se.tm = s["tm","Std. Error"], coef(f))
}
full <- fitTm(Td, y)
r38  <- fitTm(Td[Td <= 37.5], y[Td <= 37.5])
r35  <- fitTm(Td[Td <= 35],   y[Td <= 35])
rbind(full = round(full[1:2], 2), drop_37.5 = round(r38[1:2], 2), drop_35 = round(r35[1:2], 2))

          est.tm se.tm
full       42.11  0.11
drop_37.5  42.75  0.31
drop_35    42.95  0.63

The hottest measurement is 40 degrees C, yet CTmax comes out near 42.1: it is already an extrapolation 2 degrees beyond the data, and only 12 of the 52 points sit above the optimum at all. Drop the descending-limb points and the estimate reads 43.0, but the real tell is the standard error, which grows from 0.1 to 0.6 degrees. If sampling stops near the optimum, CTmax is not estimated, it is guessed. The check: report CTmax with its interval, and refuse to compare CTmax across species unless each is anchored by measurements on the descending limb.

cf <- function(cc) { g <- seq(cc["t0"], cc["tm"], 0.05); data.frame(T = g, P = cc["aa"]*g*(g-cc["t0"])*sqrt(pmax(cc["tm"]-g, 0))) }
cA <- cf(full[3:5]); cB <- cf(r35[3:5])
p_ctmax <- ggplot(data.frame(Td, y), aes(Td, y)) +
  geom_point(colour = col_ref, size = 2, alpha = 0.7) +
  geom_line(data = cA, aes(T, P), colour = col_full, linewidth = 1) +
  geom_line(data = cB, aes(T, P), colour = col_cut, linewidth = 1, linetype = "22") +
  geom_vline(xintercept = 35, colour = col_ref, linetype = "dotted") +
  annotate("text", x = 35.3, y = 0.15, label = "data dropped above here", hjust = 0, colour = col_ref, size = 3.2) +
  labs(x = "Temperature (C)", y = "Performance", title = "Full-data fit and refit without the hottest points")
p_ctmax

Figure 1: The upper limit is pinned by a handful of hot points. The full-data fit (green) reaches zero near 42 degrees; the refit without the two hottest temperatures (ochre) reaches zero further right, and its standard error is the larger of the two.

Check 2: does the answer depend on the curve family?

Because CTmax lives in the extrapolated region, the model chosen to get there matters more than anywhere else on the curve. Fit the asymmetric Briere and a symmetric Gaussian to the same data and compare where each says performance vanishes.

cb <- coef(nls(y ~ aa*Td*(Td-t0)*sqrt(pmax(tm-Td, 0)), start = list(aa = 6e-4, t0 = 6, tm = 44)))
cg <- coef(nls(y ~ A*exp(-((Td-mo)^2)/(2*s^2)),        start = list(A = 2, mo = 34, s = 8)))
c(briere_CTmax = unname(cb["tm"]), gaussian_upper = unname(cg["mo"] + cg["s"]*sqrt(-2*log(0.05))))

  briere_CTmax gaussian_upper 
      42.10962       59.23226 

The Briere model puts the upper limit at 42.1 degrees; the Gaussian, reaching five per cent of its peak, at 59.2. Inside the measured range the two curves stay close without being identical, and the symmetric family puts its optimum below the asymmetric one, pulled down by the long left shoulder it cannot represent. On the extrapolated limit they part by 17 degrees, from identical data, much the larger of the two disagreements. When a number depends this strongly on a modelling choice, fit more than one family and report the spread, rather than trusting the decimal a single model returns.

Check 3: is the temperature measured with error?

Performance is regressed on body temperature, but body temperature is hard to measure and is often replaced by air temperature. Error in the predictor distorts a nonlinear curve, the same regression dilution that attenuates a straight-line slope, and it biases every summary.

set.seed(1394)                                  # not 394: that seed would replay the noise already in y
Terr <- Td + rnorm(length(Td), 0, 2.5)          # body temperature measured with error
ce <- coef(nls(y ~ aa*Terr*(Terr-t0)*sqrt(pmax(tm-Terr, 0)),
               start = list(aa = 6e-4, t0 = 6, tm = 46), control = nls.control(warnOnly = TRUE)))
round(c(true_CTmax = 42, fitted_CTmax = unname(ce["tm"])), 1)

  true_CTmax fitted_CTmax 
        42.0         44.9 

A temperature standard deviation of 2.5 degrees pushes the apparent CTmax up to 44.9 degrees. A curve fitted to air temperature will look more heat-tolerant than the animal really is. The check: use body temperature, not a proxy, and where a proxy is unavoidable, treat the estimated breadth and optimum as approximate (measurement error and regression dilution has the general version).

ggplot() +
  geom_line(data = data.frame(T = Tg, P = briere(Tg, a0, T0, Tm)), aes(T, P), colour = col_true, linewidth = 1) +
  geom_line(data = data.frame(T = TgE, P = PgE), aes(T, P), colour = col_err, linewidth = 1, linetype = "22") +
  annotate("text", x = 14, y = Pmax0*0.96, label = "true", colour = col_true, size = 3.6) +
  annotate("text", x = 14, y = PmaxE*0.9, label = "fitted with temperature error", colour = col_err, hjust = 0, size = 3.4) +
  labs(x = "Temperature (C)", y = "Performance", title = "True curve and curve fitted with temperature error")

Figure 2: The true curve (green) and the curve fitted to error-laden temperatures (ochre): the ochre curve reaches zero further right.

The deeper caveat

These three checks bound the numbers a curve returns. They do not fix what the curve is. A performance curve is measured by holding an animal at a constant temperature and reading a rate, and the two earlier posts showed why that is not the field: performance under fluctuating temperature is not performance at the mean, and the curve itself shifts with acclimation and can be edited by behaviour (Sinclair et al. 2016 Ecology Letters 19(11):1372-1385). Run the three checks so the numbers are honest, then remember that even an honest constant-temperature curve is a laboratory object, and the animal lives outside.

References

• Kingsolver JG 2009. The American Naturalist 174(6):755-768 (10.1086/648310).
• Schoolfield RM et al. 1981. Journal of Theoretical Biology 88(4):719-731 (10.1016/0022-5193(81)90246-0).
• Sinclair BJ et al. 2016. Ecology Letters 19(11):1372-1385 (10.1111/ele.12686).

Related tutorials

• Thermal performance curves in R
• Jensen’s inequality and thermal variability
• Thermal safety margins and warming
• Measurement error and regression dilution