A statistical model of COVID-19 testing in populations: effects of sampling bias and testing errors

1. Introduction

Real-time estimation of the level of infection in a population is important for assessing the severity of an epidemic as well as for guiding mitigation strategies. Several previous studies have addressed the issue of correcting for errors and testing biases. However, inferring disease prevalence via patient testing is challenging due to testing inaccuracies, testing biases and heterogeneous and dynamically evolving populations and severity of the disease.

There are two major classes of tests that are used to detect previous and current SARS-CoV-2 infections [1]. Serological, or antibody, tests measure the concentration of antibodies in infected and recovered individuals. Since antibodies are generated as a part of the adaptive immune system response, it takes time for detectable antibody concentrations to develop. Serological tests should thus not be used as the only method to detect acute SARS-CoV-2 infections. An alternative testing method is provided by viral-load or antigen tests, such as reverse transcription polymerase chain reaction (RT-PCR), enzyme-linked immunosorbent assay (ELISA) and rapid antigen tests, which are able to identify ongoing SARS-CoV-2 infections by directly detecting SARS-CoV-2 nucleic acid or antigen.

Test results are mainly reported as binary values (0 or 1, negative or positive) and often do not include further information such as the cycle threshold ($\mathrm{Ct}$) for RT-PCR tests. The cycle threshold $\mathrm{Ct}$ defines the minimum number of PCR cycles at which amplified viral RNA becomes detectable. Large values of $\mathrm{Ct}$ indicate low viral loads in the specimen. An increase in $\mathrm{Ct}$ by a factor of about $3.3$ corresponds to a viral load that is about one order of magnitude lower [2]. Cycle threshold cutoffs are not standardized across jurisdictions and range from values between 37 and 40, making it difficult to compare RT-PCR test results [3]. Lower $\mathrm{Ct}$ cutoffs in the range of 30–35 may be more reasonable to avoid classifying individuals with insignificant viral loads as positive [3].

Further uncertainty in COVID-19 test results arises from different type I errors (false positives) and type II errors (false negatives) that are associated with different assays. Note that inherent to any test, the threshold (such as $\mathrm{Ct}$ mentioned above) may be tunable. Therefore, besides intrinsic physical limitations, binary classification of ‘continuous-valued’ readouts (e.g. viral load) may also lead to an overall error of either type [4]. In this work, we will assume that there is a standardized threshold and the test readout is binary; if any virus is detected, the test subject is positive. We will not explicitly model the underlying statistics of the errors but assume that the test readouts are binary but can be erroneous at specified rates. Some uninfected individuals will be wrongly classified as infected with rate $\mathrm{FPR}$ and some infected individuals will be wrongly classified as uninfected with rate $\mathrm{FNR}$. For serological COVID-19 tests, the estimated proportions of false positives and false negatives are relatively low, with $\mathrm{FPR}\approx 0.02-0.07$ and $\mathrm{FNR}\approx 0.02-0.16$ [5–8]. The $\mathrm{FNR}$s of RT-PCR tests depend strongly on the actual assay method [9,10] and may be significantly larger than those of serological tests. Typical values of $\mathrm{FNR}$ for RT-PCR tests lie between 0.1 and 0.3 [11,12] but might be as high as $\mathrm{FNR}\approx 0.68$ if throat swabs are used [7,12]. False-negative rates may also vary significantly depending on the time delay between initial infection and testing [8]. According to a systematic review [13] that was conducted worldwide, the initial value of $\mathrm{FNR}$ is about 0.54, underlying the importance of retesting. Similar to serological tests, reported false-positive rates of RT-PCR tests are about $\mathrm{FPR}=0.05$ [7].

Estimates of disease prevalence and other surveillance metrics [14,15] need to account for $\mathrm{FPR}$s and $\mathrm{FNR}$s, in particular if reported positive-testing rates [16] are in the few percent range and potentially dominated by type I errors. In addition to type I/II testing errors, another confounding effect is biased testing [17], that is preferential testing of individuals that are expected to carry a high viral load (e.g., symptomatic and hospitalized individuals). Biasing testing towards certain demographic and risk groups leads to additional errors in disease prevalence estimates that need to be corrected for.

In §2, we discuss related studies that developed statistical methods to correct for erroneous and biased testing. To account for type I/II errors, bias, retesting and exclusion after testing, we develop a corresponding framework for disease testing in §3. We apply our testing model to COVID-19 testing and case data in §4 and estimate testing bias by comparing random-sampling testing data [18] with officially reported, biased COVID-19 case data in §5. We conclude our study in §6.

2. Related work

Several previous studies have addressed the issue of correcting for errors and testing biases. In the random-sampling study [18], specificity and sensitivity corrected ELISA results are reported without specifying the actual statistical correction method. In another work [19], corrected case numbers for different European countries are derived based on the assumption that the infection fatality ratio ($\text{IFR}$) is independent of the geographical location. If the $\text{IFR}$ were known exactly, this method could be used to estimate the sampling bias by comparing the reported number of cases with the corresponding reported number of deaths divided by $\text{IFR}$. However, the framework in [19] does not account for false negatives and false positives. In addition, there are geographical variations of the $\text{IFR}$ that may be attributed to significant differences in incidence rates, population density, preparedness of public health systems and age structure [14,20,32]. Therefore, the assumption of a time and location-independent $\text{IFR}$ may yield inaccurate results.

In [21], a semi-Bayesian probabilistic bias analysis is used to estimate the cumulative number of SARS-CoV-2 infections in the United States. The employed corrections for erroneous testing are similar to the results that we derive in §4. Corrections for incomplete testing are based on distributions associated with random sampling studies similar to [18], which we use in §5.

One major difference between [21] and our work is that we derive the distributions with and without retesting, and explicitly account for test-type-dependent specificities and sensitivities.

3. Statistical testing model

Here, and in the following subsections, we develop a general statistical model for estimating the number of infected individuals in a jurisdiction by testing a sample population. The relevant variables and parameters to be used in our derivations are listed and defined in table 1. Suppose we randomly administer $Q$ tests within a given short time period (e.g. within 1 day or 1 week) to a total effective population of $N$ previously untested individuals. This population comprises $S$ susceptible, $I$ infected and $R$ removed (i.e. recovered or deceased) individuals, which are unknown. $S$, $I$ and $R$ can dynamically change from one testing period to another due to transmission and recovery dynamics, as well as removal from the untested pool by virtue of being tested. The total population $N=S+I+R$ can also change through intrinsic population dynamics (birth, death and immigration), but can assumed to be constant over the typical time scale of an epidemic that does not cause mass death.

We start the derivation of our statistical model by first fixing $S$, $I$ and $R$, assuming both perfect error-free testing, considering a ‘testing with replacement’ scenario, in which tested individuals can be retested within the same time window. Under these conditions, the probability that $q$ tests are returned positive and $Q^{-}=Q-q$ tests are returned negative is

Equation (3.1) describes perfect error-free and random testing. However, if there is some prior suspicion of being infected, the administration of testing may be biased. For example, certain jurisdictions focus testing primarily on hospitalized patients and people with significant symptoms [17], thus biasing the tests to those that are infected. We quantify such testing biases through a biased-testing function $K(Q^{+})\in {\mathbb{R}}_{\ge 0}$, leading to the following modification of equation (3.1):

The expected value of $f_{\text{b}}$, $\overline{\mu }(f,b)$, can be understood as a product of the true underlying infected fraction $f$ and a bias function $B(f,b)$ that depends on $f$ and the bias parameter $b$, i.e. $\overline{\mu }(f,b)=fB(f,b)$, where $B(f,b)$ is given by

Figure 1a shows the bias function (3.9) as a function of $b$ for different infection fractions $f$. For $B(f,b)>1$, the biased-testing fraction $fB$ is larger than the unbiased-testing fraction $f$. The opposite holds for $B(f,b)<1$. If only susceptible individuals are tested (i.e. $b\to -\mathrm{\infty }$), the bias function $B(f,b)$ and the expected observed positive testing fraction $\overline{\mu }(f,b)=fB(f,b)$ approach zero. For a complete bias towards infected individuals (i.e. $b\to \mathrm{\infty }$), the bias function approaches $f^{-1}$ and $\overline{\mu }(f,b)\to 1$.

• Download figure
• Open in new tab
• Download PowerPoint

The variance ${\overline{\sigma }}^2$ of the Gaussian approximation (3.7) is plotted as a function of $b$ in figure 1b and exhibits a maximum value of $1/(4Q)$ at $b^{\ast }=\ln [(1-f)/f]$.

The probabilities $P_{\text{true}}$ derived in equations (3.1) and (3.3) correspond to ‘testing with replacement’. The opposite limit is ‘testing without replacement’; once an individual is tested they are labelled as such and removed from the pool of test targets, at least within the specified testing period. This concept of sampling with and without replacement commonly arises in the measurement of diversity in ecological settings [24]. Without replacement, and still under conditions of perfect random testing, two slightly different forms for $P_{\text{true}}$ arise for the different type of tests (e.g. antibody versus PCR/viral load). For antibody tests that perfectly identify recovered (or deceased) individuals as being previously infected, equation (3.1) is replaced by

To incorporate testing bias into the probabilities $P_{\text{true}}$ for testing without replacement, we first consider equation (3.11) where $\left(\genfrac{}{}{0ex}{}{I}{q}\right)\left(\genfrac{}{}{0ex}{}{S+R}{Q-q}\right)$ can be interpreted as the number of ways of distributing $q$ positive tests among $I$ infected individuals, and $Q-q$ negative tests among $S+R$ uninfected individuals. As in the biased-testing formulation of equation (3.3), we interpret the bias as a factor $K(q)$ that assigns more weight to tests in the $I$ or $S+R$ pools

The above choice for $K(q)$ also allows us to explicitly evaluate the denominator in equation (3.12)

(a) Testing errors

The probability distributions $P_{\text{true}}$ that we derived in equation (3.3) and in equations (3.10)–(3.11) assume that testing is error-free, i.e. that the false-negative rate $\text{FNR}=1-\text{TPR}=0$ and false-positive rate $\text{FPR}=1-\text{TNR}=0$, or equivalently that the true positive rate $\text{TPR}=1$ and the true negative rate $\text{TNR}=1$. To incorporate erroneous testing, we now construct the probability distribution of error-generated deviation $P_{\text{err}}({\tilde{Q}}^{+}=q|Q^{+}=q^{\prime })$ over the number of ‘apparent’ positives $q$ from tests that carry nonzero FPRs and FNRs, given that $q^{\prime }$ positives would be recorded if the tests were perfect. If $q$ apparent positive tests are tallied, $p$ of them might have been true positives drawn from the perfect-test positives $q^{\prime }$ in $\left(\genfrac{}{}{0ex}{}{q^{\prime }}{p}\right)$ ways, while the remaining $k=q-p$ apparent positives might have been erroneously counted as positives drawn from the $Q^{-}\equiv Q-q^{\prime }$ true negatives. The remaining $q^{\prime }-p$ true positive tests might have been erroneously tallied as false negatives, while the remaining $Q^{-}-k$ negative tests might have been correctly tallied as true negatives. Assuming non-zero $\text{FPR}$ and $\text{FNR}$, we find that the probability distribution of finding $0\le q\le Q$ apparent positive tests is

$$\begin{array}{rl}& P_{\text{err}}({\tilde{Q}}^{+}=q|Q^{+}=q^{\prime },Q,\text{FPR},\text{FNR})\\ & =\sum _{p=0}^q\left(\genfrac{}{}{0ex}{}{q^{\prime }}{p}\right){(\text{TPR})}^p{(\text{FNR})}^{q^{\prime }-p}\left(\genfrac{}{}{0ex}{}{Q-q^{\prime }}{q-p}\right){(\text{FPR})}^{q-p}{(\text{TNR})}^{Q-q^{\prime }-(q-p)},\end{array}$$ 3.17

where we invoked the identities $\text{TPR}+\text{FNR}=1$ and $\text{FPR}+\text{TNR}=1$. The total error-prone distribution $P_{\text{TOT}}({\tilde{Q}}^{+}=q|Q,S,I,R,b)$ of recording $q$ positives after having administered $Q$ tests under bias and/or testing errors is given by convolving the probability $P_{\text{err}}$ of finding $q$ apparent tests given $q^{\prime }$ true positive tests with the probability $P_{\text{true}}$ of finding $q^{\prime }$ positive tests under perfect testing (equation (3.4) or equation (3.15)):

$$\begin{array}{rl}& P_{\text{TOT}}({\tilde{Q}}^{+}=q|Q,S,I,R,b,\text{FPR},\text{FNR})\\ & =\sum _{q^{\prime }=0}^Q{P}_{\text{err}}({\tilde{Q}}^{+}=q|Q^{+}=q^{\prime },Q,\text{FPR},\text{FNR})P_{\text{true}}(Q^{+}=q^{\prime }|Q,S,I,R,b).\end{array}$$ 3.18

This convolution can be further simplified by taking the Gaussian limit of the binomial distributions that appear in $P_{\text{err}}$ and in $P_{\text{true}}$, once more invoking the de Moivre–Laplace theorem [22]. To be concrete, we use $P_{\text{true}}$ as written in equation (3.3) under the testing with replacement scenario and use the approximation (3.5). The same Gaussian approximation can be used for both binomial terms in $P_{\text{err}}$ from equation (3.17), under the conditions illustrated above, which yield $\text{TPR}(1-\text{TPR})Q^{+}>9$, $\text{FPR}(1-\text{FPR})(Q-Q^{+})>9$. These will complement the condition $\overline{\mu }(1-\overline{\mu })Q>9$ that arises from using equation (3.5) in lieu of equation (3.17). Once the summand in equation (3.17) has been expressed as a product of two Gaussians, we approximate the summation over $p\in (0,q)$ by an integral over $p\in (-\mathrm{\infty },\mathrm{\infty })$. Provided $q$ is sufficiently large, we find

$$P_{\text{err}}({\tilde{Q}}^{+}=q|Q,Q^{+}=q^{\prime },\text{FPR},\text{FNR})\approx \frac{\exp [-\frac{{(q-q^{\prime }\text{TPR}-(Q\hspace{0.17em}-\hspace{0.17em}{q}^{\prime })\text{FPR})}^2}{2q^{\prime }\text{FNR}\hspace{0.17em}\text{TPR}+2(Q\hspace{0.17em}-\hspace{0.17em}{q}^{\prime })\text{FPR}\hspace{0.17em}\text{TNR}}]}{\sqrt{2\pi }\sqrt{q^{\prime }\text{FNR}\hspace{0.17em}\text{TPR}+(Q\hspace{0.17em}-\hspace{0.17em}{q}^{\prime })\text{FPR}\hspace{0.17em}\text{TNR}}}.$$ 3.19

We now convolve equation (3.5) with equation (3.19) as prescribed by equation (3.18) to find

$$P_{\text{TOT}}({\tilde{Q}}^{+}=q|Q,f,b,\text{FPR},\text{FNR})\approx \frac{1}{{\sigma }_{\text{T}}\sqrt{2\pi }}\exp [-\frac{{(q-{\mu }_{\text{T}})}^2}{2{\sigma }_{\text{T}}^2}],$$ 3.20

where

$$\begin{array}{rl}{\mu }_{\text{T}}(f,b,\text{FPR},\text{FNR})& \equiv Q[\overline{\mu }(f,b)(1-\text{FNR})+(1-\overline{\mu }(f,b))\text{FPR}],\\ {\sigma }_{\text{T}}^2(f,b,\text{FPR},\text{FNR})& \equiv Q(1-\overline{\mu }(f,b))\text{FPR}(1-\text{FPR})+Q\overline{\mu }(f,b)\text{FNR}(1-\text{FNR})\\ & +Q\overline{\mu }(f,b)(1-\overline{\mu }(f,b)){(1-\text{FNR}-\text{FPR})}^2,\end{array}\}$$ 3.21

with $\overline{\mu }(f,b)$ the mean value of the fraction of observed positives under biased, but perfect testing (see equation (3.8)). To derive equation (3.20), we replaced the summation over $q^{\prime }\in (0,Q)$ by an integral over $q^{\prime }\in (-\mathrm{\infty },\mathrm{\infty })$. This approximation is valid if the value of the two Gaussians that approximate $P_{\text{err}}$ and $P_{\text{true}}$ in equation (3.18) are negligible when evaluated at $q^{\prime }\le 0$ and $q^{\prime }\ge Q$ for all values of $q$. We have verified that this is true provided that $0<2/(Q+2)\ll \text{FPR, FNR}$, $\overline{\mu }\ll Q/(Q+2)<1$. In the limit of large numbers of tests $Q$, this condition will almost always be satisfied. Finally, to perform the integral that now appears in equation (3.18), we fixed the value of $q^{\prime }=\overline{\mu }$ in the $q^{\prime }$-dependent variance in equation (3.19). We do this to obtain an analytic expression for the integral of $P_{\text{err}}{P}_{\text{true}}$ over $q^{\prime }$, and select the value $q^{\prime }=\overline{\mu }$ for the variance of $P_{\text{err}}$ as the one that carries the most weight from $P_{\text{true}}$.

Equation (3.20) reveals that the mean number of apparent positive tests ${\mu }_{\text{T}}$ is given by the sum of the expected value of true positive tests (i.e. $Q\overline{\mu }(f,b)(1-\mathrm{FNR})=Q\overline{\mu }(f,b)\mathrm{TPR}$) and the expected value of false positive tests (i.e. $Q(1-\overline{\mu }(f,b))\mathrm{FPR}$). Based on the derived expressions for ${\mu }_{\text{T}}$ and ${\sigma }_{\text{T}}^2$, we define the random variable as the fraction of observed positive tests ${\tilde{f}}_{\text{b}}\equiv {\tilde{Q}}^{+}/Q$ under biased and error-prone testing and obtain

$$P_{\text{TOT}}({\tilde{f}}_{\text{b}}=x|Q,f,b,\text{FPR},\text{FNR})\approx \frac{{\text{e}}^{-{(x-{\overline{\mu }}_{\text{T}})}^2/(2{\overline{\sigma }}_{\text{T}}^2)}}{{\overline{\sigma }}_{\text{T}}\sqrt{2\pi }},$$ 3.22

where

$$\begin{array}{rl}& {\overline{\mu }}_{\text{T}}(f,b,\text{FPR},\text{FNR})\equiv [\overline{\mu }(f,b)(1-\text{FNR})+(1-\overline{\mu }(f,b))\text{FPR}]\\ \text{and}& Q{\overline{\sigma }}_{\text{T}}^2(f,b,\text{FPR},\text{FNR})\equiv (1-\overline{\mu }(f,b))\text{FPR}(1-\text{FPR})+\overline{\mu }(f,b)\text{FNR}(1-\text{FNR})\\ & +\overline{\mu }(f,b)(1-\overline{\mu }(f,b)){(1-\text{FNR}-\text{FPR})}^2.\end{array}\}$$ 3.23

We have numerically verified that the Gaussian approximation (3.23) is quite accurate provided that (i) the number of positive and apparent positive tests, $Q^{+}=q^{\prime }$ and ${\tilde{Q}}^{+}=q$, are sufficiently large, and (ii) the quantities $\overline{\mu }$, $\text{FNR}$ and $\text{FPR}$ are not too close to 0 or 1. This is shown in figure 2 and is in accordance with the conditions implied by using the de Moivre–Laplace theorem and by exchanging the sum with the integration when evaluating equation (3.18), as described above. Figure 2 shows the distribution of apparently infected individuals $P_{\text{TOT}}$ for different numbers of infected/recovered individuals (figure 2a), testing biases (figure 2b), and testing sensitivities (i.e. true positive rates, $\text{TPR}$s) and specificities (i.e. true negative rates, $\text{TNR}$s) (figure 2c–d). Solid light blue lines (colour available in online version) represent the Gaussian approximation (3.22) and dashed black lines and the remaining thicker coloured lines are calculated by directly evaluating equation (3.18) with replacement (equation (3.3)) and without replacement (equation (3.16)), respectively. The $\text{FNR}$s that we consider in figure 2c are chosen in accordance with reported sensitivities of serological and RT-PCR tests for SARS-CoV-2 [9–12]. We observe that an increase in the $\text{FNR}$ slightly shifts the distribution $P_{\text{TOT}}$ towards smaller values of apparently infected individuals, which is consistent with the $\mathrm{FNR}$ dependence of the mean ${\overline{\mu }}_{\text{T}}$ (equation (3.23)). For serological and RT-PCR tests, the $\mathrm{FPR}=1-\mathrm{TNR}$ is about 5%. A smaller specificity would lead to larger $\text{FPR}$s and a shift of $P_{\text{TOT}}$ towards larger values of ${\tilde{Q}}^{+}$ and ${\tilde{f}}_{\text{b}}$ (figure 2d). Our results show that $P_{\text{TOT}}$ is more affected by variations in the testing specificity than by variations in testing sensitivity.

• Download figure
• Open in new tab
• Download PowerPoint

4. Inference of prevalence and application to COVID-19 data

One often wishes to infer the evolution of $I_t+R_t$ and $S_t$, or $I_t$ and $S_t+R_t$ over a given time period from values of $b_{t,c}$, $Q_{t,c}$ and $q_{t,c}$. Since $b_{t,c}$ is difficult to independently ascertain, one may only be able to infer ${(f_{\text{b}})}_{t,c}=f(S_t,I_t,R_t,b_{t,c})$. For a single test result ${\tilde{q}}_{t,c}$ (or ${({\tilde{f}}_{\text{b}})}_{t,c}$), we can generate the maximum likelihood estimate (MLE) of the bias-modified prevalence ${({\hat{f}}_{\text{b}})}_{t,c}$ by setting the measured value ${({\tilde{f}}_{\text{b}})}_{t,c}={\overline{\mu }}_{\text{T}}[{({\hat{f}}_{\text{b}})}_{t,c}]$ to find

As an example, we collected US testing data [25] from March 2020 to March 2021. Figure 3a shows the daily number of observed positive tests ${\tilde{Q}}_t^{+}$ (red bars) and the corresponding total daily number of tests $Q_t$ (blue bars). The 7-day average of the observed positive testing rate ${({\tilde{f}}_{\text{b}})}_t={\tilde{Q}}_t^{+}/Q_t$ is indicated by the black solid line. The first drop in ${({\tilde{f}}_{\text{b}})}_t$ in March 2020 was associated with the initially very limited number of available SARS-CoV-2 testing infrastructure followed by the ramping up of testing capacity. After new cases surged by the end of March and in April 2020, different types of stay-at-home orders and distancing policies with different durations were implemented across the USA [26]. In June and July 2020, reopening plans were halted and reversed by various jurisdictions to limit the resurgence of COVID-19 [27].

• Download figure
• Open in new tab
• Download PowerPoint

In figure 3b, we show the corrected proportion of positive tests ${\hat{f}}_t$, found by numerically inverting equation (4.1) for different $\mathrm{FPR}$, $\mathrm{FNR}$ and bias combinations. We observe that a small $\mathrm{FPR}=0.04$ shifts values ${({\tilde{f}}_{\text{b}})}_t\approx 0.05$ towards zero such that the corrected positive testing rate ${\hat{f}}_t\approx 0$. Reducing the $\mathrm{FNR}$ from $0.2$ to $0.1$ has only little effect on the corrected proportion of positive tests ${\hat{f}}_t$ (solid black and dashed lines in figure 3b). Accounting for a positive testing bias of $b=1$ (i.e. preferential testing of infected and symptomatic individuals by a factor of $e$), however, markedly changes the inferred ${\hat{f}}_t$ (dashed-dotted black line in figure 3b). Since the 7-day average of the number of tests $Q$ is about ${10}^6$ in the USA (figure 3a), the variance terms ${({\hat{\overline{\sigma }}}_{\text{T}})}_t^2$ are very small compared to the values of ${\hat{f}}_t$.

5. Inference of bias b

One way to estimate the testing bias $b$ is to identify a smaller subset of control tests within a jurisdiction that is believed to be unbiased and compare it with the reported fraction of positive tests obtained via standard (potentially biased) testing procedures. Given this scenario, we can derive a rather complete methodology to estimate bias by formally comparing the statistics of two sets of tests applied to the same population. The first set of control tests with testing parameters ${\theta }_0=\{Q_0,{\text{FPR}}_0,{\text{FNR}}_0\}$ is known to be unbiased (has prior distribution $\delta (b)$), while the second set is taken with known parameters $\theta =\{Q,\text{FPR},\text{FNR}\}$, but unknown testing bias $b$. For example, the control set may consist of a smaller number $Q_0$ of tests that are administered completely randomly, while the second set may be the scaled-up set of tests with $Q>Q_0$. Since both sets of tests are applied roughly at the same time to the same overall population, the underlying positive fraction $f$ is assumed to be the same in both test sets. We can then use Bayes’ rule on the first unbiased test set to infer $f$

Of course, a simpler MLE can also be applied to data by first inferring the most likely value of $f$ from the control test set. We can use the number of positive tests in the control sample $Q_0^{+}$ to define the variable ${\tilde{f}}_0=Q_0^{+}/Q_0$. One can then maximize $P_{\text{TOT}}({\tilde{f}}_0|f,{\theta }_0,b=0)$ with respect to $f$ and use this value $\hat{f}$ in $P_{\text{TOT}}({\tilde{f}}_{\text{b}}|\hat{f},\theta ,b)$. Maximizing $P_{\text{TOT}}({\tilde{f}}_{\text{b}}|\hat{f},\theta ,b)$ with respect to $b$ then gives the MLE estimate $\hat{b}$. We can use random and unbiased sampling results obtained in the German jurisdiction of Gangelt, North Rhine–Westphalia [18]. A total of 600 adult persons with different last names were randomly selected from a population of 12 597 and asked to participate in the study together with their household members. The resulting study comprised of $Q_0=919$ subjects who underwent serological and PCR testing between 31 March and 6 April 2020. The specificity and sensitivity corrected, unbiased positive test fraction was determined to be $f=15.53\mathrm{\%}$ (95% CI 12.31–18.96%). Thus, we use this value as an estimate for the true underlying positivity rate $\hat{f}$. The larger sample taken across North Rhine–Westphalia between 30 March 30 and 5 April 2020 was measured ($Q\approx 25\hspace{0.17em}000$) to be ${\tilde{f}}_{\text{b}}\approx 0.1$ [28]. Assuming that this value is also error-corrected, an estimate of the bias $\hat{b}$ in this main testing set can be found by solving ${\tilde{f}}_{\text{b}}\approx 0.1={\overline{\mu }}_{\text{T}}(\hat{f}=0.1553,\hat{b}$, $\text{FPR}=\text{FNR}=0)=\overline{\mu }(\hat{f}=0.1553,\hat{b})$ for $\hat{b}$. We find that the difference between the unbiased positive testing rate of 15.53% and 10% corresponds to a bias of $\hat{b}=-0.50$. This negative bias likely arises because Gangelt was an infection hotspot within the entire North Rhine–Westphalia region, so the control sample was probably not unbiased. For comparison, a higher biased positive testing rate of 20% would lead to an estimated testing bias $\hat{b}=0.31$.

The number of total deaths on 6 April 2020 amounted to 7. Hence, the corresponding estimate of the $\mathrm{IFR}$, the number of disease-induced deaths $D_t$ divided by the total number of cases $N_t$ at time $t$, in this jurisdiction on 6 April 2020 was $7/(0.1553\times \text{12,597})=0.36$% (95% CI 0.29–0.45%) [18]. If only a biased estimate of the proportion of positive cases is known and not the true value $f$, we can use our framework to distinguish between the true ${\text{IFR}}_t=D_t/(N_t-S_t)=D_t/(f{N}_t)$ and the observed infection fatality ratio

• Download figure
• Open in new tab
• Download PowerPoint

6. Summary and conclusion

Radiological testing methods such as chest computed tomography are used sporadically to identify COVID-19-induced pneumonia in patients with negative tests [29]. However, the overwhelming majority of COVID-19 tests are based on serological (or antibody) tests and rapid antigen tests, ELISA and RT-PCR assay [1]. These tests are designed to subsequently output a binary signal, either infected or not. The population statistics of this output are affected by testing errors and bias. False-positive and false-negative rates of serological tests are generally smaller than those of rapid antigen tests and RT-PCR tests. However, serological tests are unable to identify early-stage infections since they are measuring antibody titres that usually develop a few days up to a few weeks after infection. In addition to the occurrence of false positives and false negatives (i.e. type I and type II errors), certain demographic groups (e.g. elderly people or those with comorbidities such as heart and lung diseases) may be overrepresented in testing statistics.

To quantify the impact of both type I/II errors and testing bias on reported COVID-19 case and death numbers, we developed a mathematical framework that describes erroneous and biased sampling (both with and without replacement) from a population of susceptible, infected and removed (i.e. recovered or deceased) individuals. We identify a positive testing bias $b>0$ with an overrepresentation of previously or currently infected individuals in the study population. Conversely, a negative testing bias $b<0$ corresponds to an overrepresentation of susceptible and/or non-infectious individuals in the study population. We derived MLEs of the testing-error and testing-bias-corrected fraction of positive tests. Our methods can be also applied to infer the full distribution of corrected positive testing rates over time and for different types of tests across different jurisdictions.

The mathematical quantity that underlies most of our analysis is the proportion of apparent positive tests. As pointed out in [30], the absolute number of positive tests may not capture the actual growth of an epidemic due to limitations in testing capacity. Still, many jurisdictions report absolute case numbers without specifying the total number of tests or additional information about test type, date of test and duplicate tests [31], rendering interpretation and application to epidemic surveillance challenging. For a reliable picture of COVID-19 case numbers, more complete testing data, including total number of tests, number of positive tests, test type and date of test, has to be reported and made publicly available at online data repositories. To correct for false positive, false negatives and testing bias in testing statistics (figure 3), it will be also important to further improve estimates of $\text{FPR}$, $\text{FNR}$, and $b$ through in field studies. In particular, estimating the testing bias $b$ requires random sampling studies similar to that carried out in [18]. Finally, while we have presented our analysis in the context of the COVID-19 pandemic, the general results presented in this paper apply to testing and estimation of severity of any infectious disease afflicting a population.

Data accessibility

Our source codes are publicly available at https://github.com/lubo93/disease-testing.

Authors' contributions

L.B., M.R.D. and T.C. contributed equally to the study design, data analyses and manuscript writing.

Competing interests

The authors declare no competing interests.

Funding

L.B. acknowledges financial support from the Swiss National Fund (P2EZP2_191888). The authors also acknowledge financial support from the Army Research Office (W911NF-18-1-0345), the NIH (R01HL146552) and the National Science Foundation (DMS-1814364, DMS-1814090).