Theoretical relationship between epidemic growth rates and Ct value statistics
To understand how biological and practical factors might affect Ct-based nowcasting performance, we first examined their impact on the theoretical relationship between Ct value statistics and epidemic growth rates. We implemented analytical convolutions of population-level infection incidence from a Susceptible-Infected-Recovered model assuming R0 = 1.5 with various within-host viral kinetics models to generate mean Ct values observed through random or symptom-driven testing over an epidemic (Supplementary Fig. 1; Supplementary Methods).
Figure 1 demonstrates the impact of viral kinetics properties (peak viral load, time from infection to peak viral load, individual-level variance) and sampling approaches (distribution of delays between infection and sampling) on the predictability of growth rates using mean Ct values. These results demonstrate several key ideas. First, in most scenarios, there is a monotonic relationship between mean Ct value and growth rate, providing a basis for predicting the latter from the former. This monotonic relationship arises because the epidemic growth rate determines the time-since-infection distribution in cross-sectional samples, which in turn determines the Ct value distribution. For instance, over the course of an outbreak with R0 = 1.5, the mean time-since-infection ranges from 8 days to 25 days (assuming a maximum time-since-infection of 35 days), which sits entirely within the viral decline phase; as such the population-mean Ct value monotonically increases with decreasing growth rate. Conversely, if the time-since-infection range sat entirely within the viral growth phase (which could arise either from a much higher epidemic growth rate or a longer viral growth phase), the relationship would be reversed. Second, the precision of growth rates estimated from Ct values depends on the shape and gradient of this relationship. Because there is uncertainty in the true population mean Ct value from a finite number of observations, a steeper gradient will correspond to a wider range of compatible growth rates for a given set of measured Ct values. Third, biological and logistical factors influence the shape of the relationship and hence prediction accuracy, but in different ways. For example, increasing the duration of the viral growth phase decreases the precision of the estimated growth rate; a sufficiently long growth phase (producing a fully symmetric viral kinetics curve) results in a non-monotonic relationship and thus bi-modal estimates (Fig. 1, second row, high scenario).
A Varying viral kinetics and sampling delay distribution assumptions. Top two rows show assumed mean Ct value over time-since-infection. Third row shows multiple assumed 95% quantiles of observed Ct values. Bottom row shows varying assumed distributions of delays between symptom onset and testing. Scenarios are colored to match (B). B Curved colored lines (blue: low scenario, orange: medium, magenta: high) show the analytical relationship between mean Ct value and daily growth rate of new infections during a single wave epidemic, observed through random testing (first column) or following symptom onset (second column). The horizontal, colored lines show the mean and 95% confidence intervals for the population mean Ct value, assuming the true mean corresponds to a growth rate of 0.02 (horizontal dashed line), a standard deviation of 3 (varied between 1, 3 and 5 in the bottom row), and sample size of 500. The vertical, colored lines show the range of epidemic growth rate values corresponding to the mean and 95% CI of the mean Ct value. Horizontal dashed lines indicate the true epidemic growth rate.
Overall, this analysis demonstrates that although a relationship between mean Ct value and epidemic growth rate is always present, the strength of this relationship, and thus the statistical power of a Ct-based inference approach, varies depending on the underlying time-since-infection distribution, viral kinetics, and sampling approach.
Modeling and nowcasting epidemic growth rates
We next examined how the Ct value-growth rate relationship can be used for epidemic nowcasting, applying the same set of basic models across several synthetic and real-world datasets. We used generalized additive models (GAM) with smoothing splines to predict growth rates of cases as a non-linear function of daily mean and skewness of observed Ct values (see Supplementary Methods). We also fit GAMs to predict the epidemic direction (i.e., whether incidence is growing or declining) using a logit link function and with the same explanatory variables and smoothing splines.
With each dataset, we assessed in-sample fits of model-predicted vs. observed growth rates and directions across the entire dataset, based on root mean squared error (RMSE) for predicted growth rate and the area under the receiver operating characteristic curve (AUC) for predicted epidemic direction. For all datasets, model residuals were normally distributed around zero (example shown in Supplementary Fig. 2). We then refit the models using separate training and testing subsets of the data. To approximate a realistic application of the Ct value-based approach in an ongoing epidemic, we fit the models using only the first 16 weeks of data and then performed rolling nowcasts with a two-week time horizon, using the fitted model to estimate the epidemic growth rate and direction daily over the next two weeks based on the Ct values reported during that time. At the end of each two-week window, we re-fit the model using all Ct values and incidence data up to that time point, then nowcast the next two-week window, and so on. As a sensitivity analysis, we also compared RMSE and AUC with a fixed train-test split date on 31 Dec 2021 (Supplementary Table 4).
Nowcasting performance on synthetic datasets
We first applied these GAMs to several synthetic Ct value datasets, generated as numerical implementations of the analytical convolutions presented in the previous section. These synthetic datasets combined a range of viral kinetics models and sampling regimes with population-level reported SARS-CoV-2 incidence curves for Massachusetts, USA to produce linelists of individual-level test results and Ct value observations (see Supplementary Methods; Supplementary Fig. 3 & 4). Using synthetic datasets in this way allowed us to incorporate or exclude the effects of confounding factors on observed population-level Ct value distributions, in addition to the effect of the epidemic trajectory itself (see Supplementary Table 1). The simulations generated approximately 22,000 positive tests under asymptomatic screening (mean of 20 observations per day, range 1–250) and 250,000 positive tests under symptom-driven testing (mean of 230 observations per day, range of 3–2700), over 1000 days.
We simulated two main synthetic datasets: (1) an ideal scenario where a strong Ct value/growth rate relationship is expected, assuming (a) highly asymmetric viral kinetics, with a very short growth phase and long clearance phase; (b) low variation in observed viral load/Ct value for a given time-since-infection; and (c) a uniform probability of testing an individual at any time up to 7 days post symptom onset or tested once at random during the simulation regardless of infection status, and (2) a realistic baseline scenario assuming (a) increased symmetry in viral kinetics, with a growth and clearance phase duration based on ancestral SARS-CoV-2 (Supplementary Fig. 5); (b) moderate variation in observed viral load/Ct value for a given time-since-infection; and (c) individuals are tested after a low-variance, gamma-distributed delay following symptom onset or once at random during the simulation (Supplementary Fig. 6). We also performed sensitivity analyses on datasets with variations of these three factors (Supplementary Fig. 6), as well as a dataset incorporating viral kinetics across epidemic waves to represent differences in viral variants (Supplementary Fig. 7). The standardized mean Ct value within each simulation and testing strategy (subtracting the overall mean Ct and dividing by the standard deviation) showed the same trend over time for all simulations, excepting the scenario assuming a symmetric viral kinetics model with equal growth and clearance phase durations (Supplementary Figs. 8 & 9).
Both ideal and realistic synthetic datasets had a negative correlation between the 7-day rolling average epidemic growth rate of cases and 7-day rolling average mean Ct value (Fig. 2A; Supplementary Fig. 4B). Ct values observed through symptom-based testing were typically lower and exhibited less variation than those observed through random testing. These results suggest that a strong relationship between mean Ct value and epidemic growth rates should be identifiable, at least in theory.
A Scatterplot showing 7-day rolling mean Ct value against 7-day rolling mean incidence growth rate. Inset text shows Pearson correlation coefficient R with two-tailed p-value. Blue line shows fitted linear regression line with 95% CI. B Model-predicted mean (black) vs. observed (blue) log incidence growth rates for the two main synthetic datasets, with model-predicted 95% confidence intervals (dark shading) and 95% prediction intervals (light shading). Predictions are from 2-week rolling nowcasts, concatenated into a single time series. Vertical dashed line denotes the end of the initial training period.
With the ideal synthetic dataset, the estimates from the GAMs closely tracked observed growth rates using Ct value means and skew, with an in-sample RMSE of 0.0184—approximately 10% of the range in observed log incidence growth rates (Supplementary Figs. 10 & 11; Table 1)—as well as accurately predicting epidemic direction (in-sample AUC = 0.926). Nowcast accuracy over the two-week windows was slightly worse than the in-sample predictive performance (mean across all nowcast windows, RMSE = 0.0192, AUC = 0.910) but still accurately tracked the epidemic over the full time period (Fig. 2B). In comparison, a null model using the mean growth rate over the preceding two weeks as the prediction for the following two-week window gave an RMSE of 0.0345 and AUC of 0.758.
Model predictive performance was worse when using the realistic baseline synthetic dataset, with an in-sample RMSE of 0.0331 and AUC of 0.748 (Supplementary Fig. 10 & 11; Table 1). Of the three factors examined individually, moving to a symmetric viral kinetics model with a longer viral growth phase had the largest impact on model performance, resulting in the greatest increase in RMSE and decrease in AUC (Supplementary Fig. 10–12, Table 1). When these models were applied to nowcasting growth rates in two-week increments, the greatest performance reduction again occurred when using either a realistic or symmetric viral kinetics model (Table 1). Nowcasting predictive performance for the realistic baseline scenario was similar to the null model using the mean growth rate of the preceding two weeks as the prediction for the following two-week window. In the sensitivity analysis assuming different variant waves, nowcasting performance was poor, as unaccounted for changes to the Ct/growth rate relationship due to changes to within-host viral kinetics resulted in biased predictions (Supplementary Fig. 13).
Prediction accuracy is increased with higher sample sizes in the training set
Despite the clear theoretical relationship between mean Ct value and the epidemic growth rate shown in Fig. 1 and Supplementary Fig. 4, prediction accuracy using synthetic Ct values under realistic assumptions was relatively poor. To understand why, we further examined the relationship between sample size and model predictive performance, evaluating outputs from a range of positive tests each day under the realistic baseline scenario collected under either random testing or symptom-based testing. We compared (a) in-sample predictive performance, fitting the GAM to the 7-day rolling average of the mean and skew of Ct values for each sample size and (b) nowcasting predictive performance, varying both the number of samples used to train the model and subsequently collected each day for prediction during the nowcasting windows.
These simulations showed three patterns. First, in-sample predictive performance (assessed by RMSE) improved substantially as the number of daily positive samples increased from 25 to 1000 (Fig. 3), reaching similar RMSE values to the ideal synthetic scenario (Table 1) with 100 positive samples collected at random per day or 500 samples through symptom-based testing. Second, both in-sample and nowcasting predictive performance was greater for the same sample size using Ct values collected at random compared to symptom-based testing (Fig. 3, Supplementary Figs. 14 & 15). The lower sample size requirement under random testing reflects the shallower gradient between mean Ct value and growth rate compared to symptom-based testing, as described in the section Theoretical relationship between epidemic growth rates and Ct value statistics. Third, increasing the number of positive samples available to predict growth rates did not compensate for having too few samples available to train the model in the first place (Supplementary Figs. 14 & 15). Model predictions were only accurate if the relationship between growth rates and Ct values was accurately characterized to begin with. These findings likely explain the limited performance under the realistic baseline scenario, where sample sizes were around 20 positive tests per day under random testing and around 250 for symptom-based testing, but with wide variation depending on incidence (approximately 10x higher at peak incidence).
Real-world relationship between observed Ct-value statistics and epidemic trajectories
Having established a baseline for model nowcasting performance using the synthetic data, we next tested the nowcasting models on two real RT-qPCR datasets: (1) routine hospital testing data from the Mass General Brigham hospital system in eastern Massachusetts (MGB), spanning Mar 2020-Feb 2023, and (2) municipal testing data from Los Angeles County, California (LAC), spanning May 2020-Jul 2021 and Jan-Sep 2022 (see Supplementary Methods). The MGB data came largely from mandatory screening testing of outpatient, inpatients and emergency room admissions, while the LAC data were primarily symptom-driven voluntary testing (Supplementary Table 3). Both datasets contained specimen collection dates and Ct values for SARS-CoV-2 positive results; LAC data also included vaccination status, symptom status, and symptom onset dates. MGB Ct values came from seven platform/assay combinations, while Ct values from LAC data came from one PCR platform with two possible assays (Supplementary Table 3).
We limited our analysis to tests reporting Ct values, using the first available recorded Ct value for each infection episode. The final analyzed sample included 104,534 (MGB) and 279,492 (LAC) Ct values. We compared these Ct values against reported COVID-19 incidence for Massachusetts and Los Angeles County, respectively. We segment the data into four variant eras based on the SARS-CoV-2 variant known or believed to be dominant in the U.S. during different approximate time periods, to allow for differences in viral kinetics by variant.
Ct value distributions from both MGB and LAC datasets showed substantial variation over the course of the pandemic (Fig. 4A& Supplementary Fig. 16A). Reported COVID-19 incidence in both locations varied over time, as well (Fig. 4B& Supplementary Fig. 16B), with large infection waves in the winters of 2020-21 and 2021-22, though the pattern of incidence was not synchronized across both settings. While absolute incidence varied widely, incidence growth rates remained largely between ±0.2 throughout the course of the pandemic (Fig. 4B& Supplementary Fig. 16B).
A Weekly Ct value quantiles over time, showing weekly median Ct value (black line) and 50/80/90/95% quantiles (grey shading, darkest to lightest). B 7-day rolling average reported incidence (grey bars), growth rate in 7-day rolling average reported incidence (grey line), and smoothed growth rate (blue line). Background is shaded by time periods of different variant dominance (grey: Wild Type, yellow: Alpha, blue: Delta; red: Omicron). Vertical dashed line demarcates the test-train split. C Incidence growth rate compared to the smoothed daily mean and the skewness of Ct value distributions. Colored lines and shaded grey regions show fitted cubic spline GAMs with 95% confidence intervals, stratified by period of variant dominance.
We found the mean and skewness of observed Ct value distributions (calculated daily over a seven-day moving window and excluding days with fewer than 10 Ct values reported) negatively correlated with the growth rate in reported incidence (Fig. 4C& Supplementary Fig. 16C). Analysis of cross-correlation functions found Ct value distributions lagged incidence growth rate in the MGB data, with strongest correlations at around 19-days lag (autocorrelation function, ACF = −0.462), and led incidence growth rates for the LAC data, with strongest correlations at around 9-days lead (ACF = −0.062) (Supplementary Figs. 17 & 18). However, for real-time nowcasting, we focused on the relationship between same-day Ct values and incidence (i.e., lag = 0 days; Fig. 4C & Supplementary Fig. 16C), which still showed significant correlation. Higher incidence growth rates corresponded with lower same-day average Ct values (Spearman’s correlation coefficient: MGB Rho = −0.43, LAC Rho = −0.22) and with positively skewed Ct distributions (MGB Rho = 0.35, LAC Rho = 0.43).
Nowcasting epidemic growth rates using Ct values in Massachusetts, USA, and Los Angeles County, USA
We next re-trained the same GAM models used with synthetic data to the MGB and LAC dataset using smooth functions of mean and skewness of Ct values to predict log incidence growth rates, with corresponding logistic models to predict epidemic direction. Model predictions were compared against observed values first in-sample across the entire dataset, over a rolling two-week prediction window with an expanding training period encompassing all preceding weeks, and with a single fixed train-test split date on 31 Dec 2021.
In both datasets, this simple model achieved in-sample prediction accuracy for incidence growth rate slightly worse than performance on the realistic synthetic data, with relatively small absolute errors (MGB RMSE = 0.0451; LAC RMSE = 0.0335, see Table 2, Supplementary Fig. 19–S22, and Supplementary Table 4). Corresponding logistic regression models successfully discriminated growing from declining incidence (Area under the curve: MGB AUC = 0.785, LAC AUC = 0.843).
The models were able to nowcast growth rates, in two-week increments with models periodically refitted to more recent data, with accuracy slightly worse than in-sample model fits (MGB RMSE = 0.0521, LAC RMSE = 0.0386) (Fig. 5A & Supplementary Fig. 23A). This level of nowcast accuracy was likewise slightly worse than nowcasting performance with realistic synthetic data. While average prediction error was relatively small, comparable to in-sample model error and to prediction error with realistic synthetic data, accuracy was highly variable from one two-week window to the next (Fig. 5B & Supplementary Fig. 23B). Nowcast accuracy was comparable to model performance over a fixed multi-month prediction window, slightly better for one dataset and worse for the other (MGB RMSE = 0.047, LAC RMSE = 0.0458; Supplementary Table 4). Nowcast predictions of epidemic direction were slightly worse than in-sample ones (MGB AUC = 0.718, LAC AUC = 0.799) and outperformed the directional predictions for realistic synthetic data. In addition, over all two-week nowcast windows combined, model-predicted growth rates correlated moderately well with observed ones (Spearman’s Rho: MGB Rho = 0.387, LAC Rho = 0.576).
A Model-predicted mean (black) vs. observed (blue) log incidence growth rates for MGB data, with 95% confidence intervals (dark shading) and 95% prediction intervals (light shading). Background is shaded by time periods of different variant dominance (grey: Wild Type, yellow: Alpha, blue: Delta; red: Omicron). Vertical dashed line demarcates the end of the initial training period. B RMSE of predicted vs. observed log incidence growth rates for each 2-week nowcasting window. “Inflection periods” refer to times when the absolute change in smoothed log incidence growth rate exceeded 0.025 over a one-week period, marked with points above each subplot.
Nowcasting state- or county-level epidemic growth rates combining Ct values and the number of new positive tests
To understand what information Ct values add over simply counting new positive tests, we evaluated GAM and logistic models which included the smoothed log incidence growth rate of positive tests within each dataset as a predictor. In all cases, we aimed to predict state- or county-wide incidence trends using data from a single reporting system (i.e., MGB hospital testing or LAC community testing). These within-dataset growth rates would be expected to correlate strongly with population-wide growth rates, provided they are a representative sample. We compared model predictions using only Ct value statistics, only growth rate of positive tests, or combined Ct value statistics and positive tests.
For both the MGB and LAC datasets, the models incorporating both positive tests and Ct value mean and skew had the greatest prediction accuracy (MGB RMSE = 0.0395, LAC RMSE = 0.0366), close to the accuracy achieved with realistic synthetic data, indicating that both the number of positive tests and Ct value distributions were informative about incidence growth rates. However, for the MGB dataset, the model using the growth rate of positive tests alone performed substantially better than the model using Ct mean and skew alone (RMSE 0.0399 vs. 0.0521). In contrast, for the LAC dataset, the model using only Ct mean and skew outperformed the model using positive test data only (RMSE 0.0386 vs. 0.0418). The same rank ordering of nowcast performance holds for epidemic direction prediction accuracy as measured by AUC.
Sensitivity analyses
We performed several sensitivity analyses, assessing nowcasting performance under various sample restrictions and stratifications. First, we restricted model performance assessment to periods of rapid change in the epidemic trajectory (defined as days when the absolute change in smoothed incidence growth rate exceeded 0.025 over a one-week period). Across both datasets, prediction error over periods of rapid change was greater than over the whole nowcast period (MGB RMSE = 0.0644 [rapid change] vs. 0.0521 [nowcast], LAC RMSE = 0.0462 [rapid change] vs. 0.0386 [nowcast]; see Table 2), while directional prediction accuracy was comparable (MGB AUC = 0.7 vs. 0.718, LAC AUC = 0.779 vs. 0.799).
Next, we assessed model performance on the same datasets stratified or restricted by covariates including (1) test setting (i.e., outpatient vs. inpatient vs. emergency room), (2) symptom status, and (3) immune history (accounting for both vaccination and past infection) (Supplementary Fig. 24). The relationship between Ct values and growth rate sometimes differed when subsetting or stratifying by these variables (Supplementary Fig. 25), but including these stratifications in the model did not always improve predictive performance. Restricting to outpatient tests only improved prediction error compared to baseline (nowcast RMSE = 0.0491 vs. 0.0521 base), whereas incorporating symptom status or immune history slightly worsened prediction error (nowcast RMSE = 0.0454 for symptom-stratified, 0.0415 for asymptomatic/no symptom status only, 0.0394 for immunologically naïve only, vs. 0.0386 base, see Supplementary Table 6).
Nowcasting performance with variable sample size and outlier removal
Finally, we assessed the sensitivity of nowcasting performance to sample size both by randomly downsampling the MGB dataset (100 random draws) and by analyzing a third, smaller dataset from Tufts Medical Center—a hospital from the same region of Massachusetts—using the same response variable (i.e., log incidence growth rates for Massachusetts) but with approximately 10% of the total sample size of the MGB data (see Supplementary Methods; Supplementary Figs. 26 & 27). In most cases, prediction accuracy for incidence growth rate was comparable with the downsampled datasets and the equivalent full datasets (Fig. 6; see also Supplementary Table 5). Only with 10% of the full dataset (but not with the Tufts dataset) did nowcasting accuracy degrade appreciably, though with a lot of variation between draws; with 50–75% downsampling or a daily maximum of 25 positive samples, accuracy improved slightly compared to baseline. Likewise, directional prediction accuracy was generally similar between downsampled and full datasets, with substantially worse accuracy only for the 10% and 25% downsamples. Improved accuracy may reflect reduced influence of outliers – downsampling the full dataset tends to exclude the days with the smallest sample sizes, which are otherwise given equal weight in model training to days with more observations, while sub-sampling each day’s observations reduces the impact of outliers on each day’s observed Ct value distribution. To test this, we examined model performance with trimming of outlier Ct values from each day’s observed data. Trimming outliers reduced prediction error with 2.5%, 5%, and 10% trims (Fig. 6, Supplementary Table 5), while 2.5% and 5% trims also improved directional prediction accuracy.
Baseline comparison metrics (red) are re-calculated for only the days included in each downsampled dataset’s nowcast. For proportional and daily max downsampling, both downsampled and baseline performance are averaged over 100 random draws (and their corresponding days included). Trim percentages indicate quantiles trimmed from each end of daily Ct value distributions (i.e., 5% trim yields the 5-95 percentile range of Ct values). Grey bars indicate the mean number of days of data included for the various downsampling analyses; the random downsampling can reduce some days below the 10-sample inclusion threshold, hence the variation.
We also repeated these analyses on the models using within-dataset positive test counts instead of, or in combination with, Ct value distributions. In both instances, but particularly for the models using only within-dataset positive tests, reducing sample size substantially worsened model predictive performance for both growth rate and epidemic direction (Fig. S28 & 29). Greater reductions resulted in greater performance degradation, indicating that models relying on the growth rates of case counts in available data are sensitive to sample size. Similarly, the models relying on only case counts within the Tufts dataset performed substantially worse than models using the full dataset.





