Estimating County-Level COVID-19 Exponential Growth Rates Using Generalized Random Forests

10/31/2020 ∙ by Zhaowei She, et al. ∙ Georgia Institute of Technology Harvard University 0

Rapid and accurate detection of community outbreaks is critical to address the threat of resurgent waves of COVID-19. A practical challenge in outbreak detection is balancing accuracy vs. speed. In particular, while accuracy improves with estimations based on longer fitting windows, speed degrades. This paper presents a machine learning framework to balance this tradeoff using generalized random forests (GRF), and applies it to detect county level COVID-19 outbreaks. This algorithm chooses an adaptive fitting window size for each county based on relevant features affecting the disease spread, such as changes in social distancing policies. Experiment results show that our method outperforms any non-adaptive window size choices in 7-day ahead COVID-19 outbreak case number predictions.

READ FULL TEXT VIEW PDF
POST COMMENT

Comments

There are no comments yet.

Authors

page 1

page 2

page 3

page 4

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

Early and accurate detection of community outbreaks is critical to address the threat of resurgent waves of COVID-19. Specifically, an epidemic outbreak is confirmed when its incident cases are estimated to grow exponentially. Furthermore, the potential impact of an outbreak is also measured by its exponential growth rate, as higher rates indicate more rapid disease spread. Last but not least, for a fixed epidemiology model (e.g. SIR, SEIR), there is an one-to-one correspondence between the exponential growth rate of an epidemic outbreak and its basic reproduction number , a common measure of intensity of epidemic outbreaks (lipsitch2003transmission). Therefore, the exponential growth rate of an epidemic outbreak’s incident cases is the most important ”model-free” parameter to estimate for detecting the outbreak (chowell2003sars).

It remains an epidemiological challenge to obtain accurate exponential growth rate estimates for disease outbreaks (ma2014estimating)

. Specifically, it is difficult to choose the fitting window size for the exponential growth rate estimation to balance the speed and accuracy of outbreak detection. On one hand, a longer fitting window is preferable as larger sample size would reduce variance of the exponential growth rate estimates of outbreaks. On the other hand, shorter fitting windows are better at detecting early-stage outbreaks, especially if these outbreaks were driven by recent policy changes such as school reopening. In the current practice, this fitting window size is treated as a hyperparameter that is either directly specified by the user (c.f.

Mel) or determined by some cross validation methods (c.f. chowell2007comparative).

This paper develops a machine learning framework that balances the speed-accuracy tradeoff of outbreak detection via dedicated feature engineering and GRF (c.f. athey2019generalized), and apply it to detect county-level COVID-19 outbreaks. Specifically, the algorithm chooses an adaptive fitting window size for each county based on a rich set of features that affect the disease spread, such as face mask mandates, social distancing policies, the CDC’s Social Vulnerability Index, changes in tests performed and rate of positive tests. Furthermore, for counties with insufficient data to capture the recent policy changes, the algorithm pools together all relevant incident case growth trends across U.S. counties and throughout the COVID-19 pandemic history to adjust for these policy changes.

2 Background

2.1 Exponential Growth Model and Exponential Growth Rate

During an epidemic outbreak, incident case number at county on day , i,e. , is governed by an exponential growth model,

(1)

(ma2014estimating).111See Appendix B for the epidemiological definition of incident case number and how we compute . Here is the incident cases exponential growth rate of outbreaks for all counties , while captures the initial incident case numbers at the beginning of the exponential case growth at county . To obtain the most recent exponential growth rate of county-level COVID-19 incident case number, we estimate the instantaneous counterpart of (1),

(2)

where the dependent variable is the log-linearized incident case number ; independent variable is day . The parameter of interest is the COVID-19 incident case exponential growth rate of county at day , . The intercept term, , is a parameter capturing the log-linearized initial incident case number of county at the beginning of the outbreak. is the error term.

Notably, the exponential growth rate in (2) varies in both day and county . In other words, we are interested in estimating the instantaneous county-level exponential growth rate of COVID-19 incident cases. Specifically, since COVID-19 related regulations are changing every day in the U.S. (c.f. CUSP), the county-level exponential growth rates, affected by these policies, also change from day to day. Therefore, in order to detect recent outbreaks, we need to estimate the most recent incident cases exponential growth rate. This instantaneous exponential growth rate , similar to the instantaneous reproduction number commonly used in the epidemiology literature (c.f. fraser2007estimating; cori2013new), captures the projected exponential growth rate of incident cases in county should the future COVID-19 regulations remain the same as those in day .

2.2 Relevant Features Affecting COVID-19 Disease Spread

The exponential growth rate of COVID-19 incident cases can be affected by many factors ranging from day-to-day changes in COVID-19 regulations to difference in population density and healthcare resources between counties. Hence, in order to estimate the instantaneous county-level exponential growth rate defined in (2), we need to control for these day-level and county-level heterogeneity. Provided that we have relevant features that captures these aforementioned factors affecting the instantaneous county-level exponential growth rate ,222Refer to Appendix C for the list of relevant features used in this study. we can identify the conditional average partial treatment effect as defined by wooldridge2010econometric. When data of these relevant features affecting COVID-19 disease spread is available, we can rewrite (2) as

(3)

following the common “redundant” assumption (c.f. wooldridge2010econometric), i.e.

3 Model Estimation

This section discusses how we estimate the exponential growth rate of COVID-19 incident cases defined in (3). Specifically, we first formulate the estimation problem in §3.1 and then present our estimation algorithm in §3.2.

3.1 Problem Formulation

First, we need the following “unconfoudness” assumption (c.f. rosenbaum1983central):

That is, day

is independent of all unobservable heterogeneity conditional on the feature vector

. Under Assumption 3.1

, we can derive the following moment equations for (

3): ,

(4)

through the instrumental variable , where is a Boolean function equal 1 if .

However, we note that the moment equations (4) cannot identify the exponential growth rate . Specifically, the moment equation system (4) is underidentified as the number of unknown parameters is exactly two times the number of moment equations. To address this problem, we assume that the exponential growth rate of day equal the local average causal response (c.f. abadie2003semiparametric; angrist1995two) from day to day , i.e.

Parameters are exactly identified by the moment equations (4) under Assumption 3.1.

3.2 Estimation Algorithms

Given the above assumptions, the exponential growth function of COVID-19 incident cases (3) can be estimated using the GRF algorithm by athey2019generalized. Specifically, Assumption 3.1 implies that the exponential growth rate of a county on day is identified by a targeted data “block” . Hence, when estimating , we can partition the panel dataset into data blocks, and pool those data blocks “similar” to the targeted block to construct an adaptive window size for this estimation. Particularly, the similarity measure is provided by the GRF algorithm. In Appendix A, Algorithm LABEL:alg:block provides the pseudocode to construct these data blocks. Algorithm LABEL:alg:grf explains how we feed these data blocks into the GRF algorithm to obtain the estimates of interest.

4 Performance Evaluation

To benchmark our method against non-adaptive window size choices, we compare the Mean Absolute Percentage Error (MAPE) of these methods’ 7-day ahead predictions. Specifically, we use the NYTimes COVID-19 Dataset (c.f. NYTimes2020) as the source of daily reported cases per county. We note that there are many data quality issues with these reported case counts, especially during the first few months of this pandemic (c.f. NYTimesIssue). As an example, bergman2020oscillations found that parts of the oscillation in U.S. COVID-19 case counts can be explained by increased testing during some week days and test backlogs during the weekend. To alleviate these issues, this section preprocesses the raw case data using 4-day moving average and only compares the median MAPEs. Additional results for performance evaluation are available at Appendix D.

As shown in Figure LABEL:fig:MAPE4 and Table LABEL:tab:4, our method outperforms methods with fixed 2-, 4-, 8- or 16-day fitting window sizes. Specifically, Figure LABEL:fig:MAPE4 shows that our method provides a uniformly better performance roughly 100 days after the first recorded COVID-19 case in the NYTimes COVID-19 Dataset, when there were enough historical data for the GRF algorithm to conduct meaningful partition. Furthermore, Table LABEL:tab:4 demonstrates that even if the early-day MAPEs are included in the comparison, our method still has the best median MAPE among the 4 methods. Last but not least, we note that when only comparing the performance of non-adaptive window size choices, there are no obvious best choice. While choosing shorter fitting window sizes could in general lead to lower median MAPEs (c.f. Table LABEL:tab:4), it was still frequently outperformed by longer fitting window size choices (c.f. Figure LABEL:fig:MAPE4). This inherent difficulty in non-adaptive window size choices suggests that treating fitting window size as a hyperparameter and choosing it through cross validation methods (c.f. chowell2007comparative) cannot balance the the speed-accuracy tradeoff in outbreak detection.

fig:MAPE4

Figure 1: MAPE Plot (4-day Moving Average)

tab:4 Method MAPE OLS.wsize=16 0.0585 OLS.wsize=8 0.0529 OLS.wsize=4 0.0483 OLS.wsize=2 0.0455 GRF 0.0423

Table 1: Median RMSE and MAPE (4-day Moving Average)

5 Conclusion and Future Work

In this work, we developed a novel framework in feature engineering that allows GRF to match observations across time and space, and adequately balance the the speed-accuracy tradeoff in COVID-19 outbreak detection. This estimation framework can be readily extended to other panel data estimation problems in epidemiology where fitting window size choices play an important role, e.g. estimating the county-level instantaneous reproduction number (c.f. fraser2007estimating; cori2013new).

References

Appendix A Pseudocode for Estimation Algorithms

alg:blockInput:
Output: ,
for  to  do
        for  to  do
              
              Normalize the incident case numbers for each block
              Generate the initial rough OLS estimates for each block
              Append the feature data for each block
        end for
       
end for
Algorithm 1 Block Transformation
alg:grfInput: ,
Output:
Compute the congruence classes
for day
Assign the feature variable
for day
where
(5)
Assign the outcome variable
for day
where
(6)
Assign the treatment variable
for day
where
(7)
Feed (X,Y,W) into the GRF algorithm
Algorithm 2 GRF Training for Day

Appendix B Epidemiological Definitions

For this work, where we wish to provide useful estimates of case trajectories for epidemiologists and decision makers in the government, we ideally wish to measure the growth rate of active case numbers

(8)

where is the cumulative number of cases so far, is the cumulative number of deaths, and is the cumulative number of recovered cases. This intuitively captures the number of still ”infectious” cases, as we can assume for COVID-19 that those who have recovered or died are no longer capable of infecting others. Unfortunately, as the number of recovered cases are no longer reported, we are unable to directly calculate the number of active cases. While other works (c.f. Mel) approximate the number of recovered cases with

(9)

Here the underlying assumption is that those who were infected but not dead after 22 days have recovered. However, upon testing this assumption at the county-level, we find that this assumption does not hold as there will be some days where this newly approximated active case number becomes negative.

We hence rely upon a commonly used proxy: the incident case count, , defined as

(10)

which is the number of new cases within a 22 day time period and is a useful proxy for infectious cases.

Appendix C Feature Data

All feature data we used, i.e. , are publicly available online. In this section, we briefly describe what these datasets are and how we incorporated them in our work.

c.1 2019 US Census Gazetteer Files

The 2019 United States Census Gazetteer Files (c.f. USCensus2019) were used to obtain the geographic locations (latitude-longitude centroids) of each county officially registered by the United States Census Bureau. This provides a spatial feature space for the GRF to further split upon.

c.2 Centers for Disease Control and Prevention Social Vulnerability Index 2018 Database

The Social Vulnerability Index (SVI) database (c.f. CDCSVI2018) is a compilation of socio-economic factors such as unemployment rate, poverty rate, education attainment level etc. at the county level. These features are also included in our methodology.

c.3 COVID-19 US state policy database (CUSP)

The CUSP database (c.f. CUSP) tracks when each state implemented and ended policies such as mask mandates, lockdowns, economic policies in response to the COVID-19 pandemic. As such, it is a vector of features capturing daily policy changes in each state. As these policies is state-wide, we naturally extend them to the respective county level.

c.4 The COVID Tracking Project

From the COVID Tracking Project (c.f. COVIDTracking), we obtained the he daily numbers of PCR, Antibody and Antigen tests performed and their positivity rates at each state. These are used as features in our framework as well.

Appendix D Additional Results for Performance Evaluation

d.1 MAPE Plots

fig:MAPE3

Figure 2: MAPE Plot (3-day Moving Average)

fig:MAPE4L

Figure 3: MAPE Plot (4-day Moving Average)

fig:MAPE5

Figure 4: MAPE Plot (5-day Moving Average)

fig:MAPE6

Figure 5: MAPE Plot (6-day Moving Average)

fig:MAPE7

Figure 6: MAPE Plot (7-day Moving Average)

d.2 Median RMSEs and MAPEs

tab:3 Method RMSE MAPE OLS.wsize=16 0.3579 0.0584 OLS.wsize=8 0.3403 0.0532 OLS.wsize=4 0.3391 0.0510 OLS.wsize=2 0.3447 0.0511 GRF 0.2913 0.0468

Table 2: Median RMSE and MAPE (3-day Moving Average)

tab:4all Method RMSE MAPE OLS.wsize=16 0.3574 0.0585 OLS.wsize=8 0.3374 0.0529 OLS.wsize=4 0.3120 0.0483 OLS.wsize=2 0.3011 0.0455 GRF 0.2600 0.0423

Table 3: Median RMSE and MAPE (4-day Moving Average)

tab:5 Method RMSE MAPE OLS.wsize=16 0.3499 0.0576 OLS.wsize=8 0.3256 0.0509 OLS.wsize=4 0.2867 0.0450 OLS.wsize=2 0.2689 0.0410 GRF 0.2429 0.0395

Table 4: Median RMSE and MAPE (5-day Moving Average)

tab:6 Method RMSE MAPE OLS.wsize=16 0.3426 0.0570 OLS.wsize=8 0.3120 0.0490 OLS.wsize=4 0.2672 0.0420 OLS.wsize=2 0.2460 0.0378 GRF 0.2241 0.0361

Table 5: Median RMSE and MAPE (6-day Moving Average)

tab:7 Method RMSE MAPE OLS.wsize=16 0.3342 0.0559 OLS.wsize=8 0.2956 0.0467 OLS.wsize=4 0.2480 0.0392 OLS.wsize=2 0.2231 0.0346 GRF 0.2083 0.0336

Table 6: Median RMSE and MAPE (7-day Moving Average)

d.3 RMSE Plots

fig:RMSE3

Figure 7: RMSE Plot (3-day Moving Average)

fig:RMSE4

Figure 8: RMSE Plot (4-day Moving Average)

fig:RMSE5

Figure 9: RMSE Plot (5-day Moving Average)

fig:RMSE6

Figure 10: RMSE Plot (6-day Moving Average)

fig:RMSE7

Figure 11: RMSE Plot (7-day Moving Average)