key: cord-0804958-z6b5ci3r
authors: Ben-Nun, M.; Riley, P.; Turtle, J. A.; Riley, S.
title: Consistent Pattern of Epidemic Slowing Across Many Geographies Led to Longer, Flatter Initial Waves of the COVID-19 Pandemic
date: 2022-04-01
journal: nan
DOI: 10.1101/2022.03.31.22273267
sha: c0dde3318b8ae25490ac62f63af7134c5bc436f5
doc_id: 804958
cord_uid: z6b5ci3r

To define appropriate planning scenarios for future pandemics of respiratory pathogens, it is important to understand the initial transmission dynamics of COVID-19 during 2020. Here, we fit an age-stratified compartmental model with a flexible underlying transmission term to daily COVID-19 death data from states in the contiguous U.S. and to national and sub-national data from around the world. The daily death data of the first months of the COVID-19 pandemic was categorized into one of four main types: "spring single-peak profile", "summer single-peak profile", "spring/summer two-peak profile" and "broad with shoulder profile". We estimated a reproduction number R as a function of calendar time tc and as a function of time since the first death reported in that population (local pandemic time, tp). Contrary to the multiple categories and range of magnitudes in death incidence profiles, the R(tp) profiles were much more homogeneous. We find that in both the contiguous U.S. and globally, the initial value of both R(tc) and R(tp) was substantial: at or above two. However, during the early months, pandemic time R(tp) decreased exponentially to a value that hovered around one. This decrease was accompanied by a reduction in the variance of R(tp). For calendar time R(tc), the decrease in magnitude was slower and non-exponential, with a smaller reduction in variance. Intriguingly, similar trends of exponential decrease and reduced variance were not observed in raw death data. Our findings suggest that the combination of specific government responses and spontaneous changes in behaviour ensured that transmissibility dropped, rather than remaining constant, during the initial phases of a pandemic. Future pandemic planning scenarios should be based on models that assume similar decreases in transmissibility, which lead to longer epidemics with lower peaks when compared with models based on constant transmissibility.

The roll out of effective vaccines [1, 2] and the emergence of more transmissible [3] [4] [5] 2 and antigenically distinct [6] lineages of SARS-Cov-2 virus [7] marked the end of the 3 global first wave of the COVID-19 pandemic. Over the first year, the COVID-19 4 pandemic has negatively impacted the health and well being of almost every population 5 around the world. In the absence of an effective vaccine, most countries implemented 6 non-pharmaceutical interventions (NPIs), including travel restrictions, school and work 7 closures, social distancing, contact tracing, quarantining and mask requirements [8] [9] [10] . 8 However, the degree of compliance of the population with these measures, and their 9 effectiveness, varied greatly from one setting to another and, even now, is not fully 10 understood [11] [12] [13] [14] [15] . 11 The transmission of the SARS-Cov-2 virus is often quantified using the time-varying 12 reproduction number, R(t), which represents the mean number of secondary cases that 13 a single index case will infect. Many studies have focused on estimating the impact of 14 different interventions on R(t) under the implicit assumption that interventions are 15 similar between different populations [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] . However, the analytical approach in these 16 studies conditions on the assumption that the interventions as measured are the main 17 drivers of changes in R(t), with transmissibility assumed to be constant otherwise. 18 In this study, we develop an age-stratified compartmental model [26] with a 19 smoothly varying reproduction number and use it to study the epidemic, from January 20 to October of 2020, in 49 jurisdictions in the contiguous U.S. and 89 locations globally. 21 Our model allows for multiple values of R(t) and is an extension of our previous work 22 which used a smoothly varying two-value functional form [28] [29] [30] . We consider the idea 23 that, whereas epidemic patterns vary from one location to another, trends in R(t) are 24 consistent if measured relative to "pandemic time", defined as the time elapsed since 25 the first reported death in a location. Trends in pandemic time R(t p ) are compared to 26 those of calendar time R(t c ). Unlike earlier studies, our model fits the inferred daily 27 death and we use a Markov Chain Monte Carlo (MCMC) procedure [31] to fit R(t) to 28 an increasing number of pandemic days. For each fitting time-window, we analyze the 29 value of R(t) for the prior two weeks and discuss the results for the U.S. and globally 30 without attempting to correlate any changes in R(t p ) with specific NPIs. 31 

Data Selection 33 Many data streams can be used as a measure for the spread of the COVID-19 34 pandemic [11] [12] [13] ; however, daily confirmed number of cases, confirmed 35 hospitalizations [1] and cumulative deaths are the most commonly used, and, in 36 particular, these datasets as reported by the Johns Hopkins University (JHU) team [32] . 37 Because of likely biases in the confirmed number of cases (e.g., the large change in test 38 availability over time and the transition to rapid home testing) and the limited 39 availability of hospitalization data, we use the reported confirmed deaths as the most 40 March 31, 2022 2/24 accurate and least biased measure of the pandemic. We note that the death data are 41 also imperfect and may underestimate the true toll of the pandemic (Viglione et al. [33] given an equal weight of one. 48 The data were split into two groups: US-states and global. The US dataset included 49 all 48 contiguous states and the District of Columbia. Global data locations were chosen 50 from both country level and administrative level-one divisions 51 (state/province/territory/etc). Sub-country locations include those that appeared in the 52 JHU database on April 31, 2020 (Canada, Australia), but exclude European island were taken from census data [36, 37] , and converted to the same combined-sex and 61 decadal ages format.

Model and Fits 63 We use an age-stratified SE[I] 4 RX compartmental model (where X refers to death and 64 four levels of severity for the infectious compartment are included: asymptomatic, flu 65 like, mild and severe) to fit the inferred daily reported death for each location 66 separately. Given the relatively low level of mobility between states/countries during 67 the first nine months of the pandemic, this is a reasonable approximation. The 68 age-specific parameters of the model are taken from [38] , and the age specific contact 69 matrix is directly derived from the work of Walker et al. [39] using the"squire" R 70 package available on github: https://github.com/mrc-ide/squire.

In previous studies on influenza and influenza-like-illness [28-30] we used a smoothly 72 varying two-value functional form to describe R(t). Here we extend this model to an 73 arbitrary number of values:

(a n − a n−1 ) tanh

This produces a smooth curve where at roughly time t n−1 , the value of R(t) 75 transitions to a n with an approximate transition time of a ≈ 2L days.

For each location, we determine the joint posterior distribution for the model 77 parameters by the fitting the inferred daily reported death using an adaptive step size 78 MCMC [31] procedure with 10 6 steps. Only the parameters that govern the time We use a simulate-and-recover procedure to validate the model. A known synthetic 88 profile for R(t) is used to generate synthetic daily death incidence data. The model is 89 fit to synthetic data and the recovered incidence and R(t profile is compared with the 90 known input (S1 Fig in Supporting Information). The model is able to recover a large 91 variety of synthetic profiles with high accuracy.

To investigate how R(t) evolved during the course of the pandemic we introduce the 93 concept of "local pandemic time", defined as the number of days elapsed since the first 94 reported death in a location. For each location, the fitting procedure is repeated using 95 an increasing number of local pandemic days: we start with 30 days of data since the 96 first reported death and increase it to 45, 60, 75, 90 and 105. We also fit and report the 97 reproduction number as a function of calendar time R(t c ). Starting with data only until 98 mid-April 2020 and increasing it in five increments of 15 days. This analysis was first 

The U.S. outbreak was first detected in the state of Washington in late February [41] .

The next six months of the pandemic can be thought of as a sequence of four state-level 110 archetypal epidemic curves of reported deaths (Fig 1) . The first appears as a "spring 111 single-peak profile". This north-east wave spreads from New York and New Jersey to 112 neighboring states (e.g., Connecticut and Massachusetts), and to the entire north-east 113 corridor (e.g., Rhode Island, New Hampshire, Virginia, the District of Columbia, 114 Maryland, Delaware and Pennsylvania). Overall, during this period, most states in HHS 115 regions 1-3 exhibit the north-east profile with a large peak in daily deaths. We describe 116 the second typical shape we observe as a "summer single-peak profile". It is observed in 117 a few states that avoided a spring peak but saw their first peak in the summer (e.g.

South Carolina, Tennessee, Florida, Texas, Arizona, Arkansas and Idaho). The third 119 typical shape exhibits two peaks, in both the spring and the summer. For example 120 Georgia, Louisiana and Nevada. The fourth typical profile, "broad with shoulder 121 profile", comes from states that do not exhibit a clear sharp peak (e.g North Carolina, 122 Kentucky, Oregon, Utah, and California) but rather a broad long shoulder-like profile 123 with increases in deaths appearing only in the summer/early fall. Finally, we note that 124 sparsely populated states (such as North and South Dakota), which largely avoided the 125 spring and summer peaks, did not start to show increases in deaths until September. represent the different data profiles described above (summer peak, extended shoulder 129 with summer peak, spring peak, and spring and summer peaks). Similarly, the global 130 locations were chosen to be representative of four continents and all data profiles. The 131 different epidemic profiles lead to different patterns in the time-varying transmissibility 132 when viewed in calendar time. For example, the profile for Pennsylvania starts at R ≈ 3 133 before dropping to just below 1 in May. It maintained that value up to July before . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted April 1, 2022. jurisdictions and most of the 15 locations (exhibiting very different dynamics) the initial 144 value of R(t) is high (above 3 in Italy, Iran and Sweden) and it decreases to either below 145 or around 1 (dashed grey line). In all cases the model produces a smooth result for R(t), 146 and with its' flexible form we are able to fit the variety of profiles described before: 147 early spring only peak (e.g Italy and Sweden), summer only peak (Egypt), spring and The distribution of deaths per capita across states in the continental U.S. was stable 153 for the first half of the study. In contrast, the distribution of R numbers declined 154 substantially during the same period (Fig 3) . We refitted overlapping subsets of the 155 daily death data for the contiguous U.S. starting with data up to mid-April 2020 and 156 increasing in five increments of 15 days (Fig 3 panel(a) ). We found initial values of 157 R(t c ) that were large (mean/median of 2.46/2.33), started to decrease only in May and 158 approached one in June (S1 Table of Supporting Information). We then defined 159 pandemic time t p -as an alternative for calendar time -as the number of days elapsed 160 since the first reported death in a given population and estimated R(t p ) the 161 reproduction number as a function of pandemic time (a time shift for each population, 162 Fig 3 panel (b) ). Although patterns of R in calendar and pandemic time were similar, it 163 dropped more quickly in pandemic time and had a lower variance (S2 Table of . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted April 1, 2022. ; Table of SI) . In panels (b) and (d), the mean date associated with the local pandemic times is indicated above each set of results. In all four panels, to increase readability, a jitter is applied to the displayed data points. difference between R in pandemic time and calendar time was also apparent when we 165 fitted an exponential decay model to both sets of estimates: the model was a much 166 better explanation for the temporal pattern of R in pandemic time than it was for 167 calendar time (Tables S3 Table and S4 Table of Tables S7 Table and S8 Table of SI for the exponential fit parameters.

types present that were observed in the contiguous U.S. (see above). In the Americas, 177 the spring single peak profile was observed only in Canada and Ecuador whereas most 178 states in Central and South America showed either the summer or (late summer peak) 179 single peak. The two countries with a large number of deaths in South America (Brazil 180 and Mexico) also showed a large wide peak that extended over more than three months. 181 The data for Asia-Pacific showed an increase in death in most places only in June (e.g., 182 the single summer peak profile). However, Iran was a notable exception, showing the 183 spring-summer double peak profile (with the first peak having already occurred in 184 April) followed by a third resurgence in September. In Australia (and particularly the 185 state of Victoria) we see the single summer peak profile (with the peak in death 186 occurring in August/early September). Examples for the summer single peak profile are 187 Bangladesh and Saudi Arabia, which largely avoided any excess death during the spring. 188 The data for Europe showed clear regional grouping with Italy and Spain leading the 189 spring single peak followed closely by most of the larger European countries.

For both calendar and pandemic time we found similar trends to those observed in 191 the contiguous U.S. (Fig 4 and S5 the initial values in mid-April were above two and had a large variance (mean/median 195 and standard deviation of 2.48/2.30 and 0.76). As with the U.S. populations, the 196 decline in magnitude and variance of the calendar R was slower than the exponential 197 decline of R in pandemic time R, which hovered around a value of one for many days. 198 However, the rate of decline in transmissibility in pandemic time was slower than that 199 for the contiguous U.S. (compare the top right panels of 3 and 4, and S2 . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted April 1, 2022. ; https://doi.org/10.1101/2022.03.31.22273267 doi: medRxiv preprint Table, and S4 Table and S8 Table) . We have verified that our conclusions for the 201 calendar and pandemic R(t) do not change when more global locations are included in 202 the analysis (S3 Fig in Supporting information) . For the contiguous U.S., here, too, we 203 found that the daily number of deaths per capita remained far more stable than 204 transmissibility when viewed in both calendar and pandemic times (Fig 4 lower panels) . 205 A histogram and heat map representation of the pandemic time evolution of R is investigated overall trends in R. This, we believe, may help guide higher-level planning 221 for the next pandemic caused by a severe respiratory pathogen. We found a common 222 pattern of a reduction in transmissibility during the initial period, rather than constant 223 transmissibility, of which there were no obvious examples. Therefore, regardless of the 224 specific intervention plans in place for individual countries, we suggest that initial 225 planning for future similar pandemics should not be based on assumptions of prolonged 226 constant transmissibility driving a rapid peak and the development of population 227 immunity.

Our study relies on a number of potentially important assumptions, approximations 229 and limitations. First, we chose the daily deaths as the dataset to fit to the model, 230 inferring it from cumulative confirmed deaths, as opposed to using other data such as 231 cases , which have been associated with known and considerably larger biases. However, 232 while the deaths are likely a more reliable measure of the pandemic than cases they are 233 far from perfect. Different states and countries use different criteria when registering 234 deaths (for example, some report only confirmed deaths while others report both 235 probable and confirmed deaths) and have different delays in reporting. Additionally, the 236 reported numbers have been shown to be lower than the true toll of the pandemic (see 237 e.g [33, 34, [46] [47] [48] [49] ). local travel were significantly disrupted in 2020, they were responsible for the initial 241 spread of COVID-19, and continue to play a role during latter periods. In this work, we 242 ignored the initial introduction of cases to a location and focused on the dynamics of 243 the virus following its introduction to a population. Also, we assumed the same quality 244 of care over time at all locations. In practical terms, in the model, we assume that the 245 probability of death from COVID-19 depended only on age during this period prior to 246 wide spread vaccination. In reality, treatments for severe patients have improved over . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted April 1, 2022. ; https://doi.org/10.1101/2022.03.31.22273267 doi: medRxiv preprint and likely had little impact in most places.) 251 We note that a large number of studies have used available databases of 252 interventions [8] [9] [10] , coupled with statistical methods, to estimate the impact of 253 different interventions on the time-varying reproduction number [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [42] [43] [44] [45] . In an 254 extensive study, Liu et al. [16] applied statistical methods to 13 categories of NPIs using 255 confirmed cases data from January to June of 2020 from 130 countries. They concluded 256 that there is a strong association between school closure and internal movement 257 restrictions on the reduction in R eff (t), which was taken from EpiForecast [51] and 258 based on the Cori method [52] . Another major statistical study by Haug et al. [15] used 259 a large database of interventions [8] and multiple statistical methods to quantify 6,068 260 NPIs in 79 locations on R eff (t), concluding that less intrusive and costly NPIs can be as 261 effective as more intrusive ones. A Bayesian hierarchical model was used by Brauner et 262 al.

[21] to study 8 NPIs in 41 countries using national case and death counts al. [44] studied the effect of social distancing measures using early (March-April 2020) 280 case-count data from British Columbia and five other jurisdictions. They developed an 281 SEIR model with physical distancing compartments that contribute a reduced amount 282 to the force of infection. In this model, the fractional reduction to the force of infection 283 is due to increase in physical distancing and it is allowed to linearly change, over the 284 course of one week, from 1 (no physical distancing) to a final value which they estimate 285 using a Bayesian statistical model. Duque et al. [45] describe an age and risk stratified 286 SEIR-style model with a transmission parameter that is reduced during 287 stay-at-home/work-safe-order time periods and use it as part of an optimization 288 framework that seeks to minimize lockdown days while preventing healthcare surges.

Weitz et al.

[53] developed an SEIR model with "short-term awareness of risk", that 290 depends on the death rate, and showed that this feedback mechanism can result in 291 highly asymmetric epidemic curves and death time-series with extended plateau periods. 292 However, this model was not used to fit any data. Implicit in all of these studies is the 293 assumption that interventions are similar between different populations and that they 294 can be correlated to changes in R(t). The main differences between our work and these 295 studies are that we did not attempt to correlate changes in R(t) with any specific NPIs 296 and we used an age-stratified compartmental model with a flexible, smoothly varying, 297 time-dependent reproduction number.

Whereas the conclusions of our study do not directly contradict any of these 299 previously published results [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [42] [43] [44] [45] , our approach of using the dynamics and September 2020 [7] , with the former estimated to be greater than 50% more 

General properties of the initial wave of the COVID-19 pandemic will likely remain a 335 topic of considerable interest for many years [66] . In future studies, we plan to extend 336 the framework developed here. Should other datasets (e.g., case counts and/or 337 hospitalizations) become less biased and more available, we will incorporate them into 338 the objective function we fit. Our model can also be made more flexible by including 339 differences in quality of healthcare over time and location, which will allow for coupling 340 between geographic regions [29] . Our previous work on forecasting ILI in the 341 U.S. [28] [29] [30] highlights the important role that spatial coupling can play in respiratory 342 disease transmission, and which we anticipate will become increasingly more important 343 as travel restrictions have been relaxed and movement between states, countries, and 344 continents significantly increases. . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity.

The copyright holder for this preprint this version posted April 1, 2022. . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted April 1, 2022. ; https://doi.org/10.1101/2022.03.31.22273267 doi: medRxiv preprint S5 . CC-BY-ND 4.0 International license It is made available under a is the author/funder, who has granted medRxiv a license to display the preprint in perpetuity. (which was not certified by peer review)

The copyright holder for this preprint this version posted April 1, 2022. ; https://doi.org/10.1101/2022.03.31.22273267 doi: medRxiv preprint

COVID-19) Vaccinations

Estimated transmissibility and impact of SARS-CoV-2 lineage B.1.1.7 in 357

Assessing transmissibility of SARS-CoV-2 lineage B.1.1.7 in England

Emergence of a new SARS-CoV-2 variant in 362 the UK

African SARS-CoV-2 variant 501Y.V2 into the UK

A structured open dataset of government interventions in response to 368

Tracking Public Health and Social Measures

Effectiveness 374 of isolation, testing, contact tracing and physical distancing on reducing 375 transmission of SARS-CoV-2 in different settings

Mitigation Strategies and Compliance in the Covid-19 378 Fight; How Much Compliance is Enough?

Feasibility of controlling COVID-19 outbreaks by isolation of cases and contacts

Effect of early application of social distancing interventions on 384

COVID-19 mortality over the first pandemic wave: An analysis of longitudinal 385

Markov Chain Monte Carlo 438 in Practice

COVID-19 Map -Johns Hopkins Coronavirus Resource Center

How many people has the coronavirus killed?

Excess Deaths Associated with COVID-19

Load US Census Boundary and Attribute

American Community Survey 2014-2018 5-Year Estimates

Estimates of the severity of coronavirus disease 2019: a model-based analysis

456 The impact of COVID-19 and strategies for mitigation and suppression in low-457 and middle-income countries

A new look at the statistical model identification

First 461 Case of 2019 Novel Coronavirus in the United States

The reproduction number of COVID-19 and its 464 correlation with public health interventions

A SEIR-like model with a time-dependent contagion factor describes 466 the dynamics of the Covid-19 pandemic. medRxiv

Iyaniwura 468 SA, et al. Quantifying the impact of COVID-19 control measures using a 469

Bayesian model of physical distancing

distancing to avert unmanageable COVID-19 hospital surges

The pandemic's true death toll: millions more than official counts. 474 Nature

Tracking excess mortality across countries during the 476

COVID-19 pandemic with the World Mortality Dataset. eLife

Excess mortality during the COVID-19 pandemic

Aden governorate

Supporting information 532 S1 Fig. Simulate and Recover. Sample fits to synthetic data generated using known R(t) profiles. The synthetic daily death (red circles and left y-axis) is generated using known synthetic R(t) profiles (red line and right y-axis). The median result of the fit is shown in blue and the shaded light-blue is the 95% CI. The recovered, median, R(t) is also shown in blue along with the 95% CI in light-grey.