Temporal autocorrelation and effective sample size

How autocorrelation shrinks the effective sample size of an ecological time series, why the naive correlation test between two series is anticonservative, and how to correct it in base R.

A field ecologist rarely gets independent observations through time. Counts, temperatures, water levels and phenological dates recorded month after month carry over: a warm April makes a warm May more likely, a high vole year props up the next. That carry-over is temporal autocorrelation, and ignoring it does more than annoy a reviewer. It quietly shrinks how much information the series actually holds, and it makes ordinary significance tests reject far too often. This post builds the autocorrelation and partial autocorrelation functions by hand, shows how they translate into an effective sample size, and demonstrates the single most common trap: two unrelated series that appear correlated simply because both are autocorrelated.

We stay in base R throughout. Everything here is stats: arima.sim, acf, pacf, lm, cor.

An autocorrelated series and its correlogram

A first-order autoregressive process, AR(1), is the simplest useful model of carry-over. Each value is a fraction of the previous one plus fresh noise: y[t] = phi * y[t-1] + e[t]. The single parameter phi sets the persistence. We simulate ten years of monthly values with strong persistence.

set.seed(4240)
phi <- 0.7; n <- 120
y <- as.numeric(arima.sim(list(ar = phi), n = n))

K <- 20
acf_h <- acf_hand(y, K)                 # by hand
acf_b <- as.numeric(acf(y, lag.max = K, plot = FALSE)$acf)
pacf_h <- pacf_dl(acf_h, K)
pacf_b <- as.numeric(pacf(y, lag.max = K, plot = FALSE)$acf)

acf_match  <- max(abs(acf_h - acf_b))
pacf_match <- max(abs(pacf_h - pacf_b))

The hand computation reproduces stats::acf to 1.1^{-16} and stats::pacf to 3.8^{-16}, so we understand exactly what the built-in functions return. The estimator is the standard biased one: the lag-k autocovariance is divided by the lag-0 sum, not by n - k.

The two functions read a series in complementary ways. The autocorrelation function (ACF) measures the correlation between the series and a copy of itself shifted by k steps. For an AR(1) process it should decay geometrically as phi^k. The lag-1 sample value here is 0.695, close to the true phi of 0.7. The partial autocorrelation function (PACF) removes the intervening lags, asking what lag k adds once lags 1 to k - 1 are accounted for. For an AR(1) process the PACF has a single spike at lag 1 and is near zero afterwards: here lag 1 is 0.695, lag 2 is -0.012 and lag 3 is 0.106. That contrast, a slow ACF decay against a sharp PACF cutoff, is the signature that tells you the order of an autoregressive model. The next post uses it directly.

Figure 1: Sample ACF and PACF of an AR(1) series. The ACF decays like phi^k; the PACF spikes once at lag 1 then falls silent.

How much information does the series really hold?

The dashed band in the ACF panel is the interval a naive analysis would use to decide whether an autocorrelation is real: plus or minus 1.96 / sqrt(n), built on the assumption that the data are independent. That assumption is exactly what fails here. When observations carry over, each new value repeats part of what the last one already told you, so the series behaves as if it had fewer independent points than its length suggests.

For estimating a mean from an AR(1) series, the variance of the sample mean is inflated relative to the independent case by a factor that depends only on phi:

vif_mean <- (1 + phi) / (1 - phi)
n_eff <- n * (1 - phi) / (1 + phi)

With phi of 0.7, the variance inflation factor is 5.67 and the effective sample size drops from 120 to about 21. Ten years of monthly counts carry roughly the information of 21 independent observations. A short simulation confirms it: we draw many AR(1) series, take each sample mean, and compare its true spread against the naive standard error sd / sqrt(n).

B <- 4000
means <- naive_se <- numeric(B)
set.seed(4240)
for (b in seq_len(B)) {
  z <- as.numeric(arima.sim(list(ar = phi), n = n))
  means[b] <- mean(z)
  naive_se[b] <- sd(z) / sqrt(n)
}
emp_sd <- sd(means)
ratio_se <- emp_sd / mean(naive_se)
cover_naive <- mean(abs(means) <= 1.96 * naive_se)
cover_eff <- mean(abs(means) <= 1.96 * naive_se * sqrt(vif_mean))

The true spread of the sample mean is 0.302, but the naive standard error averages 0.125. The naive figure is too small by a factor of 2.42, which matches sqrt(VIF) of 2.38. The practical damage: a naive 95% confidence interval for the mean covers the truth only 58% of the time, not 95%. Rescaling the standard error by sqrt(VIF), equivalently using the effective sample size, restores coverage to 94%.

The classic trap: two series that are not related

The most costly consequence is not about a single mean but about comparing two series. Suppose you record two monitoring variables, each persistent through time, and ask whether they move together. Regressing one on the other with lm, or testing their correlation, uses a standard error built for independent data. Because both series drift and linger, spurious agreement is common.

set.seed(4240)
x1 <- as.numeric(arima.sim(list(ar = phi), n = n))
y1 <- as.numeric(arima.sim(list(ar = phi), n = n))   # generated independently of x1
r_demo <- cor(x1, y1)
p_ols <- summary(lm(y1 ~ x1))$coefficients[2, 4]
p_corr <- mod_t_p(x1, y1)

These two series share nothing by construction, yet their sample correlation is 0.16 and the ordinary slope test returns a p-value of 0.083, close to conventional significance. The Dutilleul (1993) modified test, which replaces the sample size with an effective one derived from both series’ autocorrelation, returns 0.27, comfortably non-significant. To see that this is systematic rather than one unlucky draw, we repeat the experiment many times and count how often each test rejects the (true) null of no relationship.

B <- 3000
p_naive <- p_dutil <- p_white <- numeric(B)
set.seed(4241)
for (b in seq_len(B)) {
  x <- as.numeric(arima.sim(list(ar = phi), n = n))
  y <- as.numeric(arima.sim(list(ar = phi), n = n))
  p_naive[b] <- summary(lm(y ~ x))$coefficients[2, 4]
  p_dutil[b] <- mod_t_p(x, y)
  xw <- rnorm(n); yw <- rnorm(n)                       # white-noise control
  p_white[b] <- summary(lm(yw ~ xw))$coefficients[2, 4]
}
rej_naive <- mean(p_naive < 0.05)
rej_dutil <- mean(p_dutil < 0.05)
rej_white <- mean(p_white < 0.05)

set.seed(4242)
r_ar <- r_wn <- numeric(B)
for (b in seq_len(B)) {
  r_ar[b] <- cor(as.numeric(arima.sim(list(ar = phi), n = n)),
                 as.numeric(arima.sim(list(ar = phi), n = n)))
  r_wn[b] <- cor(rnorm(n), rnorm(n))
}

The naive slope test rejects 23% of the time when it should reject 5%, an error rate inflated roughly 5-fold. The white-noise control sits at 6%, confirming that the test itself is fine and the fault is the autocorrelation. The effective-df correction pulls the rate back to 5%. The left panel below shows why: the null distribution of the sample correlation is far wider for autocorrelated series than for white noise. A correlation of 0.2 between two independent series happens 20% of the time when both are AR(1), against 3% for white noise.

Figure 2: Left: the null distribution of the sample correlation between two independent series is much wider under autocorrelation. Right: the naive t-test rejects the true null far too often; an effective-degrees-of-freedom correction restores the nominal rate.

What to carry away

Temporal autocorrelation is not a nuisance to be swept into the residuals. It sets how much a series actually tells you. A persistent monthly record can hold the information of a couple of dozen independent points, and any test that assumes independence, whether a confidence interval for a mean, a correlation, or an ordinary regression slope, will be too confident. The effective-degrees-of-freedom correction shown here is the temporal cousin of the spatial adjustment used when neighbouring sites are not independent: same logic, same fix. Detecting that a real ecological relationship exists between two time series takes an honest accounting of the autocorrelation, not just a small p-value.

Two habits follow. First, always look at the correlogram before testing anything on a time series; the ACF and PACF tell you both how strong the carry-over is and what model might describe it. Second, when you compare two series, correct for their autocorrelation, whether through an effective sample size, a generalised least squares fit with an autoregressive error, or the prewhitening approach we return to when checking a time series model. The next post reads the ACF and PACF the other way round, using them to identify and fit an autoregressive-moving-average model.

Related tutorials

• Fitting ARIMA models to ecological time series
• Spectral analysis of population cycles
• Checking a time series model
• Regression with spatially correlated errors (GLS)

References

• Bartlett 1946 Supplement to the Journal of the Royal Statistical Society 8(1):27-41 (10.2307/2983611)
• Box, Jenkins, Reinsel & Ljung 2015 Time Series Analysis: Forecasting and Control, 5th ed, Wiley, ISBN 978-1-118-67502-1
• Clifford, Richardson & Hemon 1989 Biometrics 45(1):123-134 (10.2307/2532039)
• Dutilleul 1993 Biometrics 49(1):305-314 (10.2307/2532625)
• Legendre & Legendre 2012 Numerical Ecology, 3rd ed, Elsevier, ISBN 978-0-444-53868-0