Introduction

The USA is in the midst of an opioid crisis that first emerged in the 1990s but has reached epidemic proportions within the last decade. Fueled by the advent of synthetic opioids, the number of annual overdose deaths almost doubled between 2010 (38,329 deaths) and 2019 (70,360 deaths) in the USA (Center for Disease Control and Prevention: CDC). While death totals remained relatively steady between 2017 and 2018, recent estimates place the 2020 national overdose death total at 93,331—an alarming 32% increase from 2019 [1]. Research indicates that the onset of the pandemic was associated with an increase in overdose deaths nationally [2, 3], regionally [4,5,6,7], and locally [8,9,10,11]. The increase in opioid overdose deaths (OODs) during the pandemic has been attributed to various factors: reduced access to interventions, increased levels of stress due to isolation and loss of mental health support, and changes in the types/combinations and purity of drugs and patterns of drug use. The impacts have severely affected racial/ethnic minority communities, lowering life expectancies in Black and Hispanic Americans. The National Center for Health Statistics notes that the drop in life expectancy among Hispanic Americans is 3 years; Black Americans is 2.9 years, while that of the White population is 1.2 years. Simultaneously, the rise in OODs among Black Americans (2013–2018) is outpacing that of White Americans [12].

To understand the effects of the pandemic on the opioid crisis, studies have primarily examined changes in OODs over time. However, the locational or spatial variations in OODs have not been addressed. Our study examines the spatiotemporal changes in OODs in Milwaukee County, Wisconsin, a racially diverse segregated metropolitan region with dramatic racialized health disparities [13]. We obtained our OOD dataset from the Milwaukee County Medical Examiner’s Office. This dataset is comprised of those deaths where an opioid was identified in the attached toxicology reports. These include overdose deaths where the solely identified drug was listed as an opioid in the toxicology report (293 cases) and “mixed-drug” overdoses where an opioid was identified along with other drugs in the report (1188 cases).

To assess the impacts of the pandemic, we examined OOD variations over time (January 1, 2017, to December 31, 2020) and space (census tracts within Milwaukee County). First, we conducted interrupted time series analysis, which uses autocorrelational analysis to address seasonality. Next, we examined the spatial changes over time in OODs through spatial modeling in geographic information systems (GIS), as it offers a powerful analytical approach to the study of the opioid crisis [14,15,16]. Our georeferenced OOD data, used in combination with Census data in GIS, enabled us to examine the demographic shifts in OODs over time and place.

Our study offers rich empirical insights and innovative methodological approach to examining the opioid crisis during the COVID-19 pandemic. It is attentive to place-based differences and racialized health disparities and offers unique empirical insights for effective policy responses. The study extends our past research on Milwaukee’s pre-pandemic opioid crisis [17], which highlighted both the gravity of the crisis among low-income Black and Hispanic neighborhoods and the ineffectiveness of harm reduction policy responses within such neighborhoods. This paper offers new empirical and methodological contributions to assist in the formulation of data-guided policy interventions, so that lives in vulnerable communities can be saved.

Study Site, Data, Methodology

Milwaukee County, with an estimated population of 945,726, is the most populous county in Wisconsin. It is racially diverse, with 58.9% White, 26.3% Black, and 15.6% Hispanic or Latino origin (U.S. Census Bureau, 2020). The county is also racially and economically highly segregated, leading to dramatic racialized health disparities. The hyper-segregated City of Milwaukee (population 590,157) exhibits concentrated urban poverty in neighborhoods occupied primarily by Black and Hispanic communities [18]. Fatal overdose rates are related to such racial and economic segregation [17], making Milwaukee an appropriate case study for examining the effects of the pandemic on OODs.

To examine the impact of COVID-19 on opioid overdoses in Milwaukee County, we collected data on opioid-related overdose deaths between January 1, 2017, and December 31, 2020 (we exclude overdose deaths directly attributed to non-opioids). Wisconsin issued its state-wide Stay at Home Order on March 23, 2020. Thus, for the current study, we established this as the “intervention” date and defined the preceding time period from January 1, 2017, to March 23, 2020, as the pre-intervention/pre-pandemic stage and the time period between March 24, 2020, and December 31, 2020, as the post-intervention/post-pandemic stage.

We obtained OOD data from detailed death certificates provided by the Milwaukee County Medical Examiners’ Office. The death certificates included toxicology reports listing the drug(s) involved in the overdose, allowing us to identify and extract only opioid overdose deaths data. This data contained addresses noting the location of each overdose death. Through address-based geocoding, we joined this dataset to the administrative boundary shapefile of census tracts collected from the TIGER/Line database (www.census.gov) using ArcGIS Desktop 10.7. We also linked the locations of overdoses with census tract level demographic data obtained from the U.S. Census Bureau (2010 Census data from data.census.gov platform). Other relevant data sets were obtained from the US Department of Health & Human Services and from Policy Map (incarceration rates, health insurance), a repository of policy-oriented data sets (Policymap.com). We conduct statistical analysis and spatial modeling. Two significant methods are applied: interrupted time series analysis and spatial empirical Bayesian smoothing.

Interrupted Time Series

Epidemiologists have long used interrupted time series analysis to assess the effectiveness of interventions aimed at reducing OODs [19,20,21,22,23]. This technique has also proven useful for examining the impact of the COVID-19 pandemic and related responses on OODs [8, 24].

Social phenomena exhibit autocorrelation and seasonal effects that can cause a temporary change in the rate of events; interrupted time series analysis captures the effects of autocorrelation and seasonality to enable us to assess the impact of a particular policy. This is done by analyzing the change in the level and slope of the time series after an intervention has been applied and comparing to the structure of the temporal dynamic before the intervention. More recently, scientists have developed a framework for evaluating the causal influence of interventions relying on Bayesian statistics. Here we investigate the effect of COVID-19 containment measures, specifically the “stay-at-home” order, on OOD trends in Milwaukee using Bayesian structural time series (BSTS) models [25]. A diffusion-regression state-space approach was used to forecast counterfactual patterns in trends that would have existed had there been no intervention in a hypothetical scenario. This method helps us to measure the effect and statistical importance of a specific incident on our variable of interest, OOD rate. Bayesian structural time series (BSTS) models were implemented using the “bsts” [26] and “Causal Impact” [27] packages in R.

Spatial Empirical Bayesian Smoothing

To control for the effects of population size, we computed OOD rates by census tract in Milwaukee County; these rates exhibit intrinsic variance instability, suggesting that a rate’s precision as a measure of underlying risk is inversely proportional to the size of the at-risk population leading to significant standard error in rates calculated from small populations. In other words, an area with a smaller population will naturally have fewer OODs; this area could easily be misclassified as a low-risk area or vice versa, if the population is not controlled for [28]. For example, an area with 100 residents and 10 cases could be classified as high risk, whereas an area with 1,000,000 residents and 75,000 cases is low risk (if a 10% rate is used as a boundary). To reduce variation caused by the size of the population, spatial empirical Bayes smoothing techniques were applied, thereby improving the precision of the crude rate by borrowing strength from other observations. In spatial empirical Bayes, rather than taking the same overall reference rate for all, the reference rate is computed for a spatial window for each individual observation. The spatial empirical Bayes smoothed rate is calculated as a weighted average of the crude rate and the prior or reference rate estimated from the spatial window surrounding a given observation, consisting of the observation and its neighbors [29].

Results

To investigate the impact of the pandemic on OODs, we examined monthly OODs from January 1, 2017, to December 31, 2020. The stay-at-home order issue date, March 23, 2020, was used to divide the study into pre-intervention (pre-pandemic) and post-intervention (pandemic) periods. The average monthly OODs increased by approximately 12. The peak monthly overdose death total during the pre-intervention stage was 37. However, in the post-intervention period, Milwaukee County experienced a peak of 57 OODs per month and a minimum of 23 monthly OODs. These differences suggest that factors related to the pandemic increased OODs in Milwaukee.

In Table 1, the numbers and percentages of OODs prior to and after the March 23, 2020 “stay-at-home” directive, disaggregated by mode of death, toxicology/cause of death, gender, race, geographic discordances (when overdose incident location differs from the residential address), and age, are presented. In line with OOD patterns in other U.S. cities [8], we found that OODs increased significantly in the post-intervention stage in Milwaukee County. Simultaneously, we find that the principal factors (race, gender, age, mode of death) shaping OODs remained consistent in the pre- and post-intervention stages. For example, the primary mode of death in both periods is accidental, but the percentage shifted slightly from the pre (96%) to the post (95%) stage. Males are the primary victims in both periods, but the female OODs slightly declined from 32 to 29%. In terms of race/ethnicity, the OOD rate was highest among White community members (66%) followed by Black (21%) and Hispanic (9%) community members during the pre-intervention stage. These numbers slightly increased for Black community members (25%) in the post-intervention period but slightly declined among the White (63%) and Hispanic community members (8%). Age groups between 30–40 (26%) and 50–60 (23%) were most affected by OODs in the pre-intervention period, followed by the 40–50 (21%) age group. These rates remained steady in the post-intervention stage, but OODs in the 0–30 age group showed a slight decline from 20 to 13%. We also calculated geographic discordance to see if the county’s OODs primarily affected Milwaukee County residents or outsiders who travel into the county. In both stages, OODs primarily affected county residents.

Table 1 Opioid overdose deaths in Milwaukee County disaggregated by mode, drug evident in the toxicology report, cause, gender, race, geographic discordance, and age
Full size table

In terms of drug involvement in overdose cases, fentanyl-related OODs increased by 25 and heroin-related OODs decreased by 23. As depicted in Fig. 1, fentanyl-related deaths comprised a greater proportion of overall monthly OODs, and heroin-related deaths declined.

Fig. 1
figure 1

Monthly overall, fentanyl-related, heroin-related, and other opioid-related OODs

Full size image

Next, we tested for dependence among monthly opioid overdose death counts using the Box-Ljung test. Using a standard value of the risk of α (confidence interval) = 0.05 (95% confidence interval), the Box-Ljung test was conducted to investigate whether autocorrelations of OOD time series are different from zero (null hypothesis: “the monthly OODs are independently distributed”; alternative hypothesis: “the monthly OODs exhibit serial correlation). Our results (X-squared = 13.614, df = 1, p-value = 0.0002245) indicate that monthly opioid OODs were dependent on each other. Because of this dependency, we next applied the BSTS model and forecasted a monthly OOD rate for the post-intervention period if pandemic had not taken place and compared the model-estimated mortality rate trend with the actual observed OOD rate in the post-intervention period.

As seen in Fig. 2, during the post-intervention period, the average number of monthly OODs was 38.00. In the absence of the pandemic and “lockdown,” we would have expected average monthly OODs of 26.29 (95% confidence interval of 22.13, 30.28). Subtracting this prediction from the observed monthly OOD yields an estimate of the causal effect of various factors arising from the pandemic had on the monthly OOD. This effect is 11.71 with a 95% interval of 7.72, 15.87. Summing up the individual data points during the post-intervention period (March 23, 2020, to December 31, 2020), there were 342 OODs. Had the pandemic not occurred, the predicted number would be 236.59 with a 95% confidence interval of 199.21 to 272.52. Thus, the monthly opioid overdose death showed an increase of + 45% (the 95% confidence interval is + 29% and + 60%). These results show that the positive trend observed after March 23, 2020, is statistically significant. The probability of obtaining this effect by chance is very small (Bayesian one-sided tail-area probability p = 0.001).

Fig. 2
figure 2

Monthly opioid overdose deaths, vertical line indicates the stay-at-home order date (March 23, 2020). Blue shaded region indicates 95% confidence interval for BSTS forecast

Full size image

On May 13, 2020, the Wisconsin Supreme Court struck down “stay-at-home” orders as unconstitutional [30]. Between the stay order and the supreme court ruling, there was an increase in monthly OODs. However, this increase persisted well after the court ruling, indicating the high impact of the pandemic upon the increase in substance abuse. In fact, when we examine the time period between the stay-at-home order and the end of 2020, we find that for several months, monthly OODs remained above the 99% confidence bounds of the BSTS forecast (Fig. 2), highlighting the importance and impacts of unobserved confounding factors (e.g., pandemic induced socioeconomic impacts and socio-political changes).

Next, we examined spatial variation in OODs in Milwaukee County from 2017 to 2020. Milwaukee’s reputation as one of the most racially segregated metropolitan areas in the country is a defining factor influencing the region’s socioeconomic, political, and cultural environments [18]. Figure 3 displays a map of segregated census tracts in Milwaukee County.

Fig. 3
figure 3

Racial segregation map displaying majority racial group in each census tracts

Full size image

To obtain OOD estimate rates and map the underlying risk for OODs within each census tract, we used a spatial empirical Bayes smoothing technique. The parameters of the prior distribution were obtained empirically from the Milwaukee County Medical Examiner Office dataset. The posterior distribution was obtained by combining the Poisson likelihood with the prior distribution. This resulted in a distribution of possible values of the underlying risks of the census tracts that were conditional upon the observed number of OODs. The empirical Bayes estimate of the true risk was taken as the mean of the posterior distribution.

Figure 4 displays the distribution of the underlying OOD risks at the census tract level in Milwaukee County from 2017 to 2020. Of the 296 census tracts examined, a persistently high OOD rate was observed in 6% of them; these are urban census tracts with predominantly young (31 years), Black (39%), and Hispanic (24%) populations, with low educational attainment (10.63% hold a bachelor’s degree or higher) and low median household income ($36,351). The incarceration rate on average is higher (5%) and 10% of residents do not have health insurance (Policymap.com). Only 69% of households in these census tracts have internet subscriptions, limiting their access to information. About two-thirds (63%) of the census tracts with persistently high OOD rates are designated as medically underserved areas (US Department of Health & Human Services, 2021), as characterized by high infant mortality, high poverty, and few primary care providers. A comparison with Fig. 3 shows that the OOD distribution aligns with racial segregation patterns. About half (53%) of these communities have experienced concentrated persistent poverty for decades; the health and well-being of residents have been persistently low.

Fig. 4
figure 4

Spatial distribution of the underlying OOD risks at the census tracts level from 2017 to 2020

Full size image

In contrast, 3% of census tracts exhibited persistently low OOD rates. These census tracts are suburban and more affluent, occupied by a predominantly White (83%) and middle-aged (41 median age) population with a higher median household income ($87,079), higher educational attainment (34% hold bachelor’s degree or higher), and relatively high access to health-related resources (99% of the residents have access to health insurance and health care). Other measures of social well-being are evident: 85% of households have an internet subscription and the average incarceration rate is low (3%).

Overall, pre-pandemic (2017–2019) OOD patterns in Milwaukee County indicate that the opioid crisis has disproportionately impacted historically marginalized Black and Hispanic communities living in central Milwaukee neighborhoods (Fig. 4). This corresponds with findings from our previous study [17] that the policies, resources, and interventions that are being designed and implemented in Milwaukee County are primarily benefiting White communities, not communities of color. The low impact of policy interventions in Black and Hispanic communities has further widened their pre-existing health inequalities.

The 2020 map in Fig. 4 indicates that the pandemic’s effect on OODs has been highest among Milwaukee’s Black and Hispanic communities. But, the map also shows a spike in OODs in suburban tracts, which were previously minimally affected by OODs.

In line with prior studies [31, 32], our interrupted time series analysis showed that OODs have significantly increased during the pandemic. To examine the geographic impact of the pandemic on the overdose crisis in Milwaukee County, differences between the 2020 OOD rate and the average of OOD rate from 2017 to 2019 were calculated for each census tract in Milwaukee County and mapped. Figure 5 shows that 199 census tracts (2/3 of all tracts) in Milwaukee County experienced an increase in the number of OODs relative to pre-intervention levels (2017 to 2019). Geographically, these census tracts are distributed throughout the county, indicating that the pandemic impacted a variety of communities in Milwaukee County.

Fig. 5
figure 5

Increment (yes/no) and increment percentage in the numbers of opioid overdose deaths in Milwaukee County, March 23 to December 31, 2020

Full size image

For all census tracts, OOD increment percentages were calculated to examine the incremental changes in OODs after the pandemic in comparison to the 5-year average number of cases. Results show that OOD rates in two areas were particularly impacted: (1) predominantly poor, Black neighborhoods in the inner city and (2) predominantly affluent White suburban census tracts. It is notable that there are areas where the death rates remained steady, showing little percentage change over time (this includes tracts with previously high OOD rates). Importantly, due to relatively low pre-intervention OOD frequency, many suburban tracts displayed sizeable percent increases that do not reflect large increases in the total number of OODs.

The urban tracts that experienced a 10,280% increment in OODs suffer from racial and economic segregation and experience concentrated, persistent poverty. The population in these tracts is 72.46% Black (11.52% White, 16.02% other minorities) with a median age of 25.77, a 34% unemployment rate, and only 5% holding a bachelor’s degree or higher with a low median household income ($31,192), 69% internet subscription rate and relatively high (7%) incarceration rate. The effect of the pandemic upon the residents of the urban tracts has further widened pre-existing socioeconomic disparities.

In contrast, the suburban census tracts that experienced an 11,600% increment are affluent and educated. The population is 83.63% White (6.29% Black, 10.08% other minorities) with a median age of 33.53, 11% unemployed, 45% holding a bachelor’s degree or higher, and a median household income is $75,959. Other social well-being indicators include a 90% internet subscription rate and low (2%) incarceration rate. Despite high economic and social well-being, suburban tracts have nonetheless experienced pandemic-related duress, likely contributing to OODs.

Figure 6 shows the spatial patterns of increment changes in fentanyl-, heroin-, and other opioid-related OODs. In Milwaukee County, there was a disproportionate increase in OODs involving fentanyl among urban Black and Hispanic community members. However, deaths involving heroin declined significantly, perhaps due to supply chain disruptions [33]. These shifts have occurred nationwide; Black Americans are experiencing the highest rate of increase in fentanyl-related OODs [34,35,36].

Fig. 6
figure 6

Spatial patterns of fentanyl-, heroin-, and other opioid-related deaths, March 23 to December 31, 2020 (increments relative to pre-intervention numbers)

Full size image

Discussion

The opioid crisis in the USA has profoundly worsened during the pandemic. Over 2/3 of census tracts in Milwaukee County experienced significant increases in OODs. The maximum number of monthly overdose deaths increased from 37 prior to the March 23, 2020, state-wide “stay-at-home” order in Wisconsin to 57 afterwards. This increase was consistent across spatially defined communities, across demographic variables (sex, age, race) and across mode of death, indicating that OODs have not only intensified in areas where they were previously evident but have spread to communities that were previously immune to the opioid crisis. Fentanyl has emerged as the prevalent drug involved in OODs, while heroin-related OODs have declined. Overall, the involvement of other drugs (e.g., cocaine, ethanol, amphetamines, sedatives, other opioids) has remained steady, although some variation is evident at the neighborhood scale.

From a spatial perspective, the number of deaths is higher in inner city tracts, where the opioid crisis has been prevalent for many years. The worst-affected areas are primarily Black/Hispanic communities shaped by decades of racial and economic segregation, with concentrated poverty and health inequalities. These vulnerable neighborhoods, with little access to health care, have limited resources to address opioid abuse, leading to a rise in OODs. However, the pandemic has also led to a surge in OODs in more affluent and educated White suburbs with historically high access to health care. We speculate that many of the resources that these communities previously enjoyed were diminished during the pandemic, reducing support for those with opioid use disorder and increased risk for OODs. Some spatial variation in drug involvement in OODs is evident. Fentanyl is dominant in OODs in poor urban tracts, where there has been a decline in heroin-related OODs. In contrast, suburban OODs more commonly involve fentanyl, heroin, and a mix of drugs. Further community-based research is needed to investigate the impact of the pandemic on different communities in Milwaukee County so that appropriate policy interventions can be implemented.

Converging factors likely contributed to the increase in OODs during the pandemic. First, adherence to social distancing guidelines likely led to large changes in drug use patterns, including higher rates of drug use in isolation, that increased the risk for fatal overdose. Social distancing likely also limited access to and utilization of critical community resources and support networks. Second, the pandemic led to significant shifts in treatment practice that may have introduced new risks to patients [40,41,42], including a transition from in-person/in-patient interventions to telehealth and the offering of large amounts of take-home medications, such as methadone. Third, in response to the pandemic, public health departments [43, 44] were compelled to redirect their resources to address COVID-19 infections, thereby limiting their capacity to implement key harm reduction interventions, including syringe exchange and naloxone distribution programs. Fourth, the pandemic also disrupted drug supply chains internationally and nationally, creating disruptions in preferred drug availability, leading to the use of drugs from unfamiliar sources. The consequent changes in drug use patterns and the questionable quality and purity of drugs being sold on the street [33] likely increased the risk of overdoses. Finally, the pandemic likely exacerbated key determinants of opioid use, such as increase in stress, anxiety, and related mental health conditions, pre-existing health disparities, and socioeconomic uncertainties, thus causing an increase in drug use.

As in any study, there are limitations. First, our data collection did not include nonfatal OODs; studies have shown that the number of survived opioid overdoses also increased [10, 37] after the onset of the pandemic. Future research must examine the impact of the pandemic on nonfatal overdoses across various communities and identify any disparities. Second, interrupted time series analyses can be susceptible to history bias due to other interventions or events occurring around the period of analysis [38, 39]. Finally, there is a need to examine the differential impacts of the COVID-19 pandemic on OODs in suburban and rural communities.

Conclusion

The impact of the COVID-19 pandemic on OODs in Milwaukee has been tremendous and has crossed demographic, socioeconomic, and neighborhood boundaries in this diverse and highly segregated metropolitan area. A better understanding of contributing factors is needed to guide interventions at the local, regional, and national scales.