Population density and basic reproductive number of COVID-19 across United States counties.

The basic reproductive number (R0) is a function of contact rates among individuals, transmission probability, and duration of infectiousness. We sought to determine the association between population density and R0 of SARS-CoV-2 across U.S. counties. We conducted a cross-sectional analysis using linear mixed models with random intercept and fixed slopes to assess the association of population density and R0, and controlled for state-level effects using random intercepts. We also assessed whether the association was differential across county-level main mode of transportation percentage as a proxy for transportation accessibility, and adjusted for median household income. The median R0 among the United States counties was 1.66 (IQR: 1.35-2.11). A population density threshold of 22 people/km2 was needed to sustain an outbreak. Counties with greater population density have greater rates of transmission of SARS-CoV-2, likely due to increased contact rates in areas with greater density. An increase in one unit of log population density increased R0 by 0.16 (95% CI: 0.13 to 0.19). This association remained when adjusted for main mode of transportation and household income. The effect of population density on R0 was not modified by transportation mode. Our findings suggest that dense areas increase contact rates necessary for disease transmission. SARS-CoV-2 R0 estimates need to consider this geographic variability for proper planning and resource allocation, particularly as epidemics newly emerge and old outbreaks resurge.

Introduction

The COVID-19 pandemic has infected millions of people globally, and there are over 400 thousand reported deaths and 7 million confirmed cases of COVID-19 worldwide [1]. Transmission of airborne and directly transmitted pathogens, such as SARS-CoV-2 (the causative agent of COVID-19), have been previously shown to be density-dependent [2-4]. Population density facilitates transmission of disease via close person-to-person contact [5-8], and may support sustained disease transmission due to increased contact rates [9-11]. Large urban areas have more opportunities for disease transmission, and hotspots of SARS-CoV-2 have been mostly concentrated in cities [12].

The basic reproductive number (R0) describes the contagiousness and transmissibility of pathogens, and is a function of contact rates among individuals, transmission probability (probability of transmission per contact), and duration of infectiousness [13]. This is in contrast with the time-varying reproductive number (Rt), defined as the number of people in a population who were infected by an infectious individual at a given point in time, which reflects the changing levels of immunity in the population and the impact of control measures limiting transmission [14]. Control measures include implementation of non-pharmaceutical interventions (NPIs), such as face coverings and social distancing. Thus, R0 estimates of COVID-19 are not exclusively determined by the pathogen, and variability in R0 depends on local sociobehavioral and environmental settings, including population density [12].

During the initial phase of the outbreak, or the exponential growth period, we hypothesize that spatial heterogeneity in R0 occurs primarily due to geographic variability in contact rates, since transmission probability and duration of infectiousness remain constant across settings. During this time frame, transmission probabilities across localities are equivalent, because exponential growth occurs prior to the implementation of NPIs which affect likelihood of transmission during contact. Moreover, contact networks are also affected by transportation systems that facilitate disease spread due to increased interconnectivity and mobility between different geographic areas [15, 16]; thus, we also hypothesize that the association of population density and R0 may be differential depending on transportation accessibility, and areas that lack access to efficient modes of transportation would not have the same SARS-CoV-2 growth rate, even in high density areas.

In the current COVID-19 pandemic, evidence of the association of population density and disease transmission have been conflicting [17-20]. The estimation of differential R0 using these area-level factors can assist in more accurate predictions of the rate of spread of SARS-CoV-2 in geographic settings with potential resurgence, where cases have been steadily increasing. In this study, we examine the association of population density with R0 of COVID-19 across United States counties.

Methods

Data

We obtained publicly available daily COVID-19 cases among United States counties from the New York Times [21]. For each county, we assumed that the exponential growth period began one week prior to the second daily increase in cases. We assumed that the period of exponential growth lasted approximately 18 days, as previous research have shown the COVID-19 exponential period to be around 20–24 days in New York City [22]. These assumptions calibrated the period of exponential growth accordingly and created reasonable curves that approximated exponential growth across the counties for case data (Fig A in S1 Appendix). The algorithm ensured that the virus had taken hold in the area and allowed a sufficient number of days to estimate the exponential growth rate, as R0 cannot be estimated accurately with sparse data, since it would be uncertain if the county was experiencing a sustained outbreak with community transmission. We restricted our analyses to counties that met a certain threshold of cumulative case counts at the end of the exponential growth period, but including counties with less than 25 cases included daily incidence counts that were insufficient to calculate R0 and yielded computational errors. Our final analytic sample included counties with 25 cases or greater at the end of the exponential growth period.

Data on the primary mode of transportation to work and median household income were obtained from the most recent 5-year American Community Survey (ACS) 2014–2018 survey estimates from the United States Census Bureau [23]. Main commuting mode was operationalized as the total percentage of people that use private transportation to work, such as those that own private vehicles (car, truck, van, motorcycle) or use a private taxi. Private transportation is a proxy for transportation accessibility in each county, as we hypothesize that areas that that lack efficient modes of transportation may not have as fast of COVID-19 exponential growth, and we want to assess whether the effect of density on R0 would differ due to ease of transportation accessibility. Population and land area were obtained from the 2010 census, and density was calculated by population divided by total square km. All census data were extracted using the R package tidycensus [24].

Statistical analysis

We first compared the densities of counties included in the final analytical sample to those that did not have sufficient case counts with a two-sample Wilcoxon test. R0 was estimated from the exponential growth rate method developed by Wallinga and Lipstich [25] and implemented by the R package R0 [26], assuming a generation interval with a gamma distribution of mean 4.7 and standard deviation of 2.9 [27]. We then conducted a cross-sectional analysis using linear mixed models with random intercepts for each state and fixed slopes for the counties to assess the association of population density and R0. The linear mixed models allow the intercept to vary among states, which accounts for non-independence among counties within each state, potentially due to variable resource allocation and differing health systems across states.

We also adjusted for county-level main mode of transportation to work percentage and median household income to control for any potential confounding between the association of density and R0. Therefore, we fit 4 models with R0 as the outcome and the following factors as covariates: Model 1: population density; Model 2: population density and the percent of individuals reporting private transportation as their main mode of transportation to work; Model 3: population density and the percent of individuals reporting private transportation, and median household income; Model 4: population density, percent of individuals reporting private transportation, median household income, and the interaction of private transportation use with population density. The associated linear mixed model equation for model 4 is where β0 = γ00+u0, for the i-th county for the j-th state

Sensitivity analyses

We conducted three sensitivity analyses to address the limitations of our approach and assess the robustness of our results. First, we conducted a sensitivity analysis using death counts from the New York Times [21] to estimate R0 to limit bias due to differential availability of testing by geographic location. We used the same exponential period as the cases, but with a lag of 14 days to account for the delays from symptom onset to deaths of cases [28-30]. Moreover, the analysis of deaths was restricted to counties with greater than 10 deaths and more than 5 daily increases in incident deaths, in order to appropriately estimate R0 in counties with sufficient death counts. The daily death data created curves that approximated exponential growth across the counties (Fig B in S1 Appendix). Second, we excluded counties within a radius of 15 miles, the average commuting miles in the United States [31], from counties with densities greater than the 75th percentile. Removing these adjacent counties would demonstrate the extent of biases due to individuals commuting from surrounding counties to cities. If cases are imported from more densely population (i.e. cities) to less dense counties, we could potentially be biasing our estimates downwards. Lastly, we conducted an analysis excluding influential counties with a Cook’s distance measure over 4/N for each model. Cook’s distance is a commonly used indicator of influence, which measures the extent data points impact the regression parameter estimates [32]. Our sensitivity analysis excluding influential counties ensures that our findings were not driven by these highly influential observations, and that the association holds for other counties.

All analyses was conducted in R version 4.0.0 [33]. The figure and removal of adjacent counties in the sensitivity analyses were done with ArcGIS [34].

Results

The United States has 3,221 counties and county equivalents. When restricting to counties with greater than 25 cases, 1,151 (35.73%) counties were included (Fig 1). The median R0 among the counties was 1.66 (IQR: 1.35–2.11). The median start and end of the exponential growth period was March 25, 2020 and April 9, 2020 respectively, and differed for each county depending on the start of disease transmission. The median density in counties included and not included in the analysis were 53.8 people/km2 (IQR: 21.24–144.05) and 11.3 people/km2 (IQR: 3. 56–23.88) respectively, and the difference was statistically significant (p <0.0001). A population density threshold of approximately 22 people/km2 was needed to sustain an outbreak (Fig 2).

Fig 1

Basic reproductive number (R0) estimates across United States counties.

Larger R0 indicates greater transmission during the initial phase of the outbreak, or the exponential growth period. We restricted calculation of R0 to counties with greater than 25 cases at the end of the exponential growth period (n = 1,151), as R0 cannot be estimated accurately with sparse data and it would be uncertain if the county was experiencing a sustained outbreak with community transmission.

Fig 2

Population density threshold required to establish a sustained outbreak in United States counties.

A population density of approximately 22 people/km2 was needed to sustain an outbreak, which is approximately equal to the lower IQR of the counties with established COVID-19 outbreaks [median = 53.8 population/km2; IQR = 21.24, 144.05], and a slightly less than the upper IQR of the counties with no COVID-19 outbreaks [median = 11.3 population/km2; IQR = 3.56, 23.88]. The grey circles represent the individual county population densities.

An increase in one unit of log population density increased R0 by 0.16 (95% CI = 0.13 to 0.19) (Model 1; Table 1), or the doubling of population density increased the R0 on average by 0.11 (95% CI = 0.09 to 0.13). When adjusted for percent of private transportation and median household income, the association of log population density and R0 remained unchanged (Model 3; Table 1). There was no significant interaction, and the effect of population density on R0 was the same among counties with a larger percentage with private vehicles as their transportation to work, but was no longer significant (Model 4; Table 1). Therefore, model 4 should be not be interpreted, since the interaction does not contribute to the model and only serves to decrease the precision (larger confidence intervals) of the measures of association of interest. R0 decreased by 0.12 (95% CI = -0.02 to -0.04) with an 10% increase in private transportation as the main commute mode, accounting for population density and median household income (Model 3; Table 1).

Table 1

Linear mixed models (random intercept, fixed slope) evaluating the association between population density and basic reproductive number (R0) among United States counties.

	Model 1	Model 2	Model 3	Model 4
Log of population density	0.16*** (0.13–0.19)	0.15*** (0.12, 0.18)	0.15*** (0.12, 0.18)	0.20 (-0.06, 0.47)
Percent of private transportation (10%)		-0.12*** (-0.19, -0.04)	-0.12** (-0.2, -0.04)	-0.08 (-0.27, 0.11)
Median household income ($10,000)			0.00 (-0.03, 0.03)	0.00 (-0.03, 0.03)
Interaction of population density and transportation				-0.01 (-0.04, 0.02)

pvalue

* < 0.05

** < 0.01

*** < 0.001.

Model 1- unadjusted association of population density and R0.

Model 2 –association of population density adjusted for the percent of individuals reporting private transportation.

Model 3 –association of population density, percent of individuals reporting private transportation, and median household income.

Model 4 –association of population density, percent of individuals reporting private transportation, median household income, and the interaction of population density and transportation.

Estimates for each model is a slope (beta) with a null of 0; a positive slope indicates that an increase in the log of population density increases R0 by the beta estimate for the log of population density. The interaction term indicates that the association of population density and R0 differs depending on the percentage of people using the private transportation for work.

In all three sensitivity analyses, population density remained positively associated with R0, demonstrating the robustness of our main analysis. First, death data was used to calculate R0 from 301 counties. The median R0 among the counties that had sufficient death counts was 1.40 (IQR: 1.05–1.78). The unadjusted association between population density and R0 remained consistent (β = 0.18, 95% CI = 0.14 to 0.23) (Model 1a; Table 2), and there were no significant interactions (Model 4a; Table 2). Next, there were 288 counties above the 75th percentile, and 414 counties that were within 15 miles of these counties high-density counties. We removed these 414 counties from the sample, and using the subsample of 737 counties, our findings remained consistent (Table 2). Influential counties were also not driving the association of population density and R0, and our results remained robust (Table 2). For the two sensitivity analyses excluding counties adjacent to highly dense counties and excluding high influence counties, however, the association of private transportation usage and R0 did not remain (Models 2a, 3a, 2b, 3b; Table 2).

Table 2

Sensitivity analysis of linear mixed models (random intercept, fixed slope) using (a) deaths only, (b) removing counties within 15 miles of high-density counties, and (c) removing high influence counties.

	Model 1	Model 2	Model 3	Model 4
Deaths only
Log of population density	0.18*** (0.14, 0.23)	0.14*** (0.09, 0.19)	0.12*** (0.07, 0.17)	0.48* (0.01, 0.94)
Percent of private transportation (10%)		-0.16** (-0.26, -0.06)	-0.15** (-0.24, -0.05)	0.14 (-0.24, 0.53)
Median household income ($10,000)			0.06** (0.03, 0.1)	0.07*** (0.03, 0.11)
Interaction of population density and transportation				-0.04 (-0.09, 0.01)
No counties adjacent to high density counties
Log of population density	0.17*** (0.14, 0.2)	0.17*** (0.13, 0.2)	0.16*** (0.13, 0.2)	0.02 (-0.13, 0.17)
Percent of private transportation (10%)		0.01 (-0.07, 0.1)	-0.01 (-0.12, 0.1)	-0.12 (-0.28, 0.03)
Median household income ($10,000)			0.02 (-0.03, 0.06)	0.01 (-0.03, 0.05)
Interaction of population density and transportation				0.03 (0, 0.06)
Removing high influence counties
Log of population density	0.18*** (0.15, 0.2)	0.18*** (0.16, 0.21)	0.18*** (0.15, 0.2)	0.38* (0.04, 0.72)
Percent of private transportation (10%)		-0.01 (-0.09, 0.07)	-0.02 (-0.1, 0.07)	0.09 (-0.1, 0.29)
Median household income ($10,000)			0 (-0.03, 0.02)	0 (-0.02, 0.03)
Interaction of population density and transportation				-0.02 (-0.06, 0.02)

pvalue

* < 0.05

** < 0.01

*** < 0.001.

Model 1– unadjusted association of population density and R0.

Model 2 –association of population density adjusted for the percent of individuals reporting private transportation.

Model 3 –association of population density, percent of individuals reporting private transportation, and median household income.

Discussion

Our findings show that the basic reproductive number (R0) is associated with population density, even when percent of individuals that use private transportation and median income were accounted for. In these settings, greater population density may potentially facilitate interactions between susceptible and infectious individuals in densely-populated networks, which sustain continued transmission and spread of COVID-19. Moreover, we see that population density continues to have an important impact on disease transmission regardless of transportation accessibility and median income, suggesting that the opportunity for effective contacts are mostly driven by crowding in denser areas, increasing the contact rates necessary for disease spread. However, we did not see that density-dependence is differential across transportation accessibility, as demonstrated by the non-significant interaction of population density and transportation. Our findings are consistent with previous research that have demonstrated a strong relationship between population density and other infectious diseases [4, 6, 7, 9, 35]. In the current SARS-CoV-2 pandemic, recent literature has been conflicting, where some research also suggests a density-dependence of COVID-19 transmission [17, 36] and other measures of the severity of the outbreak [19, 37], while other research suggests that there are other factors that can better explain the pandemic [18, 38]. However, to our knowledge, our results are one of the first to show that population density is an important driver of COVID-19 transmission, even in areas where residents rely more on private modes of transportation. Moreover, even though transmission is less in lower density areas (i.e. rural areas), rural settings may eventually disproportionately be more vulnerable to COVID-19 morbidity and mortality. Individuals in rural areas are generally older, have more underlying conditions, have less access to care, and have fewer ICU beds, ventilators, and facilities needed for severe COVID-19 treatment [39-42]. Further research is needed on the overall burden of COVID-19 across the spectrum of population density.

Geographic estimates of R0 of SARS-CoV-2 need to take into account the specific area’s population density, since the R0 estimate is dependent on both the pathogenicity of the virus as well as environmental influences. In countries where cases are only on the starting to climb, such as countries in Latin America and Southern Africa [1, 43], or there is a resurgence of cases, such as India, Iraq, and Israel [44], area-specific density can assist in predictions of R0, which is important because epidemiological forecasts and predictive models are sensitive to small changes in R0 inputs. Accurate estimation of R0 consequently leads to more precise approximations of the epidemic size, so that governments can appropriately allocate resources and coordinate mitigation strategies. Moreover, as cities and states reopen in the United States, and if there is a second-wave of infections, areas with higher density accessibility will likely have greater SARS-CoV-2 resurgence.

Our study has a number of limitations. While we demonstrate that population density is associated with R0, the estimation of R0 can be biased depending on the data and assumptions adopted. However, our main aim in this analysis was to evaluate the association between population density and R0, and not to accurately estimate R0. Thus, any biases in estimation of R0 due to underlying assumptions would likely be non-differential across counties, and would still yield similar results. In addition, we estimated R0 based on the number of reported cases; therefore, the incidence of COVID-19 across US counties may be underestimated at varying rates due to differential testing. Testing data at the county-level currently do not exist, and we were unable to adjust for the number of tests performed. Confounding of true epidemic growth by increase in testing could also be a potential constraint to the robustness of the analysis. To mitigate this limitation, we included a random intercept term to adjust for state-level effects, and thus differential testing across states were accounted by our model. Differential testing by local governments within states are less likely to strongly impact our findings, as most funding and budgets for COVID-19 is distributed at the state-level [45, 46]. We also conducted a sensitivity analysis using death data which demonstrates the robustness of our findings. Furthermore, we utilized a number of assumptions based on previous findings to calibrate the exponential growth period, which ensured that the virus had taken hold and allowed a sufficient number of days and case counts to estimate exponential growth. There are potential for biases in our method; for example, there is the possibility that some NPIs were introduced in the initial outbreak stage of COVID-19 in some counties; however, if this was the case, then case counts and subsequently R0 would even higher than we calculated, and thus our associations of density and R0 was underestimated. However, we implemented numerous ways to limit the biases. The exponential growth period was restricted to approximately 14 days at the start of the epidemic, where we would expect limited increases in testing and thus would not affect R0 substantially. Moreover, we plotted the calibrated exponential growth curves of all the counties included in our analysis, which gave us reasonable curves that approximated exponential growth for case and death data. Another limitation is that we had to only include counties that had sufficient case data in order to accurately estimate R0; however, if we included all counties, the true association between population density and R0 would likely be greater than what we report in our analysis given our findings that the counties excluded in the analysis had a significantly lower density and expected very low R0 due to lack of cases. Another limitation is that our model also assumes homogenous mixing, which may can be an oversimplification of the heterogeneity in contact patterns within populations [4, 47]. However, previous research has shown that population structure only changes R0 estimates slightly [48], and assumptions of well-mixed populations are valid in small-to-medium spatial scales [17]. Moreover, our method loses spatial granularity in assessing R0 in counties, especially in counties with spatially heterogenous clustering. The aim of our study, however, was to provide a generalizable estimate of the association between population density and R0, in order to appropriately estimate potential for disease transmission, rather than a microspatial estimate that may not be generalizable to other settings. Finally, an important confounder that we were unable to adjust for is the number of importations of SARS-CoV-2 in these counties, as more urbanized areas are more likely to have links with countries and other states where the virus could have originated from. Even so, we still see that once an area is seeded with COVID-19, the growth rate is greater in denser areas during the time period prior to implementation of NPIs.

Conclusions

Counties with greater population density have greater rates of transmission of SARS-CoV-2, likely due to increased contact rates in areas with greater population density. Population density affects the network of contacts necessary for disease transmission, and SARS-CoV-2 R0 estimates need to consider this variability for proper planning and resource allocation, particularly as epidemics newly emerge and old outbreaks resurge.

Population density and basic reproductive number of COVID-19 across United States counties.

(DOCX)

Click here for additional data file.

2 Feb 2021

PONE-D-20-40744

Population density and basic reproductive number of COVID-19 across United States counties

PLOS ONE

Dear Dr. Nichols,

Thank you for submitting your manuscript to PLOS ONE. After careful consideration, we feel that it has merit but does not fully meet PLOS ONE’s publication criteria as it currently stands. Therefore, we invite you to submit a revised version of the manuscript that addresses the points raised during the review process.

Please submit your revised manuscript by Mar 18 2021 11:59PM. If you will need more time than this to complete your revisions, please reply to this message or contact the journal office at plosone@plos.org. When you're ready to submit your revision, log on to https://www.editorialmanager.com/pone/ and select the 'Submissions Needing Revision' folder to locate your manuscript file.

Please include the following items when submitting your revised manuscript:

A rebuttal letter that responds to each point raised by the academic editor and reviewer(s). You should upload this letter as a separate file labeled 'Response to Reviewers'.

A marked-up copy of your manuscript that highlights changes made to the original version. You should upload this as a separate file labeled 'Revised Manuscript with Track Changes'.

An unmarked version of your revised paper without tracked changes. You should upload this as a separate file labeled 'Manuscript'.

If you would like to make changes to your financial disclosure, please include your updated statement in your cover letter. Guidelines for resubmitting your figure files are available below the reviewer comments at the end of this letter.

If applicable, we recommend that you deposit your laboratory protocols in protocols.io to enhance the reproducibility of your results. Protocols.io assigns your protocol its own identifier (DOI) so that it can be cited independently in the future. For instructions see: http://journals.plos.org/plosone/s/submission-guidelines#loc-laboratory-protocols

We look forward to receiving your revised manuscript.

Kind regards,

Martial L Ndeffo Mbah, Ph.D

Academic Editor

PLOS ONE

Additional Editor Comments (if provided):

Reviewer #2 has very constructive comments which should help greatly improve the quality of the manuscript and robustness of its results.

Journal Requirements

When submitting your revision, we need you to address these additional requirements.

1) Please ensure that your manuscript meets PLOS ONE's style requirements, including those for file naming. The PLOS ONE style templates can be found at

https://journals.plos.org/plosone/s/file?id=wjVg/PLOSOne_formatting_sample_main_body.pdf and

https://journals.plos.org/plosone/s/file?id=ba62/PLOSOne_formatting_sample_title_authors_affiliations.pdf

[Note: HTML markup is below. Please do not edit.]

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Partly

**********

2. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: Yes

**********

3. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: Yes

**********

4. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

5. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: This paper is highly informative and timely. My main concern was about differential testing across counties during the exponential growth period, but the authors acknowledge this limitation and do their best to mitigate with random intercept terms. Aside from this, please note a few suggestions/points of clarification below.

• Some brief text on the difference between R0 versus R(t) might be helpful.

• Instead of saying “22 population/km2” perhaps say “22 people/km2?”

• Lines 62 and 63 read strangely, consider rephrasing.

• Consider changing “case” to “cases” in line 76.

• I have trouble understanding line 77. Does it mean to say the exponential growth period began one week prior to the second daily increase in cases?

• Some more explanation about “influential counties” around lines 131-133 would be helpful.

• Consider changing “densely-population” to “densely-populated” in line 196.

Reviewer #2: Comments

This study examined whether cross-county variation in R0 could be explained by variation in population density. The results indicated that an increase in population density increased R0. The method used was linear mixed models with random intercept and fixed slopes. A merit of this study is that it utilized subnational data that contained rich information on local settings. The following are my suggestions for making improvement of the manuscript.

1. The basic reproductive number has many limitations. It must be estimated and varies with the mathematical model and assumptions adopted. R0 was rarely observed in the real world. Authors may consider using other measures of COVID-19 transmission that are more straightforward, such as case number or doubling/halving time, as sensitivity analysis or alternative models. Analysis of those measures may better serve public policy.

2. Literature review should be elaborated, especially studies relating to population density and disease transmission.

3. It is not clear how the period of exponential growth was calibrated, and the reason for assuming the exponential growth period was one week prior to the second daily increases in cases.

4. It is not clear how mixed linear models fit into the cross-sectional data structure, as those models are usually for multilevel/hierarchical, longitudinal, or correlated data. Does each county have repeated measurements? Are measurements made on clusters?

5. Figures seem to be poorly presented.

6. It was possible that some NPIs were introduced in the initial outbreak stage of COVID-19 in some counties, such as public information campaigns. This would make the study assumption invalid.

7. Authors are advised to provide explanations for insignificance of the coefficient estimate of population density in Model 4.

**********

6. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

[NOTE: If reviewer comments were submitted as an attachment file, they will be attached to this email and accessible via the submission site. Please log into your account, locate the manuscript record, and check for the action link "View Attachments". If this link does not appear, there are no attachment files.]

While revising your submission, please upload your figure files to the Preflight Analysis and Conversion Engine (PACE) digital diagnostic tool, https://pacev2.apexcovantage.com/. PACE helps ensure that figures meet PLOS requirements. To use PACE, you must first register as a user. Registration is free. Then, login and navigate to the UPLOAD tab, where you will find detailed instructions on how to use the tool. If you encounter any issues or have any questions when using PACE, please email PLOS at figures@plos.org. Please note that Supporting Information files do not need this step.

23 Feb 2021

Response to Reviewers

Reviewer #1

This paper is highly informative and timely. My main concern was about differential testing across counties during the exponential growth period, but the authors acknowledge this limitation and do their best to mitigate with random intercept terms. Aside from this, please note a few suggestions/points of clarification below.

Comment 1

Some brief text on the difference between R0 versus R(t) might be helpful.

Response 1

We thank the reviewer for the positive and helpful assessment of the manuscript. We agree with the reviewer that briefly describing the difference between these two measures would add additional context. On page 3, paragraph 2, we added the following text:

“This is in contrast with the time-varying reproductive number (Rt), defined as the number of people in a population who were infected by an infectious individual at a given point in time, which reflects the changing levels of immunity in the population and the impact of control measures limiting transmission.”

Comment 2

• Instead of saying “22 population/km2” perhaps say “22 people/km2?”

Response 2

We thank the reviewer for the suggestion, and have changed the text accordingly. We have also changed the text of Figure 2 for consistency.

Comment 3

• Lines 62 and 63 read strangely, consider rephrasing.

Response 3

We thank the reviewer for the thoughtful suggestion. We rephrased both this sentence and the sentence before it (now lines 64 to 66), in order to improve the transition from the previous sentence.

“During the initial phase of the outbreak, or the exponential growth period, we hypothesize that spatial heterogeneity in R0 occurs primarily due to geographic variability in contact rates, since transmission probability and population size remain constant across settings. During this time frame, transmission probabilities across localities are equivalent, because exponential growth occurs prior to the implementation of NPIs which affect likelihood of transmission during contact.”

Comment 4

• Consider changing “case” to “cases” in line 76.

Response 4

We have modified the text from “case” to “cases” in line 76.

Comment 5

• I have trouble understanding line 77. Does it mean to say the exponential growth period began one week prior to the second daily increase in cases?

Response 5

We thank the reviewer for realizing our typographical error. Yes, the reviewer is right, we meant to say that “we assumed that the exponential growth period began one week prior to the second daily increase in cases”. We have changed the text accordingly (now line 82). In the same paragraph, we state that “The algorithm ensured that the virus had taken hold in the area and allowed a sufficient number of days to estimate the exponential growth rate, as R0 cannot be estimated accurately with sparse data, since it would be uncertain if the county was experiencing a sustained outbreak with community transmission.”

Comment 6

• Some more explanation about “influential counties” around lines 131-133 would be helpful.

Response 6

We thank the reviewer, and have added additional information and a citation on “influential counties”, Cook’s D, and influence (now lines 142 to 146).

“Lastly, we conducted an analysis excluding influential counties with a Cook’s distance measure over 4/N for each model. Cook’s distance is a commonly used indicator of influence, which measures the extent data points impact the regression parameter estimates [31]. Our sensitivity analysis excluding influential counties ensures that our findings were not driven by these highly influential observations, and that the association holds for other counties.”

Comment 7

• Consider changing “densely-population” to “densely-populated” in line 196.

Response 7

We thank the review for noticing this typographical error, and have changed “densely-population” to “densely-populated” in line 196 (now line 213).

Reviewer #2

This study examined whether cross-county variation in R0 could be explained by variation in population density. The results indicated that an increase in population density increased R0. The method used was linear mixed models with random intercept and fixed slopes. A merit of this study is that it utilized subnational data that contained rich information on local settings. The following are my suggestions for making improvement of the manuscript.

Comment 1

1. The basic reproductive number has many limitations. It must be estimated and varies with the mathematical model and assumptions adopted. R0 was rarely observed in the real world. Authors may consider using other measures of COVID-19 transmission that are more straightforward, such as case number or doubling/halving time, as sensitivity analysis or alternative models. Analysis of those measures may better serve public policy.

Response 1

We thank the reviewer for the positive assessment of our manuscript and suggestions for improvement. We agree that the estimation of the basic reproductive number can be biased depending on the data and assumptions adopted. However, our main aim in this analysis was to evaluate the association between population density and R0, and not to accurately estimate R0. Thus, any biases in estimation of R0 due these underlying assumptions would likely be non-differential across counties, and would yield similar results to our findings. We further outline these limitations in the Discussion section (p. 242)

“While we demonstrate that population density is associated with R0, the estimation of R0 can be biased depending on the data and assumptions adopted. However, our main aim in this analysis was to evaluate the association between population density and R0, and not to accurately estimate R0. Thus, any biases in estimation of R0 due to underlying assumptions would likely be non-differential across counties, and would still yield similar results. In addition, we estimated R0 based on the number of reported cases; therefore, the incidence of COVID-19 across US counties may be underestimated at varying rates due to differential testing. Testing data at the county-level currently do not exist, and we were unable to adjust for the number of tests performed. Confounding of true epidemic growth by increase in testing could also be a potential constraint to the robustness of the analysis.”

We also outlined numerous analytics methods and sensitivity analyses to mitigate these limitations (p. 250).

“To mitigate this limitation, we included a random intercept term to adjust for state-level effects, and thus differential testing across states were accounted by our model. Differential testing by local governments within states are less likely to strongly impact our findings, as most funding and budgets for COVID-19 is distributed at the state-level [45, 46]. We also conducted a sensitivity analysis using death data which demonstrates the robustness of our findings.”

As our main aim in this study is the association of population density and transmission of SARS-CoV-2 as measured by R0, we believe that using a different a measure of COVID-19 severity such as case counts or doubling/halving time would ask a different research question, and measures may still not be free of biases. Evaluating the association of population density and other COVID-19 measures is important, but is not a focus of this paper. In both the Introduction (p. 3) and the Discussion section (p. 12), we reiterate why we chose to focus on R0 as the main outcome.

“During the initial phase of the outbreak, or the exponential growth period, we hypothesize that spatial heterogeneity in R0 occurs primarily due to geographic variability in contact rates, since transmission probability and population size remain constant across settings.”

“Geographic estimates of R0 of SARS-CoV-2 need to take into account the specific area’s population density, since the R0 estimate is dependent on both the pathogenicity of the virus as well as environmental influences… Accurate estimation of R0 consequently leads to more precise approximations of the epidemic size, so that governments can appropriately allocate resources and coordinate mitigation strategies.”

Comment 2

2. Literature review should be elaborated, especially studies relating to population density and disease transmission.

Response 2

We appreciate the helpful suggestion by the reviewer. We have added additional literature review outlining studies relating population density and disease transmission of other infectious diseases. Moreover, we have updated the literature review on studies specifically relating population density and COVID-19 transmission (p. 12).

Our findings are consistent with previous research that have demonstrated a strong relationship between population density and other infectious diseases [4, 6, 7, 9, 35]. In the current SARS-CoV-2 pandemic, recent literature has been conflicting, where some research also suggests a density-dependence of COVID-19 transmission [17, 36] and other measures of the severity of the outbreak [19, 37], while other research suggests that there are other factors that can better explain the pandemic [18, 38]. However, to our knowledge, our results are one of the first to show that population density is an important driver of COVID-19 transmission, even in areas where residents rely more on private modes of transportation.”

Comment 3

3. It is not clear how the period of exponential growth was calibrated, and the reason for assuming the exponential growth period was one week prior to the second daily increases in cases.

Response 3

We agree with the reviewer that clarifying our method of determining the timing of the exponential growth period is very important. Our method actually assumed that the exponential growth period began one week prior to the second daily increase in cases. We have modified the text accordingly. In the same paragraph, we state the reasoning behind this assumption.

“The algorithm ensured that the virus had taken hold in the area and allowed a sufficient number of days to estimate the exponential growth rate, as R0 cannot be estimated accurately with sparse data, since it would be uncertain if the county was experiencing a sustained outbreak with community transmission.”

We understand that this method may have biases, and we added additional information in the Discussion section as a limitation (p.13)

“Furthermore, we utilized a number of assumptions based on previous findings to calibrate the exponential growth period, which ensured that the virus had taken hold and allowed a sufficient number of days and case counts to estimate exponential growth. There are potential for biases in our method… However, we implemented numerous ways to limit the biases. The exponential growth period was restricted to approximately 14 days at the start of the epidemic, where we would expect limited increases in testing and thus would not affect R0 substantially. Moreover, we plotted the calibrated exponential growth curves of all the counties included in our analysis, which gave us reasonable curves that approximated exponential growth for case and death data.”

Despite the limitations of our method, we believe there is value in using a relatively straightforward method to assess exponential growth period, that can potentially be used in other research studies of COVID-19 or future pandemics.

Comment 4

4. It is not clear how mixed linear models fit into the cross-sectional data structure, as those models are usually for multilevel/hierarchical, longitudinal, or correlated data. Does each county have repeated measurements? Are measurements made on clusters?

Response 4

We thank the reviewer for the valuable suggestion, and agree that clarification on the model specification in needed. As the reviewer correctly stated, linear mixed models are usually used for correlated data (non-independence), which can arise from a hierarchical structure. Our data has a hierarchical structure, as COVID-19 case numbers among counties within states are non-independent, potentially due to variable resource allocation and differing health systems across states. Thus, including a random intercept term for the states, and a fixed slope for the counties in the linear mixed models would be able to account for state-level correlation among counties. We have added text in the Methods section to further clarify this (p. 5):

“We then conducted a cross-sectional analysis using linear mixed models with random intercepts for each state and fixed slopes for the counties to assess the association of population density and R0. The linear mixed models allow the intercept to vary among states, which accounts for non-independence among counties within each state, potentially due to variable resource allocation and differing health systems across states.”

Moreover, in order to clarify our model specification in the manuscript, we have included the linear mixed model equation (p. 6)

“The associated linear mixed model equation for model 4 is

Y_ij=β_0j+β_1j l〖n(density1 unit increase)〗_ij+β_2j 〖% private transportation〗_ij+β_3j 〖median household income〗_ij +β_4j 〖% private transportation〗_ij*〖median household income〗_ij+e_ij

where β_0j=γ_00+u_0j, for the i-th county for the j-th state”

Comment 5

5. Figures seem to be poorly presented.

Response 5

We have changed the color scale of Figure 1 to improve the distinction between the R0 values. We also excluded Alaska and Hawaii for presentability, as is commonly done in spatial maps of the United States.

Moreover, we updated Figure 2 to add in grey points to represent the individual county population densities, which add valuable information to the Figure and increase presentability of the figure.

We have run our figures through the PACE tool, which ensures that all our figures meet PLOS One’s technical requirements.

Comment 6

6. It was possible that some NPIs were introduced in the initial outbreak stage of COVID-19 in some counties, such as public information campaigns. This would make the study assumption invalid.

Response 6

We thank the reviewer for noting the need to clarify this. There is the possibility that some NPIs were introduced in the initial outbreak stage of COVID-19 in some counties; however, if this was the case, then case counts and subsequently R0 would even higher than we calculated, and thus our associations of density and R0 was underestimated. Moreover, we plotted the calibrated exponential growth curves of all the counties included in our analysis, which gave us reasonable curves that approximated exponential growth for case and death data. This gives further validity to our calibration method of the exponential growth period (S1_Appendix).

We have added this as a we have added this as a potential limitation in our study, as well as the measures we put in place to limit these potential biases (p. 13).

“There are potential for biases in our method; for example, there is the possibility that some NPIs were introduced in the initial outbreak stage of COVID-19 in some counties; however, if this was the case, then case counts and subsequently R0 would even higher than we calculated, and thus our associations of density and R0 was underestimated. However, we implemented numerous ways to limit the biases. The exponential growth period was restricted to approximately 14 days at the start of the epidemic, where we would expect limited increases in testing and thus would not affect R0 substantially. Moreover, we plotted the calibrated exponential growth curves of all the counties included in our analysis, which gave us reasonable curves that approximated exponential growth for case and death data.”

Comment 7

7. Authors are advised to provide explanations for insignificance of the coefficient estimate of population density in Model 4.

Response 7

We thank the reviewer for the valuable suggestion for the insignificance of the coefficient estimate of population density in Model 4. Since the interaction between population density and transportation was non-significant, this demonstrates that density-dependence is not differential across transportation accessibility, and thus Model 4 results should not be interpreted. We attempted to make this clearer in the Discussion section (p. 11, last paragraph)

“However, we did not see that density-dependence is differential across transportation accessibility, as demonstrated by the non-significant interaction of population density and transportation. Therefore, model 4 should be not be interpreted, since the interaction does not contribute to the model and only serves to decrease the precision (larger confidence intervals) of the measures of association of interest.”

Submitted filename: Nichols_USDensityCOVID_ResponsetoReviewers.docx

Click here for additional data file.

16 Mar 2021

Population density and basic reproductive number of COVID-19 across United States counties

PONE-D-20-40744R1

Dear Dr. Nichols,

We’re pleased to inform you that your manuscript has been judged scientifically suitable for publication and will be formally accepted for publication once it meets all outstanding technical requirements.

Within one week, you’ll receive an e-mail detailing the required amendments. When these have been addressed, you’ll receive a formal acceptance letter and your manuscript will be scheduled for publication.

An invoice for payment will follow shortly after the formal acceptance. To ensure an efficient process, please log into Editorial Manager at http://www.editorialmanager.com/pone/, click the 'Update My Information' link at the top of the page, and double check that your user information is up-to-date. If you have any billing related questions, please contact our Author Billing department directly at authorbilling@plos.org.

If your institution or institutions have a press office, please notify them about your upcoming paper to help maximize its impact. If they’ll be preparing press materials, please inform our press team as soon as possible -- no later than 48 hours after receiving the formal acceptance. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information, please contact onepress@plos.org.

Kind regards,

Martial L Ndeffo Mbah, Ph.D

Academic Editor

PLOS ONE

Additional Editor Comments (optional):

Reviewers' comments:

Reviewer's Responses to Questions

Comments to the Author

1. If the authors have adequately addressed your comments raised in a previous round of review and you feel that this manuscript is now acceptable for publication, you may indicate that here to bypass the “Comments to the Author” section, enter your conflict of interest statement in the “Confidential to Editor” section, and submit your "Accept" recommendation.

Reviewer #1: All comments have been addressed

Reviewer #2: All comments have been addressed

**********

2. Is the manuscript technically sound, and do the data support the conclusions?

The manuscript must describe a technically sound piece of scientific research with data that supports the conclusions. Experiments must have been conducted rigorously, with appropriate controls, replication, and sample sizes. The conclusions must be drawn appropriately based on the data presented.

Reviewer #1: Yes

Reviewer #2: Yes

**********

3. Has the statistical analysis been performed appropriately and rigorously?

Reviewer #1: Yes

Reviewer #2: Yes

**********

4. Have the authors made all data underlying the findings in their manuscript fully available?

The PLOS Data policy requires authors to make all data underlying the findings described in their manuscript fully available without restriction, with rare exception (please refer to the Data Availability Statement in the manuscript PDF file). The data should be provided as part of the manuscript or its supporting information, or deposited to a public repository. For example, in addition to summary statistics, the data points behind means, medians and variance measures should be available. If there are restrictions on publicly sharing data—e.g. participant privacy or use of data from a third party—those must be specified.

Reviewer #1: Yes

Reviewer #2: Yes

**********

5. Is the manuscript presented in an intelligible fashion and written in standard English?

PLOS ONE does not copyedit accepted manuscripts, so the language in submitted articles must be clear, correct, and unambiguous. Any typographical or grammatical errors should be corrected at revision, so please note any specific errors here.

Reviewer #1: Yes

Reviewer #2: Yes

**********

6. Review Comments to the Author

Please use the space provided to explain your answers to the questions above. You may also include additional comments for the author, including concerns about dual publication, research ethics, or publication ethics. (Please upload your review as an attachment if it exceeds 20,000 characters)

Reviewer #1: Thank you for addressing my comments. While R0 is a complex measure to calculate, and is easily misrepresented, misinterpreted, and misapplied, it is a fundamental metric for the study of infectious disease. This research, therefore, adds to the scientific literature on COVID-19 in a meaningful way.

Reviewer #2: The authors have addressed issues raised in the reviewer report. I hope the readers will benefit from this work.

**********

7. PLOS authors have the option to publish the peer review history of their article (what does this mean?). If published, this will include your full peer review and any attached files.

If you choose “no”, your identity will remain anonymous but your review may still be made public.

Do you want your identity to be public for this peer review? For information about this choice, including consent withdrawal, please see our Privacy Policy.

Reviewer #1: No

Reviewer #2: No

31 Mar 2021

PONE-D-20-40744R1

Population density and basic reproductive number of COVID-19 across United States counties

Dear Dr. Nichols:

I'm pleased to inform you that your manuscript has been deemed suitable for publication in PLOS ONE. Congratulations! Your manuscript is now with our production department.

If your institution or institutions have a press office, please let them know about your upcoming paper now to help maximize its impact. If they'll be preparing press materials, please inform our press team within the next 48 hours. Your manuscript will remain under strict press embargo until 2 pm Eastern Time on the date of publication. For more information please contact onepress@plos.org.

If we can help with anything else, please email us at plosone@plos.org.

Thank you for submitting your work to PLOS ONE and supporting open access.

Kind regards,

PLOS ONE Editorial Office Staff

on behalf of

Dr. Martial L Ndeffo Mbah

Academic Editor

PLOS ONE