Autoregressive Correlation

Autoregressive Correlation and AR(1) Structures for Trial Costs, Utilities and Health Spending

Trial-based economic evaluations record utility and resource use at several visits for each patient, and health systems report spending for the same countries or regions year after year. Measurements on the same unit are not independent, and an autoregressive structure captures that dependence with a single parameter that decays with distance in time. The AR(1) form is used in generalised estimating equations and mixed models, where it competes with exchangeable and unstructured alternatives. In an illustrative worked example the choice changes the variance of a QALY total and of a change from baseline, and the same structure appears as serial correlation in analyses of health spending.

Why repeated costs and utilities move together within a patient

A patient who reports poor health at one visit tends to report poor health at the next, because the condition, comorbidities and reporting habits persist. Resource use behaves in the same way, since a patient on an expensive care pathway in one quarter is often still on it in the next. In a cost-utility analysis alongside a trial, these repeated health state utility values and costs are summed into totals for each patient, so the dependence between visits passes straight into the precision of those totals. A model that treats each visit record as an independent observation misstates the standard error of any effect built from those records, while a model that represents the correlation can use the data more efficiently.

The AR(1) structure and its continuous-time form

The first-order autoregressive structure uses a single parameter for the correlation between any two visits. The correlation is highest for adjacent visits and falls geometrically with the number of steps between them.

$$\mathrm{Corr}(Y_{ij}, Y_{ik}) = \rho^{|j-k|}$$

where $Y_{ij}$ and $Y_{ik}$ are the outcomes (for example utility or cost) of patient $i$ at visits $j$ and $k$, $|j-k|$ is the number of visit steps between them, and $\rho$ is the correlation between adjacent visits, with $-1 < \rho < 1$. For four equally spaced visits the correlation matrix is as follows.

$$R(\rho) = \begin{pmatrix} 1 & \rho & \rho^2 & \rho^3 \\ \rho & 1 & \rho & \rho^2 \\ \rho^2 & \rho & 1 & \rho \\ \rho^3 & \rho^2 & \rho & 1 \end{pmatrix}$$

where the entry in row $j$ and column $k$ is the correlation between visits $j$ and $k$. With a common variance $\sigma^2$ at every visit, the covariance between visits is $\sigma^2 \rho^{|j-k|}$. In repeated health measures $\rho$ is normally positive, so measurements close in time are more alike than measurements far apart, and the correlation tends towards zero as the gap grows.

Trial visits are often unequally spaced, for example at 0, 3, 6 and 12 months, and counting visit steps would treat a 3-month gap and a 6-month gap as the same distance. The continuous-time version uses elapsed time instead, as follows.

$$\mathrm{Corr}(Y_{ij}, Y_{ik}) = \phi^{|t_{ij} - t_{ik}|}$$

where $t_{ij}$ and $t_{ik}$ are the times of visits $j$ and $k$ for patient $i$ in a chosen unit such as months, and $\phi$ is the correlation between two measurements one time unit apart. Liang and Zeger gave this form as a working model in their 1986 paper on generalised estimating equations, noting that it accommodates an arbitrary number and spacing of observations.

How the structure enters GEE working correlation and mixed-model residuals

Autoregressive correlation is used in two families of models for repeated health economic outcomes. They treat the structure differently, which changes how much rests on choosing it well.

In generalised estimating equations, the AR(1) matrix is one choice of working correlation, a stated guess at the within-patient correlation. Liang and Zeger showed that the regression estimates and their variance estimate stay consistent whatever working correlation is chosen, provided the mean model is correct; the variance comes from the sandwich estimator, and a working structure close to the truth gives more efficient estimates. The protection has a condition. When the working correlation is wrong, the authors required data to be missing completely at random, which trial dropout rarely satisfies.

In a mixed effects model, the structure is part of a full likelihood. AR(1) is usually placed on the within-patient residuals, often with a random intercept for each patient. A random intercept alone implies an exchangeable pattern, so the combination gives correlation that decays towards a floor set by stable differences between patients, a model set out in the longitudinal texts of Diggle and colleagues and of Fitzmaurice, Laird and Ware. Gabrio, Plumpton, Banerjee and Leurent describe linear mixed models for trial-based cost-effectiveness analysis that give correct inferences when data are missing at random, with QALYs and total costs recovered as weighted linear combinations of the estimated visit means. With incomplete follow-up, the assumed covariance links missing visits to observed ones, so a structure that misdescribes the decay of correlation can shift the estimated QALY difference as well as its standard error.

Choosing between exchangeable, AR(1) and unstructured patterns

The main alternatives differ in how correlation changes with distance between visits and in how many parameters they need. The table summarises the three most used for trial outcomes measured at $T$ visits.

Liang and Zeger noted that the unstructured choice gives the most efficient estimator but is only useful when there are few observation times; with five visits it needs 10 correlation parameters. The mixed model for repeated measures used in many trials adopts an unstructured covariance, and Gabrio and colleagues note that when such a model does not converge, simpler covariance structures could be considered. AR(1) is one such option, with one correlation parameter whatever the number of visits.

A table of residual correlations by pair of visits is the first guide. Correlation that falls with the gap points to AR(1), a flat pattern to exchangeable, and a fall that levels off to a random intercept plus AR(1) residuals. For mixed models, the Akaike information criterion can compare covariance structures fitted to the same mean model. GEE has no likelihood, so Pan proposed a quasi-likelihood version, QIC; the geepack documentation for R states that QIC is used to select a correlation structure, with CIC as an alternative.

Worked example with one QALY total and one change from baseline

The figures in this example are illustrative. A trial measures utility at baseline, 6 months and 12 months, with a standard deviation of 0.20 at each visit, so $\sigma^2 = 0.04$. Each patient's 12-month QALY is the area under the utility curve by the trapezium rule, as follows.

$$Q_i = 0.25\,U_{i1} + 0.5\,U_{i2} + 0.25\,U_{i3}$$

where $Q_i$ is the 12-month QALY total for patient $i$, and $U_{i1}$, $U_{i2}$ and $U_{i3}$ are that patient's utilities at baseline, 6 and 12 months; the weights are in years, with half of each 6-month interval attached to each end. The variance of this weighted sum is given below.

$$\mathrm{Var}(Q_i) = \sigma^2 \sum_{j=1}^{3} \sum_{k=1}^{3} w_j w_k r_{jk}$$

where $w_j$ and $w_k$ are the visit weights (0.25, 0.5 and 0.25) and $r_{jk}$ is the correlation between visits $j$ and $k$, equal to 1 when $j = k$.

1. Independence. With zero correlation between different visits only the squared weights remain.

$$0.04 \times (0.0625 + 0.25 + 0.0625) = 0.04 \times 0.375 = 0.015$$

where 0.015 is the variance of the QALY total, a standard deviation of about 0.122.

2. AR(1) with adjacent-visit correlation 0.7. Baseline and 12 months are two steps apart, so their correlation is 0.7 × 0.7 = 0.49. The cross-products are 0.25 × 0.5 × 0.7 = 0.0875 for each adjacent pair and 0.25 × 0.25 × 0.49 = 0.030625 for the outer pair, each counted twice, which gives the following.

$$0.04 \times (0.375 + 2 \times (0.0875 + 0.0875 + 0.030625)) = 0.04 \times 0.78625 = 0.03145$$

where 0.03145 is the variance under AR(1), a standard deviation of about 0.177.

3. Exchangeable with an assumed correlation of 0.7. With the same correlation of 0.7 assumed for every pair, the outer cross-product is 0.25 × 0.25 × 0.7 = 0.04375 and the variance is then as follows.

$$0.04 \times (0.375 + 2 \times (0.0875 + 0.0875 + 0.04375)) = 0.04 \times 0.8125 = 0.0325$$

where 0.0325 is the variance under the exchangeable pattern, a standard deviation of about 0.180.

4. Between-arm standard error. With 100 patients per arm and equal variances, the standard error of the difference in mean QALYs is the square root of twice the variance divided by 100 and comes to about 0.0173 under independence and about 0.0251 under AR(1). An analysis that models visit records as independent and uses model-based (naive) standard errors would understate the standard error by about 31 per cent, giving a confidence interval about 70 per cent of its proper width. An independence working model in GEE with sandwich standard errors does not have this problem.

5. Change from baseline to 12 months. The variance of a difference between two visits is twice $\sigma^2(1 - r_{13})$, where $r_{13}$ is the baseline to 12-month correlation, as follows.

$$2 \times 0.04 \times (1 - 0.49) = 0.0408$$

where 0.0408 is the variance of the change under AR(1). An exchangeable structure assumed to share the same correlation of 0.7 gives 2 × 0.04 × 0.3 = 0.024, about 41 per cent below the AR(1) value. A fitted exchangeable model would not use 0.7. With data following AR(1), it would estimate roughly the average pairwise correlation, (0.7 + 0.7 + 0.49) / 3 = 0.63, giving 2 × 0.04 × 0.37 = 0.0296, still about 27 per cent below the AR(1) variance.

The two correlated structures nearly agree on the QALY total (0.03145 against 0.0325 at an assumed 0.7, or 0.03075 with the fitted 0.63) but differ clearly on the 12-month change. The choice matters most for contrasts between distant visits, such as a change from baseline or a baseline-adjusted effect on final utility, and least for time-averaged totals. The independence result is the risk when visit records are modelled one by one; when each patient's QALY is calculated first and then compared between arms, the spread of the totals already reflects the correlation.

Serial correlation in time-series and panel analyses of health spending

The same dependence appears when the unit observed repeatedly is a country, region, hospital or practice. Spending, prescribing and admission levels tend to carry over from one period to the next, so regression residuals can be serially correlated. In time-series notation the AR(1) error process is as follows.

$$e_t = \rho\, e_{t-1} + u_t$$

where $e_t$ is the regression error in period $t$, $\rho$ is the autoregressive coefficient, and $u_t$ is a new disturbance independent of earlier periods. For a stationary process this implies a correlation of $\rho^s$ between errors $s$ periods apart, the same geometric pattern as the trial structure.

Bertrand, Duflo and Mullainathan generated placebo laws in state-level wage data and found that conventional difference-in-differences standard errors severely understated the true variability, with an effect significant at the 5 per cent level for up to 45 per cent of the placebo interventions. Corrections imposing a specific parametric form on the time-series process did not perform well. Corrections allowing a general covariance across periods worked well for moderate numbers of states, and collapsing the data into before and after periods worked well even for small numbers. A health spending panel of around 20 countries or regions is a small-cluster setting, so the collapsing approach, rather than a simple AR(1) adjustment, is the safer reference point there.

For interrupted time series evaluations, Bernal, Cummins and Gasparrini note that autocorrelation in epidemiological data is often largely explained by seasonality, and recommend residual plots, the partial autocorrelation function and, where autocorrelation remains, Prais regression or ARIMA models.

What the choice changes in a cost-effectiveness result

The structure reaches a decision through the uncertainty around incremental costs and QALYs. Those standard errors feed confidence intervals, acceptability curves and the probability that an intervention is cost-effective at a given threshold, so an understated variance makes a result look more certain than the data allow.

When uncertainty is estimated by bootstrapping, resampling whole patients rather than visit records keeps each sequence intact, so the bootstrap distribution reflects the within-patient correlation whatever structure the point-estimation model assumed. Reports should state the structure used, why it was preferred, and whether incremental results change under a plausible alternative.

Where autoregressive structures are misread

Several errors recur in applied work. Most come from treating AR(1) as a neutral default rather than a specific claim about how correlation decays.

• Assuming correlation dies out. Pure AR(1) drives the correlation between distant visits towards zero, yet stable differences between patients keep utilities correlated over long follow-up; a random intercept with AR(1) residuals often fits better.
• Confusing a correlation structure with a dynamic model. An AR(1) residual structure describes dependence left after the mean model; including the previous outcome as a covariate answers a conditional question, which Liang and Zeger distinguished from their marginal approach.
• Reading an omitted trend as serial correlation. An unmodelled trend or seasonal pattern leaves autocorrelated residuals, so the mean model needs checking before an autoregressive error is added.

Sources

• Bernal JL, Cummins S, Gasparrini A. Interrupted time series regression for the evaluation of public health interventions: a tutorial. International Journal of Epidemiology. 2017;46(1):348-355.
• Bertrand M, Duflo E, Mullainathan S. How much should we trust differences-in-differences estimates? Quarterly Journal of Economics. 2004;119(1):249-275.
• Diggle PJ, Heagerty P, Liang KY, Zeger SL. Analysis of Longitudinal Data. 2nd ed. Oxford: Oxford University Press; 2002.
• Fitzmaurice GM, Laird NM, Ware JH. Applied Longitudinal Analysis. 2nd ed. Hoboken, NJ: Wiley; 2011.
• Gabrio A, Plumpton C, Banerjee S, Leurent B. Linear mixed models to handle missing at random data in trial-based economic evaluations. Health Economics. 2022;31(6):1276-1287.
• Højsgaard S, Halekoh U, Yan J. geepack: Generalized Estimating Equation Package. R package version 1.3.11, QIC help page. 2024. https://search.r-project.org/CRAN/refmans/geepack/html/QIC.html
• Liang KY, Zeger SL. Longitudinal data analysis using generalized linear models. Biometrika. 1986;73(1):13-22.
• Pan W. Akaike's information criterion in generalized estimating equations. Biometrics. 2001;57(1):120-125.