Risk score learning for COVID-19 contact tracing apps

In this paper, we show that it is possible to learn interpretable risk score models using standard machine learning tools, such as multiple instance learning and stochastic gradient descent, provided we have suitable data. However, some new techniques are also required, such as replacing hard thresholding with soft binning, and using constrained optimization methods to ensure monotonicity of the learned function. We believe these methods could be applied to learn other kinds of risk score models.

Let $\tau_{n}$ be the duration (in minutes) of the $n$’th exposure, and let $d_{n}$ be the distance (in meters) between the two people during this exposure. (For simplicity, we assume the distance is constant during the entire interval; we will relax this assumption later.) Finally, let $\sigma_{n}$ be the time since of the onset of symptoms of the index case at the time of exposure, i.e., $\sigma_{n}=t^{\text{sym}}_{i_{n}}-t^{\text{exp}}_{j_{n}}$, where $i_{n}$ is the index case for this exposure, $j_{n}$ is the user, $t^{\text{sym}}_{i_{n}}$ is the time of symptom onset for $i_{n}$, and $t^{\text{exp}}_{j_{n}}$ is the time that $j_{n}$ gets exposed to $i_{n}$. (If the index case did not show symptoms, we use the date that they tested positive, and shift it back in time by 7 days, as a crude approximation.) Note that $\sigma_{n}$ can be negative. In (Ferretti et al., 2020a), they show that $\sigma_{n}$ can be used to estimate the degree of contagiousness of the index case.

Given the above quantities, we define the ”hazard score” for this exposure as follows:

$\displaystyle s_{n}$

$\displaystyle=f_{\text{hazard}}(\tau_{n},d_{n},\sigma_{n};\phi)=\tau_{n}\times f_{\text{dist}}(d_{n};\phi)\times f_{\text{inf}}(\sigma_{n};\phi)$

(1)

where $\tau_{n}$ is the duration, $d_{n}$ is the distance, $\sigma_{n}$ is the time since symptom onset, $f_{\text{dist}}(d_{n};\phi)$ is the simulated risk given distance, $f_{\text{inf}}(\sigma_{n};\phi)$ is the simulated risk given time, and $\phi$ are parameters of the simulator (as opposed to $\psi$, which are parameters of the risk score model that the phone app needs to learn).

The ”correct” functional form for the dependence on distance is unknown. In this paper, we follow (Briers et al., 2020), who model both short-range (droplet) effects, as well as longer-range (aerosol) effects, using the following simple truncated quadratic model:

$\displaystyle f_{\text{dist}}(d;\phi)=\min(1,D^{2}_{\min}/d^{2})$

(2)

They propose to set $D^{2}_{\min}=1$ based on an argument from the physics of COVID-19 droplets.

The functional form for the dependence on symptom onset is better understood. In particular, (Ferretti et al., 2020a) consider a variety of models, and find that the one with the best fit to the empirical data was a scaled skewed logistic distribution:

$\displaystyle f_{\text{inf}}(\sigma;\phi)\propto\text{SkewedLogistic}(\sigma|\mu=-4,\sigma=1.85,\alpha=5.85,\tau=5.42)$

(3)

In order to convert the hazard score into a probability of infection (so we can generate samples), we use a standard exponential dose response model (Smieszek, 2009; Haas, 2014):

$\displaystyle p_{n}$

$\displaystyle=\Pr(y_{n}=1|x_{n};\phi)=1-e^{-\lambda s_{n}}$

(4)

where $y_{n}=1$ iff exposure $n$ resulted in an infection, and $x_{n}=(\tau_{n},d_{n},\sigma_{n})$ are the features of this exposure. The parameter $\lambda$ is a fixed constant with value $3.1\times 10^{-6}$, chosen to match the empirical attack rate reported in (Wilson et al., 2020).

A user may have multiple exposures. Following standard practice, we assume each exposure could independently infect the user, so

	$\displaystyle P_{j}$	$\displaystyle=\Pr(Y_{j}=1|\mathcal{X}_{j};\phi)=1-\prod_{n\in E_{j}}(1-p_{n})$		(5)
		$\displaystyle=1-\prod_{n\in E_{j}}e^{-\lambda s_{n}}=1-e^{-\lambda\sum_{n}s_{n}}$		(6)
		$\displaystyle=p_{1}+(1-p_{1})\times p_{2}+(1-p_{1})(1-p_{2})\times p_{3}+\cdots$		(7)

where $E_{j}$ are all of $j$’s exposure events, and $\mathcal{X}_{j}$ is the corresponding set of features.

To account for possible ”background” exposures that were not recorded by the app, we can add a $p_{0}$ term, to reflect the prior probability of an exposure for this user (e.g., based on the prevalence). This model has the form

$\displaystyle P_{j}$

$\displaystyle=1-(1-p_{0})\prod_{n\in E_{j}}(1-p_{n})$

(8)

Suppose we have a set of $J$ users. The probability that user $j$ gets infected is given by Eq. 8. However, it remains to specify the set of exposure events $E_{j}$ for each user. To do this, we create a pool of $N$ randomly generated exposure events, uniformly spaced across the 3d grid of distances (80 points from 0.1m to 5m), durations (20 points from 5 minutes to 60 minutes), and symptom onset times (-10 to 10 days). This gives us a total of $N=33600$ exposure events.

We assume that each user $j$ is exposed to a random subset of one or more of these events (corresponding to having encounters with different index cases). In particular, let $u_{jn}=1$ iff user $j$ is exposed to event $n$. Hence the probability that $j$ is infected is given by $P_{j}=1-\exp[-\lambda S_{j}]$, where $S_{j}=\sum_{n=1}^{N}u_{jn}s_{n}$ and $s_{n}=f_{\text{hazard}}(x_{n};\phi)$, where $x_{n}=(\tau_{n},d_{n},\sigma_{n})$. We can write this in a more compact way as follows. Let $\mathbf{U}=[u_{jn}]$ be an $J\times N$ binary assignment matrix, and $\mathbf{s}$ be the vector of event hazard scores. For example, suppose there are $J=3$ users and $N=5$ events, and we have the following assignments: user $j=1$ gets exposed to event $n=1$; user $j=2$ gets exposed to events $n=1$, $n=2$ and $n=3$; and user $j=3$ gets exposed to events $n=3$ and $n=5$. Thus

$\displaystyle\mathbf{U}=\begin{pmatrix}1&0&0&0&0\\ 1&1&1&0&0\\ 0&0&1&0&1\end{pmatrix}$

(9)

The corresponding vector of true infection probabilities, one per user, is then given by $\mathbf{p}=1-\exp[-\lambda\mathbf{U}\mathbf{s}]$. From this, we can sample a vector of infection labels, $\mathbf{y}$, one per user. If a user gets infected, we assume they are removed from the population, and cannot infect anyone else (i.e., the pool of $N$ index cases is fixed). This is just for simplicity. We leave investigation of more realistic population-based simulations to future work.

In the real world, a user may get exposed to events that are not recorded by their phone. We therefore make a distinction between the events that a user was actually exposed to, encoded by $u_{jn}$, and the events that are visible to the app, denoted by $v_{jn}$. We require that $v_{jn}=1\implies u_{jn}=1$, but not vice versa. For example, Eq. 9 might give rise to the following censored assignment matrix:

$\displaystyle\mathbf{V}=\begin{pmatrix}1&0&0&0&0\\ 1&1&\mathbf{0}&0&0\\ 0&0&1&0&\mathbf{0}\end{pmatrix}$

(10)

where the censored events are shown in bold. This censored version of the data is what is made available to the user’s app when estimating the exposure risk, as we discuss in Sec. 4. We assume censoring happens uniformly at random. We leave investigation of more realistic censoring simulations to future work.

DCT apps do not observe distance directly. Instead, they estimate it based on bluetooth attenuation. In order to simulate the bluetooth signal, which will be fed as input to the DCT risk score model, we use the stochastic forwards model from (Lovett et al., 2020). This models the atttenuation as a function of distance using a log-normal distribution, with a mean given by

$$\mathbb{E}\left[{a|d;\phi}\right]=e^{\phi_{\alpha}+\phi_{\beta}\log(d)}=e^{\phi_{\alpha}}d^{\phi_{\beta}}$$

(11)

Using empirical data from MIT Lincoln Labs, Lovett et al. (2020) estimate the offset to be $\phi_{\alpha}=3.92$ and the slope to be $\phi_{\beta}=0.21$.

In this paper, we assume this mapping is deterministic, even though in reality it is quite noisy, due to multipath reflections and other environmental factors (see e.g., (Leith and Farrell, 2020)). Fortunately, various methods have been developed to try to ”denoise” the bluetooth signal. For example, the England/Wales GAEN app uses an approach based on unscented Kalman smoothing (Lovett et al., 2020). We leave the study of the effect of bluetooth noise on our learning methods to future work.

The GAEN risk score model makes a piecewise constant approximation to the risk-vs-attenuation function. In particular, it defines 3 attenuation thresholds, $\theta^{\text{ble}}_{1:3}$, which partitions the real line into $4$ intervals. Each interval or bucket is assigned a weight, $w^{\text{ble}}_{1:4}$, as shown in Fig. 1.

Overall, we can view this as approximating $f_{\text{dist}}(d_{n};\phi)$ by $f_{\text{ble}}(a_{n};\psi)$, where

$\displaystyle f_{\text{ble}}(a_{n};\psi)$

$\displaystyle=\begin{cases}w^{\text{ble}}_{1}&\mbox{if $a_{n}\leq\theta^{\text{ble}}_{1}$}\\ w^{\text{ble}}_{2}&\mbox{if $\theta^{\text{ble}}_{1}<a_{n}\leq\theta^{\text{ble}}_{2}$}\\ w^{\text{ble}}_{3}&\mbox{if $\theta^{\text{ble}}_{2}<a_{n}\leq\theta^{\text{ble}}_{3}$}\\ w^{\text{ble}}_{4}&\mbox{if $\theta^{\text{ble}}_{3}>a_{n}$}\end{cases}$

(14)

The health authority can compute the days since symptom onset for each exposure, $\sigma_{n}=t^{\text{sym}}_{i_{n}}-t^{\text{exp}}_{j_{n}}$, by asking the index case when they first showed symptoms. (If they did not show any symptoms, they can use a heuristic, such as the date of their test shifted back by several days.) However, sending the value of $\sigma_{n}$ directly was considered to be a security risk by Google/Apple. Therefore, to increase privacy, the symptom onset value is mapped into one of 3 infectiousness or contagiousness levels; this mapping is defined by a lookup table as shown in Fig. 2. The three levels are called ”none”, ”standard” and ”high”. We denote this mapping by

$$\underbrace{c_{n,l}}_{\text{contagiousness level}}=\mathbb{I}\left({\theta^{\text{con}}_{l-1}<\underbrace{\sigma_{n}}_{\text{time since symptom onset}}\leq\theta^{\text{con}}_{l}}\right)$$

Each infectiousness level is associated with a weight. (However, we require that the weight for level ”none” is $w^{\text{con}}_{1}=0$, so there are only 2 weights that need to be specified.) Overall, the mapping and the weights define a piecewise constant approximation to the risk vs symptom onset function defined in Eq. 3. We can write this approximation as $f_{\text{inf}}(\sigma_{n};\psi)\approx f_{\text{con}}(\text{LUT}(\sigma_{n});\psi)$, where

$\displaystyle f_{\text{con}}(c_{n};\psi)$

$\displaystyle=\begin{cases}0&\mbox{if $c_{n}=1$}\\ w^{\text{con}}_{2}&\mbox{if $c_{n}=2$}\\ w^{\text{con}}_{3}&\mbox{if $c_{n}=3$}\end{cases}$

(15)

In Fig. 3(left), we plot the probability of infection according to our simulator as a function of symptom onset $\sigma_{n}$ and distance $d_{n}$, for a fixed duration $\tau_{n}=15$. We see that this is the convolution of the (slightly asymmetric) bell-shaped curve from the symptom onset function in Eq. 3 along the x-axis with the truncated quadratic curve from the distance function in Eq. 2 along the z-axis. We can approximate this risk surface using the GAEN piecewise constant approximation, as shown in Fig. 3(right).

We combine the risk from multiple exposure windows by adding the risk scores, reflecting the assumption of independence between exposure events. For example, if there are multiple exposures within an exposure window, each with different duration and attenuation, we can combine them as follows. First we compute the total time spent in each bucket:

$\displaystyle\underbrace{\tau_{nb}}_{\text{duration in bucket $b$}}$

$\displaystyle=\sum_{k=1}^{K_{n}}\underbrace{\tau_{nk}}_{\text{duration}}\times\underbrace{\mathbb{I}\left({\theta^{\text{ble}}_{b-1}<a_{nk}\leq\theta^{\text{ble}}_{b}}\right)}_{\text{attenuation is in bucket $b$}}$

(16)

where $K_{n}$ is the number of ”micro exposures” within the $n$’th exposure window, $\tau_{nk}$ is the time spent in the $k$’th micro exposure, and $a_{nk}$ is the corresponding attenuation. (This gives a piecewise constant approximation to the distance-vs-time curve for any given interaction between two people.) We then compute the overall risk for this exposure using the following bilinear function:

$\displaystyle\underbrace{r_{n}}_{\text{risk score}}$

$\displaystyle=\underbrace{\left[\sum_{b=1}^{N_{B}}\tau_{nb}w^{\text{ble}}_{b}\right]}_{\text{weighted exposure minutes}}\times\underbrace{\left[\sum_{\ell=1}^{N_{C}}\mathbb{I}\left({c_{n,\ell}}\right)w^{\text{con}}_{\ell}\right]}_{\text{weighted contagiousness level}}$

(17)

where $N_{B}=4$ is the number of attenuation buckets, and $N_{C}=3$ is the number of contagiousness levels.

If we have multiple exposures for a user, we sum the risk scores to get $R_{j}=\sum_{n\in E_{j}}r_{n}$, which we can convert to a probability of infection using $Q_{j}=1-e^{-\mu R_{j}}$.

We optimize the parameters by maximizing the log-likelihood, or equivalently, minimizing the binary cross entropy:

$\displaystyle\mathcal{L}(\psi)=-\sum_{j=1}^{J}Y_{j}\log Q_{j}+(1-Y_{j})\log(1-Q_{j})$

(18)

where $Y_{j}\in\{0,1\}$ is the infection label for user $j$ coming from the simulator, and

$\displaystyle Q_{j}$

$\displaystyle=1-\exp\left[-\mu\sum_{n\in E_{j}}f_{\text{risk}}(\tilde{x}_{n};\psi)\right]$

(19)

is the estimated probability of infection, which is computed using the bilinear risk score model in Eq. 17.

A major challenge in optimizing the above objective arises due to the fact that a user may encounter multiple exposure events, but we only see the final outcome from the whole set, not for the individual events. This problem is known as multi-instance learning (see e.g., (Foulds and Frank, 2010)). For example, consider Eq. 9. If we focus on user $j=2$, we see that they were exposed to a ”bag” of 3 exposure events (from index cases 1, 2, and 3). If we observe that this user gets infected (i.e., $Y_{j}=1$), we do not know which of these events was the cause. Thus the algorithm does not know which set of features to pay attention to, so the larger the bag each user is exposed to, the more challenging the learning problem. (Note, however, that if $Y_{j}=0$, then we know that all the events in the bag must be labeled as negative, since we assume (for simplicity) that the test labels are perfectly reliable.)

In addition, not all exposure events get recorded by the GAEN app. For example, the index case $i$ may not be using the app, or the exposure may be due to environmental transmission (known as fomites). This can be viewed as a form of measurement ”censoring”, as illustrated in Eq. 10, which is a sparse subset of the true infection matrix in Eq. 9. This can result in a situation in which all the visible exposure events in the bag are low risk, but the user is infected anyway, because the true cause of the user’s infection is not part of $\tilde{\mathcal{X}}_{j}$. This kind of false positive is a form of label noise which further complicates learning.

In this section, we discuss some algorithmic issues that arise when trying to optimize the objective in Eq. 18.

In this section, we describe our experimental results.