Minimizing Uncertainty in Prevalence Estimates

Estimating prevalence – the proportion of a population that has been infected by a disease – is a fundamental problem in epidemiology. Nonetheless, many core mathematical issues associated with this task have only been recently discovered and understood [1, 2]. For example, it has long been assumed that classification of samples as positive or negative is necessary to compute the prevalence. In Ref. [1] we demonstrated that this is false: unbiased estimators of prevalence can be constructed from conditional probability arguments having nothing to do with classification.

The core idea of Ref. [1] was to recognize that the probability $Q(r)$ of a diagnostic measurement outcome $r$ in a space $\Omega\subset\mathbb{R}^{N}$ is given by the convex combination

$\displaystyle Q(r)=qP(r)+(1-q)N(r),$

(1)

where $q$ is the prevalence and $P(r)$, $N(r)$ are conditional probability density functions (PDFs) for positive and negative populations. Importantly, $P(r)$ and $N(r)$ can be constructed from training data and are thus known. Taking $D$ to be an arbitrary subset of $\Omega$, one can define $q$ via

$\displaystyle q=\frac{Q_{D}-N_{D}}{P_{D}-N_{D}},$

(2)

where $S_{D}$ indicates the measure of $D$ with respect to the arbitrary distribution $S$.¹¹1It is necessary to assume that: (i) $D$ has neither zero measure nor unit measure with respect to $Q$; and (ii) the measures with respect to $P$ and $N$ are not equal. These conditions are trivial to ensure in practice. Given $M$ samples drawn at random from a population and denoted $r_{j}$, Eq. (2) implies that $q$ can be estimated by

$\displaystyle q\approx\tilde{q}=\frac{\tilde{Q}_{D}-N_{D}}{P_{D}-N_{D}},\hskip 14.22636pt\tilde{Q}_{D}=\frac{1}{M}\sum_{j=1}^{M}\mathbb{I}(r_{j}\in D),$

(3)

where $\mathbb{I}$ is the indicator function. This estimate is unbiased and converges in mean square [1].

To solve the related and important problem of optimal prevalence estimation, we minimize the variance of $\tilde{q}$, which is proportional to

$\displaystyle\sigma^{2}=\frac{Q_{D}(1-Q_{D})}{(P_{D}-N_{D})^{2}}.$

(4)

Equation (4) arises from the variance of the binomial random variable $\tilde{Q}_{D}$ but is unusual in that the denominator depends on a set that can be deformed arbitrarily. Thus, minimization is performed with respect to both the parameter $Q_{D}$ as well as the $D$. We reduce Eq. (4) to a one-dimensional (1D) problem by maximizing the denominator for fixed $Q_{D}$. We construct this maximum in terms of a bathtub principle.

Distinct sets $D$ may yield the same $Q_{D}$ but different realizations of $P_{D}-N_{D}$. This motivates treating $Q_{D}$ as an independent variable, denoted by $\hat{Q}$ to avoid confusion. By rewriting Eq. (2) in terms of $\hat{Q}$, we find a constraint

$\displaystyle\hat{Q}=q(P_{D}-N_{D})+N_{D}$

(5)

that defines the admissible $D$ corresponding to each $\hat{Q}$. Clearly the collection of sets $\{D^{\star}\}_{\hat{Q}}$ whose elements maximize $|P_{D}-N_{D}|$ ($|\cdot|$ denotes absolute value) for a fixed $\hat{Q}$ is an equivalence class. Since, each element in this class yields the same $F(\hat{Q})=(P_{D^{\star}}-N_{D^{\star}})^{2}$ as a function of $\hat{Q}$, minimizing Eq. (4) is equivalent to solving the 1D problem

$\displaystyle\hat{Q}^{\star}=\mathop{\rm argmin}_{\hat{Q}}\frac{\hat{Q}(1-\hat{Q})}{F(\hat{Q})}.$

(6)

We next construct $F(\hat{Q})$ by finding the equivalence class $\{D^{\star}\}_{\hat{Q}}$.

Equation (5) indicates that “loading” points $r$ into $D$ to increase $|P_{D}-N_{D}|$ simultaneously increases $N_{D}$. We therefore wish to construct $D$ so as to increase $|P_{D}-N_{D}|$ as much as possible and $N_{D}$ as little as possible. This motivates us to define the sets

	$\displaystyle D_{\pm}$	$\displaystyle=\{r:\pm q[P(r)-N(r)]>\pm\Delta_{\pm}N(r)\}\cup S_{\pm}$		(7)
	$\displaystyle\hat{Q}$	$\displaystyle=\int_{D_{\pm}}Q(r)dr.$		(8)

where

$\displaystyle S_{\pm}\subset\{r:q[P(r)-N(r)]=\Delta_{\pm}N(r)\}$

(9)

and Eq. (8) defines $\Delta_{\pm}$ in Eq. (7). The $S_{\pm}$ are any sets satisfying both Eqs. (9) and the constraints given by Eq. (8). An optimality proof for these sets is closely related to the bathtub principle [3]. For completeness and to facilitate comparison to classification methods [1], we present a proof adapted to the problem at hand.

Proof: First we show that Eq. (8) defines $\Delta_{\pm}$. By construction, $\Delta_{\pm}$ satisfies the inequality $-q\leq\Delta_{\pm}<\infty$. Restrict attention to $\Delta_{+}$. Let $D(\Delta)=\{r:q[P(r)-N(r)]>\Delta N(r)\}$, where $-\infty<\Delta<\infty$. Define $\displaystyle D(\infty)=\lim_{\Delta\to\infty}D(\Delta)=\{r:P(r)>0,N(r)=0\}$. Clearly,

$\displaystyle\mathcal{Q}(\Delta)=\left[\int_{D(\Delta)}Q(r)dr\right]-\hat{Q}$

(10)

is a monotone decreasing function of $\Delta$. Assume that $\mathcal{Q}(\Delta)$ crosses zero as $\Delta\to\infty$. Then $\mathcal{Q}(\Delta)$ is either continuous or discontinuous at this crossing, which defines $\Delta_{+}$ as the corresponding right-limit; that is $\Delta_{+}={\rm inf}_{\Delta}\{\Delta:\mathcal{Q}(\Delta)<0\}$. The $S_{+}$ can then be chosen as any subset of $\{q[P(r)-N(r)]=\Delta_{+}N(r)\}$ for which $\mathcal{Q}(\Delta_{+})+\int_{S+}Q(r)dr=0$. If instead $\mathcal{Q}(\Delta)$ does not vanish for any finite $\Delta$, then $\Delta$ is not finite and $S_{+}$ is a subset of $D(\infty)$. A similar argument yields existence of $\Delta_{-}$ and $S_{-}$.

Let $D$ be any set that differs from $D_{+}$ by positive measure not on $S_{+}$. Requiring that Eq. (5) also hold for $D$ yields

$\displaystyle P_{D_{+}}-N_{D_{+}}-P_{D}+N_{D}=\frac{1}{q}[N_{D/D_{+}}-N_{D_{+}/D}].$

(11)

where $/$ is the set difference operator. Combining the definition of $\Delta_{+}$ with Eq. (5), one finds that $N_{D/D_{+}}>N_{D_{+}/D}$, implying $P_{D_{+}}-N_{D_{+}}>P_{D}-N_{D}$. A similar argument can be used to show that $D_{-}$ maximizes the difference $N_{D}-P_{D}$. These differences are unique, as can be verified by the definition of $S_{\pm}$. Finally, note that $F(\hat{Q})$ is maximized by $\max_{\hat{Q}}[P_{D_{+}}-N_{D_{+}},N_{D_{-}}-P_{D_{-}}]$. \qed

The denominator of Eq. (4) can be parameterized by $\hat{Q}$. In particular, we have shown that for a fixed $\hat{Q}$ satisfying $0<\hat{Q}<1$, one of the variances

$\displaystyle\sigma^{2}_{+}(\hat{Q})=\frac{\hat{Q}(1-\hat{Q})}{(P_{D_{+}}-N_{D_{+}})^{2}},\hskip 14.22636pt\sigma^{2}_{-}(\hat{Q})=\frac{\hat{Q}(1-\hat{Q})}{(P_{D_{-}}-N_{D_{-}})^{2}}$

(12)

minimizes $\sigma^{2}$. However, we have yet to uniquely define the objective function in Eq. (6), since it not clear when to use $\sigma^{2}_{+}$ or $\sigma^{2}_{-}$. We now prove that both are equivalent.

Proof: The numerators of $\sigma^{2}_{+}$ and $\sigma^{2}_{-}$ are invariant to the transformation $\hat{Q}\to 1-\hat{Q}$. Also, for any $D_{+}$,

$\displaystyle P_{D_{+}}-N_{D_{+}}=N_{\Omega/D_{+}}-P_{\Omega/D_{+}}.$

(13)

By Eqs. (7)–(9), the set $\Omega/D_{+}$ is equal to a set $D_{-}$ for which $\Delta_{-}=\Delta_{+}$. Moreover, since $\hat{Q}$ is the $Q$-measure of domain $D_{+}$, the complement $\Omega/D_{+}$ has $Q$-measure $1-\hat{Q}$. That is, $\Omega/D_{+}\in\{D_{-}\}_{(1-\hat{Q})}$. In light of Eq. (13), one finds $\sigma^{2}_{+}(\hat{Q})=\sigma^{2}_{-}(1-\hat{Q})$. \qed

Lemma 2 proves that we may treat $\sigma^{2}_{+}$ (or $\sigma^{2}_{-}$) as the objective function in Eq. (6). We now show that this objective function has desirable properties.

Proof: By Lemma 2, it is sufficient to only consider $\sigma^{2}_{+}(\hat{Q})$. Because any point $r$ has zero measure, by the definition of $D_{\pm}$, the difference $P_{D_{+}}-N_{D_{+}}$ is continuous. Thus, $\sigma_{+}^{2}(\hat{Q})$ is continuous on $(0,1)$ provided that $P_{+}-N_{+}>0$ for every $\hat{Q}\in(0,1)$. To demonstrate this, assume that there exists a $\hat{Q}_{\infty}$ for which $P_{D_{+}}-N_{D_{+}}=0$. By assumption (ii) and the definition of $D_{+}$, the remaining points in $\Omega/D_{+}$ must have the property that $P(r)<N(r)$. Integrating over this set shows that $P_{\Omega}<N_{\Omega}$, which violates the assumption that both are probability densities. Thus $P_{D_{+}}-N_{D_{+}}$ cannot be zero in the interior of the domain, and $\sigma_{+}^{2}(\hat{Q})$ is continuous on $(0,1)$. By definition, $P_{D_{+}}\leq\hat{Q}$ and $N_{D_{+}}\leq\hat{Q}$, which implies that $\sigma^{2}_{+}(\hat{Q})\to\infty$ in the limit $\hat{Q}\to 0$. The corresponding result in the limit $\hat{Q}\to 1$ is proved in the same way by considering $\hat{Q}\to 0$ for $\sigma^{2}_{-}$ and then using Lemma 2.

Existence of the minimum follows from divergence of $\sigma_{+}^{2}(\hat{Q})$ at the boundaries and continuity in the interior $(0,1)$. Specifically, for any $0<\epsilon<1/2$, $\sigma_{+}^{2}$ is continuous on $D_{\epsilon}=[\epsilon,1-\epsilon]$. By the extreme value theorem, there exists a value $\hat{Q}_{\epsilon}\in D_{\epsilon}$ for which $\sigma_{+}^{2}(\hat{Q})$ attains a minimum. Clearly there exists an $\epsilon_{0}$ such that for all $\epsilon<\epsilon_{0}$, $\sigma_{+}^{2}(\hat{Q}_{\epsilon})$ is constant, which defines the minimum. \qed

Remark: Assumption (ii) implies that there is a set of positive measure with respect to $N$ for which $N(r)>P(r)$.

Remark: The assumption that there exists a set of positive measure for which $P(r)>N(r)$ is an important feature of diagnostic tests; it implies the ability to distinguish populations. Failure to satisfy this condition is the hallmark of a useless diagnostic.