Estimating epidemiologic dynamics from cross-sectional viral load distributions

Abstract

Estimating an epidemic's trajectory is crucial for developing public health responses to infectious diseases, but case data used for such estimation are confounded by variable testing practices. We show that the population distribution of viral loads observed under random or symptom-based surveillance-in the form of cycle threshold (Ct) values obtained from reverse transcription quantitative polymerase chain reaction testing-changes during an epidemic. Thus, Ct values from even limited numbers of random samples can provide improved estimates of an epidemic's trajectory. Combining data from multiple such samples improves the precision and robustness of this estimation. We apply our methods to Ct values from surveillance conducted during the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) pandemic in a variety of settings and offer alternative approaches for real-time estimates of epidemic trajectories for outbreak management and response.

Ct values reflect the epidemic trajectory and can be used to estimate incidence. (A and B) Whether an epidemic has rising or falling incidence will be reflected in the distribution of times since infection (A), which in turn affects the distribution of Ct values in a surveillance sample (B). (C) These values can be used to assess whether the epidemic is rising or falling and estimate the incidence curve.

Fig. 1.. The Ct value distribution reflects epidemiological dynamics over the course of an outbreak.

(A) Per capita daily incidence (histogram) and daily growth rate (blue line) of new infections in a simulated epidemic using an SEIR model. (B) Median days since infection versus daily growth rate of new infections by epidemic day. Labeled points here, and in (E) to (G), show five time points in the simulated epidemic. (C) Observed Ct value by day for 500 randomly sampled infected individuals. (D) Viral kinetics model (increasing Ct value after peak and subsequent plateau near the limit of detection), demonstrating the time course of Ct values (x axis; line shows mean, and ribbon shows 95% quantile range) against days since infection (y axis). Note that the y axis is arranged to align with (E). (E) Distribution of days since infection (violin plots and histograms) for randomly selected individuals over the course of the epidemic. Median and first and third quartiles are shown as green lines and points, respectively. (F) Skewness of observed Ct value distribution versus daily growth rate of new infections by epidemic day. (G) Distribution of observed Ct values (violin plots and histograms) among sampled infected individuals by epidemic day. Median and first and third quartile are shown as purple lines and points, respectively. (H) Time-varying effective reproductive number, R_t, derived from the SEIR simulation, plotted against median and skewness of observed Ct value distribution.

Fig. 2.. Single cross-sectional distributions of observed Ct values can be used to reconstruct epidemic trajectories in a Massachusetts long-term care facility.

(A) Estimated prevalence [faint orange lines show posterior samples, solid orange line shows posterior median, and orange ribbon shows 95% credible intervals (CrIs)] and incidence (red line shows posterior median and red ribbon shows 95% CrI) from the standard compartmental (SEEIRR) model fit to point prevalence at three sampling times (error bars show 95% binomial confidence intervals). (B) Model-predicted Ct distributions (blue) fitted to the observed Ct values (gray bars) from each of three cross-sectional samples. Shown are the posterior median (black line) and 95% CrI for the expected Ct distribution (dark blue ribbon) and 95% prediction intervals based on simulated observations (light blue ribbon). Note that prediction intervals are much wider than CrIs because they result from simulating observations with a small sample size. (C) Each panel shows results from fitting the Ct-based SEIR model separately to three cross sections of virologic data. Shown are random posterior samples (red lines) and the maximum posterior probability (MAP) trajectory (black line) for the incidence curve. (D) Fitted median (blue point) and 95% CrI (blue error bars) for the proportion of samples testing positive based on the Ct model compared with the observed proportion tested positive (gray cross). (E) Thirty-five–day (green) and 1-day (pink) average growth rates from the Ct model estimates in (C) at three time points (violin plots) compared with growth-rate estimates from the SEEIRR model in (A) (lines and shaded ribbons).

Fig. 3.. Inferring epidemic trajectory from cross-sectional surveillance samples with observed Ct values yields nearly unbiased estimates of the time-varying effective reproductive number, R_t, whereas changing testing rates lead to biased estimation using reported case counts.

(A) Number of positive tests per day by sampling time in epidemic and testing scheme for reported case counts (top row) and surveillance Ct sampling (bottom row), from a simulated SEIR epidemic. Analysis times corresponding to (B) are shown by the dashed vertical lines. (B) R_t estimates from 100 simulations for each epidemic sampling time, testing scheme, and estimation method. Each point is the posterior median from a single simulation. R_t estimates for reported case counts use EpiNow2 estimation and for surveillance Ct samples use the Ct-based likelihood for one or multiple cross sections fitted to an SEIR model. The semitransparent points at the right of the plots are those surveillance samples fit to an SEIR model using only a binary result of testing, assuming PCR positivity reflects the infectious compartment. True model-based R_t on the sampling day is indicated by the black star and dashed horizontal line, whereas an R_t of 1, indicating a flat outbreak, is indicated by the solid horizontal line.

Fig. 4.. Cross-sectional distributions of observed Ct values can estimate the complex statewide epidemic trajectory from hospital-based surveillance at Brigham and Women’s Hospital in Massachusetts.

(A) Daily confirmed new cases in Massachusetts (gray bars) and estimated time-varying effective reproductive number, R_t. (B) Estimated R_t from the case counts versus median and skewness of observed Ct value distribution by weekly sampling times. (C) Distribution (violin plots and points) and smoothed median (blue line) of observed Ct values by sampling week. Red box highlights data used to inform estimates in (D). (D) Posterior median (yellow arrow) and distribution (blue shaded area) of estimated daily growth rate of incident infections from an SEIR model fit to a single cross section of observed Ct value data from the week commencing 14 June 2020. Shading density is proportional to posterior density. Fits to all single weekly cross sections are shown in fig. S14. (E) Posterior distribution of relative probability of infection by date from a GP model fit to all observed Ct values (ribbons show 95% and 50% CrIs; line shows posterior median). Note that the y axis shows relative rather than absolute probability of infection, as the underlying incidence curve must sum to one: Only positive samples were included in the estimation, and all samples were therefore assumed to have been from infections. (F) Comparison of estimated daily growth rate of incident infections from the GP model (blue line and shaded ribbons show posterior median and 95% CrI) to that from R_t estimation using observed case counts (red and green line and shaded ribbons show posterior median and 95% CrI) by date. Note that estimates of infection incidence are made for dates before the first observed sample date of 15 April 2020, as far back as 1 March 2020, but the x axis is truncated at 1 April 2020 (fig. S19).