Meta-analysis of diagnostic test accuracy

A diagnostic study reports a pair, not a number

An ordinary meta-analysis pools a scalar: one effect measure with its variance (see effect measures and their variances). A diagnostic accuracy study does not work that way. What it hands over is a 2×2 table, and the sensitivity and specificity read off that table are two outcomes fixed by one and the same cut-off (see the 2x2 table, sensitivity and specificity and choosing a cut-off).

Move the cut-off towards calling more people positive and sensitivity rises while specificity falls; move it the other way and they swap. So when a dozen studies each use a different cut-off — or nominally the same cut-off, read more or less strictly — the two numbers travel along a curve together across studies rather than scattering independently.

The consequence: handing sensitivity to one random-effects model and specificity to another amounts to declaring that the two are unrelated. What this page does is measure the price of that declaration rather than repeat the slogan. The data are mada::AuditC — 14 studies of the AUDIT-C alcohol screening questionnaire, 18,332 people in total, of whom 2,071 had the target condition.

What these studies look like

The sensitivity and specificity columns are the uncorrected observed values — what the studies actually reported. The DOR (diagnostic odds ratio) column is computed after adding 0.5 to every cell, because 2 studies report FN = 0 and their uncorrected DOR would be infinite. The studies differ enormously in size: the smallest has 112 participants, the largest 6,955, and that one study alone carries 38% of the total sample. That fact returns in a later section.

The mechanism: move the cut-off and the two numbers move opposite ways

Before fitting anything, plot each study’s logit sensitivity against its logit specificity:

The Pearson correlation is -0.737 (95% CI -0.911 to -0.339, p = 0.0026, df = 12) and Spearman’s ρ is -0.732. The fitted line has slope -0.280 and R² 0.543.

This negative correlation is not held up by one high-leverage study, which matters enough to be checked rather than assumed. Refitting with each study dropped in turn leaves the correlation between -0.805 and -0.569 — the least negative value arises when study 8 is removed, which is the largest one and the one sitting furthest out along the logit sensitivity axis. Removing both zero-cell studies (studies 7 and 8) together leaves a correlation of -0.750.

In other words, the argument that follows rests on positive evidence you can cite, not on a disclaimer.

Where the threshold effect comes from, and one gap in these data

The most direct reading of that downward slope is that the studies used different cut-offs. In practice several sources leave the same footprint, and a reader has to keep them apart:

• Explicitly different cut-offs — the same questionnaire or biomarker, scored against a lower threshold in some studies and a higher one in others.
• The same cut-off read more or less strictly — commonest in imaging and pathology: one nominal criterion, different effective thresholds between readers.
• Different reference standards — change what counts as truly having the target condition and both measures move together.
• A different spectrum of disease — symptomatic inpatients versus community screening changes both measures at once.

The first two are threshold effects in the strict sense; the last two are not, but what they leave in the data looks much the same. The bivariate model does not distinguish the four sources. All it does is acknowledge that the two measures are correlated and build that into the variance structure; separating the sources requires meta-regression with cut-off, reference standard or population characteristics as covariates (the machinery is the same as in the heterogeneity page, with two outcomes instead of one).

Pooling separately: the summary point barely moves

Now fit two separate univariate random-effects models (REML, logit scale) and one bivariate model with mada::reitsma(), holding the continuity correction the same across both, and put the results side by side:

The summary point barely moves; the joint statement does. Separate pooling does sit above the bivariate model on both measures, the direction the textbooks describe, but by 0.19 and 0.11 percentage points — a difference nobody could spot in an abstract, and one that vanishes on rounding to two decimal places. This page will not say that separate pooling badly overstates accuracy, because on these data it does not.

Nor are these data unusually well behaved. With no covariates in the model, the marginal point estimates from a bivariate fit are expected to sit very close to those from two univariate fits: the fixed-effect parts are estimating the same thing, and the correlation enters through the variance structure rather than the location. Correlation changes the shape of the uncertainty, so the difference only becomes visible one level up, at the joint statement.

For the record, between-study heterogeneity is large on both measures: τ² = 1.215 with I² = 93.5% for sensitivity, τ² = 0.427 with I² = 98.6% for specificity. I² is a proportion rather than an amount (see heterogeneity — I², τ² and prediction intervals), but τ² that size on the logit scale is genuinely wide spread, and the prediction ellipse below turns it into something visible.

The joint statement is what breaks

The commonest shape of a published result is two sentences: pooled sensitivity 0.893 (95% CI 0.814–0.941), pooled specificity 0.781 (95% CI 0.715–0.835). A reader then does the natural thing: crosses the two intervals into a rectangle and says the test’s true performance lies somewhere in there.

That step is wrong, and it is wrong on two levels.

First: crossing two 95% intervals does not give a 95% region. Even if the two were genuinely independent, the joint coverage would be 90.25% (0.95 squared). That has nothing to do with the model; it is arithmetic.

Second, and this is the point of the page: the shape is wrong. Drawing 200,000 Monte Carlo samples from the bivariate sampling distribution, the rectangle’s actual joint coverage is 91.9%. The same draws give the bivariate confidence ellipse 94.93%, which is the check on the procedure itself — that number should come out near 95%, and it does.

The rectangle’s area is 0.0152 against 0.0117 for the ellipse, a ratio of 1.30. Roughly a third more area, and less coverage — because the extra area is in the wrong places.

Panel B above shows it at a glance. The upper left corner of the rectangle is high sensitivity together with high specificity; the lower right corner is both low at once. Those are exactly the combinations the negative correlation rules out: the ellipse excludes them and the rectangle includes them. What the rectangle leaves out instead are the two ends of the ellipse’s long axis — high sensitivity with low specificity, and low sensitivity with high specificity — which are the possibilities the data actually support.

What happens if you take a corner of the rectangle to a likelihood ratio

The problem still sounds abstract until the rectangle is pushed into a quantity clinicians use. The likelihood ratio is such a quantity, and its definition consumes both measures at once (see likelihood ratios): the positive likelihood ratio is sensitivity divided by the false positive rate, and the negative one is 1 minus sensitivity divided by specificity.

Computed from the bivariate summary point, the positive likelihood ratio is 4.04 and the negative one is 0.140. Now take the most flattering corner of the rectangle instead — sensitivity at its upper bound 0.941, specificity also at its upper bound 0.835 — and the positive likelihood ratio becomes 5.69, some 1.41 times the summary value, while the negative one becomes 0.071, about half of it.

That pair is not an optimistic estimate; it is a combination the model says will not occur together. A reader who crosses the two intervals ends up with a persuasive-looking test that no data support.

Hence a simple working rule: for anything that needs both measures at once — likelihood ratios, predictive values, decision curves — the interval cannot be built by arithmetic on two sets of bounds. It has to be computed on the model’s joint distribution: sample from the bivariate sampling distribution, compute the likelihood ratio for each draw, take quantiles; or use the delta method. The simulation snippet in the previous section becomes exactly this with a small change.

The confidence ellipse and the prediction ellipse are not the same thing

That ellipse has a more frequently misread sibling. Both are centred on the same point; they differ in what they are the uncertainty of:

• The confidence ellipse uses vcov(fixed effects) and describes where the average operating point lies. It shrinks as studies accumulate.
• The prediction ellipse uses vcov(fixed effects) + Psi and describes where the true operating point of the next study lies. It does not shrink with more studies, because Psi is real between-study spread rather than our ignorance.

⚠️ The prediction ellipse refers to the next study’s true operating point and excludes that study’s own sampling error. It is the two-dimensional analogue of a univariate prediction interval (see the heterogeneity page), not an interval for the point a next study would be observed at — predicting an observed point would require adding that study’s sampling variance, which depends on its sample size.

The two differ enormously here. The prediction ellipse covers 10.37 times the area of the confidence ellipse:

The interval around the summary sensitivity looks reassuringly tight, yet the next study’s true sensitivity could be as low as 0.291. That contrast is far stronger than the 0.19 percentage points between the two pooling methods, and it answers the question clinicians actually have: what will this test do in my clinic. A meta-analysis reporting only the confidence region is answering how precisely the average is known when the question was what happens next time.

The SROC curve: past the data it is extrapolation

The curve in that figure is the summary ROC (SROC) curve, here in the Rutter-Gatsonis parameterisation. It shows how sensitivity and specificity trade off if the cut-off is moved along the trajectory the model estimates. The area under it is AUC = 0.887, and restricted to the observed range the partial AUC is 0.861.

This AUC is not the area under a single study’s ROC curve (see ROC curves and the area under them). A single study’s ROC comes from moving the cut-off within one group of people; an SROC is across studies, each with a different population, different readers and possibly a different reference standard. They look alike and mean quite different things.

Only part of the curve has data behind it. These 14 studies have false positive rates between 0.059 and 0.520; inside that range the curve is drawn solid, outside it dotted. The dotted segments are extrapolation — their shape is determined entirely by the model’s assumptions, and no study contributed information there. When a paper draws its SROC curve out to a false positive rate of 0.01 and quotes the sensitivity at that point, that number was computed, not measured.

Two model names you will see in reports

The methods section of a diagnostic accuracy meta-analysis usually names one of the following two, and readers often take them for competing approaches:

• Bivariate random effects (Reitsma and colleagues) — a two-dimensional random-effects model on the pair (logit sensitivity, logit false positive rate), estimating the summary operating point and its joint uncertainty. This is what the page uses.
• HSROC (hierarchical summary ROC, Rutter and Gatsonis) — parameterised in terms of an accuracy component and a threshold component, estimating the curve.

With no covariates the two are the same model written two ways and can be converted into each other — which is why this page can both report a summary point (from the bivariate parameters) and draw an SROC curve (the ruttergatsonis parameterisation). They part company once covariates are added: putting a covariate on the accuracy component is not the same question as putting it on sensitivity and specificity.

So a paper saying bivariate while showing an SROC curve is not contradicting itself, and neither is one saying HSROC while reporting a summary point. The things worth asking are which level the covariates sit on, and whether the summary point was estimated directly or picked off the curve.

Why not just pool the DOR

Faced with an awkward pair, the natural move is to squash it into a scalar. The diagnostic odds ratio does that, and it is the last column of the first table on this page. Once it is a scalar, ordinary meta-analysis tooling works again — the early Moses-Littenberg approach regressed each study’s log DOR against a proxy for its threshold to trace an SROC curve.

There are two costs, and neither is small.

First, the DOR discards the distinction that matters clinically. Two tests with the same DOR can differ completely: one with high sensitivity and unremarkable specificity (good for ruling out), the other the reverse (good for ruling in). An equal DOR does not make them interchangeable, and where a test belongs in a diagnostic pathway is precisely what has to be decided.

Second, the DOR is unstable on extreme tables. The 2 zero-cell studies here end up with the two most extreme DORs in the set after correction, and with the two largest standard errors of log DOR among all 14 studies (both near 1.42, while the rest fall between 0.17 and 0.87). Across the set the DOR runs from 8.9 to 1388.9, two orders of magnitude, and the largest two are entirely an artefact of the zero cells.

So the DOR’s place now is descriptive rather than the target of pooling: useful in a table to show how much overall discrimination varies between studies, unsuitable as the headline result.

One point, or the whole curve

The stronger the threshold effect, the less a statement of the form pooled sensitivity is X means on its own — it averages across studies operating at different thresholds and lands on a point at which possibly no study operated. The between-study spread here is genuinely wide: Psi’s two standard deviations on the logit scale are 1.175 for sensitivity and 0.638 for the false positive rate, which is why the prediction ellipse is so much larger.

Three usual ways of handling this; which one applies is decided by the data, not by taste:

1. Report the whole SROC curve together with the summary point and the prediction region. This page’s approach, and the honest presentation when the cut-offs really do differ: the curve says roughly how good the test is, the prediction region says where a next study might land.
2. Pull out the studies using one common threshold and pool those separately, leaving the rest for a sensitivity analysis. Only then does the summary point mean something definite — performance at that threshold. This requires the thresholds to be obtainable.
3. Meta-regression on the cut-off. The most informative option and the most demanding; it runs into the gap described earlier.

A different mistake entirely: adding up the 2x2 tables

Everything above is about the model. There is a more basic error that remains common in practice: adding the four cells of all 14 tables and computing sensitivity and specificity from the totals.

The three approaches agree well enough on sensitivity. They do not agree on specificity: adding up and dividing gives 0.691, about 8.9 percentage points below both meta-analytic estimates.

The cause has nothing to do with thresholds and everything to do with the unit of analysis. Adding up treats each participant as the unit, so the study with 6,955 participants — 38% of the sample — very nearly decides the answer on its own, and its specificity is the lowest of all 14. In a meta-analysis the unit is the study, not the person; random-effects weighting pulls the weights towards each other precisely so that one study cannot drown out a dozen others (see fixed-effect and random-effects models).

This error and the separate-pooling one are two different things and should not be run together: one ignores the link between the two measures, the other ignores that studies are the unit of analysis. A paper can commit either without the other.

Continuity correction: hold it fixed or the comparison is not clean

The side-by-side table above rests on one condition: both approaches used the same continuity correction — 0.5 added to every cell, which is what mada::reitsma() does internally. Only then does the remaining difference reflect the model.

Switch the univariate side to metafor’s default instead (adding 0.5 only where a cell is zero) and the gap widens to 0.98 and 0.20 percentage points:

The extra difference is the correction, not the model. A paper comparing zero-cell-only univariate results against a bivariate fit and attributing the whole gap to separate pooling has attributed it wrongly. Zero-cell handling has no settled answer of its own, and the costs of the various choices are covered in effect measures and their variances; the requirement here is only this: when comparing, hold fixed everything that is not what you are comparing.

Run it yourself

library(mada)
library(metafor)
library(ellipse)

data(AuditC, package = "mada")     # 14 studies, columns TP / FN / FP / TN
cc <- AuditC + 0.5                 # 0.5 on every cell, as reitsma() does

# -- 1. look at the mechanism before fitting anything ------------------------
madad(AuditC, correction = 0.5)$cor_sens_fpr   # sens vs FPR (positive)
cor.test(qlogis(cc$TP / (cc$TP + cc$FN)),      # sens vs spec (negative)
         qlogis(cc$TN / (cc$TN + cc$FP)))

# -- 2. separate pooling: two independent univariate models ------------------
es_sens <- escalc(measure = "PLO", xi = cc$TP, ni = cc$TP + cc$FN)
es_spec <- escalc(measure = "PLO", xi = cc$TN, ni = cc$TN + cc$FP)
m_sens  <- rma(yi, vi, data = es_sens, method = "REML")
m_spec  <- rma(yi, vi, data = es_spec, method = "REML")
plogis(c(m_sens$beta, m_sens$ci.lb, m_sens$ci.ub))   # pooled sens and 95% CI
plogis(c(m_spec$beta, m_spec$ci.lb, m_spec$ci.ub))   # pooled spec and 95% CI

# -- 3. bivariate: one model, with a between-study covariance ----------------
fit <- reitsma(AuditC, correction = 0.5, correction.control = "all")
summary(fit)          # pooled sensitivity and false pos. rate with 95% CIs
fit$Psi               # between-study covariance; psi12 > 0 IS that correlation
AUC(fit)              # area under the SROC curve (not a single study's ROC)

# -- 4. two ellipses: mada draws the confidence one, add Psi for the other ---
b <- as.numeric(coef(fit))                 # (logit sens, logit fpr)
V <- fit$vcov                              # sampling covariance of the fixed effects
e_conf <- ellipse(V[c(2, 1), c(2, 1)],            centre = b[c(2, 1)], level = 0.95)
e_pred <- ellipse((V + fit$Psi)[c(2, 1), c(2, 1)], centre = b[c(2, 1)], level = 0.95)

plot(fit, sroclwd = 2, xlim = c(0, 0.62), ylim = c(0.25, 1))
points(fpr(AuditC), sens(AuditC), pch = 21)
lines(plogis(e_conf[, 1]), plogis(e_conf[, 2]), lwd = 2)
lines(plogis(e_pred[, 1]), plogis(e_pred[, 2]), lty = 2)   # the much larger one

# -- 5. what does that rectangle actually cover? simulating beats reading ----
set.seed(20260823)
z <- matrix(rnorm(2 * 2e5), ncol = 2) %*% chol(V)
sim_sens <- plogis(b[1] + z[, 1])
sim_spec <- 1 - plogis(b[2] + z[, 2])
mean(sim_sens >= plogis(m_sens$ci.lb) & sim_sens <= plogis(m_sens$ci.ub) &
     sim_spec >= plogis(m_spec$ci.lb) & sim_spec <= plogis(m_spec$ci.ub))

import numpy as np
import pandas as pd
from scipy.special import logit
from scipy.stats import pearsonr, spearmanr

# The per-study 2x2 counts of mada::AuditC, 14 studies.
tab = pd.DataFrame({
    "TP": [47, 126, 19, 36, 130, 84, 68, 752, 59, 142, 137, 57, 34, 152],
    "FN": [9, 51, 10, 3, 19, 2, 0, 0, 5, 50, 24, 3, 1, 51],
    "FP": [101, 272, 12, 78, 211, 68, 112, 3226, 55, 571, 107, 103, 21, 88],
    "TN": [738, 1543, 192, 276, 959, 89, 423, 2977, 136, 2788, 358, 437, 56, 264],
})

# TP / FN / FP / TN per study; transcribed just above
cc = tab + 0.5                      # the same continuity correction as in R
sens = cc.TP / (cc.TP + cc.FN)
spec = cc.TN / (cc.TN + cc.FP)

print(pearsonr(logit(sens), logit(spec)))     # negative = threshold effect
print(spearmanr(logit(sens), logit(spec)))

# Drop each study in turn: is the correlation carried by one leverage point?
r_loo = [np.corrcoef(np.delete(logit(sens).values, i),
                     np.delete(logit(spec).values, i))[0, 1]
         for i in range(len(cc))]
print(min(r_loo), max(r_loo))

# And that is as far as Python goes: no off-the-shelf bivariate / HSROC fit.
# Do not quietly settle for two univariate models instead -- that is the very
# thing this page is about.

Reading the report

1. Whether the model is named. Bivariate random effects (Reitsma) or HSROC (Rutter-Gatsonis) are acceptable answers; pooled sensitivity and specificity with no model named usually means two univariate models.
2. Whether a between-study covariance or its correlation is given. That is where the threshold effect appears in a report. Two τ² values and nothing else means the model treated the measures as independent.
3. Whether the two intervals are used crossed. An abstract giving sensitivity X (CI) and specificity Y (CI), with the body making a joint claim from them, is the rectangle from the section above.
4. Whether the region reported is a confidence or a prediction region. A meta-analysis with only the former cannot answer what the test will do in your population.
5. How far the SROC curve is drawn. If the observed range of false positive rates is not marked, treat the ends of the curve as extrapolation.
6. How zero cells and the continuity correction were handled. This matters most when two approaches are being compared: unless the correction is the same, the attribution of the difference is wrong.
7. Whether an add-up-and-divide result is being presented as a meta-analysis. Easy to check: if the pooled value equals total TP divided by total (TP + FN), that is what it is.
8. Whether the model converged, and the residual degrees of freedom. A bivariate fit estimates five parameters (two means, two variances, one covariance) and can converge poorly, or push the covariance to a boundary, when studies are few. This fit converged, with 23 residual degrees of freedom and a log-likelihood of 31.56. When convergence is not mentioned, at least check that Psi’s correlation is not exactly plus or minus one — usually the sign of a boundary estimate.

Common misuses

Reproducing every number on this page

/opt/homebrew/bin/Rscript figures/scripts/B7-08-diagnostic-meta.R

Report a problem on this page