The Dynamics of Infectious Diseases Associated With Injection Drug Use in Lawrence and Lowell, Massachusetts.

BACKGROUND: There are a wide variety of infectious complications of injection drug use. Understanding the trajectory of these complications might inform the development of an early warning system for human immunodeficiency virus (HIV) outbreaks that occur regularly among people who inject drugs (PWID).
METHODS: A distributed lag Poisson regression model in the Bayesian setting was used to examine temporal patterns in the incidence of injection-associated infectious diseases and their association with HIV cases in Lawrence and Lowell, Massachusetts between 2005 and 2018.
RESULTS: Current-month HIV counts are associated with fatal overdoses approximately 8 months prior, cases of infective endocarditis 10 months prior, and cases of skin and soft tissue infections and incision and drainage procedures associated with these infections 12 months prior.
CONCLUSIONS: Collecting data on these other complications associated with injection drug use by public health departments may be important to consider because these complications may serve as input to a sentinel system to trigger early intervention and avert potential outbreaks of HIV.

Outbreaks of human immunodeficiency virus (HIV) infection among people who inject drugs (PWID) have been a regular occurrence in the United States since the Scott County, Indiana epidemic in 2014–2015, with outbreaks happening in the following states: Ohio and West Virginia in 2017; Miami, Florida in 2018; and Seattle, Washington in 2018 [1-4]. From 2015 through 2018, the communities of Lawrence and Lowell, Massachusetts (MA) experienced an outbreak of HIV infection among PWID [5]. By June 30, 2018, 129 new infections met the case definition of outbreak-associated HIV infection used by the MA Department of Public Health (MDPH) and the US Centers for Disease Control and Prevention (CDC) in their joint investigation; by June 30, 2020, a total of 180 new infections have been linked to the outbreak. Until 2017 and 2018, respectively, Lawrence and Lowell did not have local syringe exchange programs in operation [6].

Many investigators have noted the wide variety of infectious complications of injection drug use and their epidemiological and economic impact. Other investigators have suggested that an increase in at least one of these conditions may be an early warning for the potential for outbreaks of HIV among PWID [7, 8]. These conditions can include the (1) viral hepatitides, (2) skin and soft tissue infections such as cellulitis and abscesses, and (3) bacteremia resulting from injection site contamination, which can lead to sepsis, endocarditis, and osteomyelitis from hematogenous spread [8, 9]. However, no studies to date have sought to look for a timeline; that is, a sequential trajectory of these conditions among populations of PWID. Before the 2015–2018 HIV outbreak, other infectious diseases were on the rise in the state, particularly among young PWID, with hepatitis C virus (HCV) infection almost doubling among young people 15–24 years of age from 65 to 113 cases per 100 000 people, between 2002 and 2009 [10]. By 2015, more than 2000 new cases of HCV infection were being identified annually in the under-30 age group in MA [11]. The MDPH effectively used reported HCV infection in people under 30 years old as a surrogate for injection drug use to prioritize follow-up of cases and guide locating public health interventions directed toward PWID. If there is a temporal pattern in the occurrence of infectious diseases among PWID, this exploratory analysis of their temporal occurrence may serve as a first step in the development of a sentinel, early warning system that could prevent HIV outbreaks among this population by directing the deployment of prevention resources to emergent hotspots. To better understand the trajectory of infectious diseases among PWID in the state, we examined the monthly hospital discharge data for 6 hospitals in the Lawrence and Lowell catchment area for the period 2005–2018. We specifically examined injection-related infectious disease diagnoses and associated medical procedures. We also examined (1) HCV and HIV infections newly reported to the MDPH Bureau of Infectious Disease and Laboratory Sciences and (2) fatal overdoses monitored by the MDPH Office of Population Health, for the period 2005–2018.

In this study, we seek to identify the temporal relationship among these other infectious diseases, overdose, and HIV diagnoses. In particular, we investigate whether case counts of any of these comorbid conditions in prior time periods have a lagged association with current-month HIV counts in examining the data from 2005 to 2018. We are not making short-term predictions about the Lawrence and Lowell outbreak in this study, but we seek to identify a pattern of relationships among these health events in these communities. However, if there are associations among these health events in Lawrence and Lowell, they can be explored in other geographical areas to determine whether they could indeed act as sentinel infections for increased HIV transmission among PWID.

METHODS

Data Sources

We obtained monthly aggregate counts of International Classification of Diseases, Ninth Revision and Tenth Revision Editions (ICD-9 and ICD-10), diagnostic codes for infective endocarditis as well as skin, soft tissue, and musculoskeletal infections among patients 15–50 years of age from hospital discharge data from 2005 to 2018 across 6 hospitals that serve the areas affected by the 2015–2018 HIV outbreak: (1) Anna Jaques, (2) Lawrence General, (3) Lowell General, (4) Lowell General, Saints Campus, (5) Holy Family, Haverhill, and (6) Holy Family, Methuen (Supplementary Table 1). These hospitals, with almost 50 miles between the 2 most distant from each other (Lowell General in the west and Anna Jaques in the east), are the primary hospitals for this entire region of the state and were chosen to capture the vast majority of the local health events occurring from 2005 to 2018 and considered in this study. As a robustness check, we also included a sensitivity analysis in which we examined data from diagnoses occurring by zip code of residence—restricted to 9 zip codes in Lawrence and Lowell only—rather than using the 6 hospitals as the sole setting for our study. This allowed us to capture people who may have resided in Lawrence and Lowell during our study period, but who traveled or were referred elsewhere in the state for medical care.

The source of the patient data by hospital and zip code is the Massachusetts Center for Health Information and Analysis (CHIA) (https://www.chiamass.gov/case-mix-data). Data on inpatient discharges are among the other health encounter data that CHIA receives and makes available for researchers. For each of these patient discharges, hospitals submit detailed information, including the following: patient demographics, admission and discharge information, diagnostic and procedural coding, provider details, and detailed charge information. The CHIA conducts routine quality assurance to assess validity and completeness and communicates with hospitals when supplemental data are needed to improve quality. This specific set of diagnoses was included because of their association with injection drug use in numerous studies in the United States and abroad [7, 9, 11–15]. More importantly, ICD-9 and ICD-10 codes corresponding to these diagnoses have been used previously to identify injection-associated infections in hospital data in the same manner in which we are doing here [7]. We restricted all of the data collected on infectious disease outcomes and associated clinical procedures for the 15–50 age group to reduce inclusion of many cases of these conditions that might occur outside the context of injection drug use (eg, endocarditis associated with intracardiac devices) and to mirror the age range of the cohort of PWID represented in the Lawrence and Lowell outbreak [16-19]. As additional proxies for these medical conditions or the level of suspicion for them among clinicians, we obtained monthly aggregate counts of (1) echocardiograms and (2) abscess incision and drainage and joint washout procedures from 2005 to 2018 at these same hospitals. To ensure that all diagnoses and procedures were consistent across ICD versions, ICD-9 and ICD-10 codes were converted using the Centers for Medicare and Medicaid Services’ General Equivalence Mappings (GEM) and confirmed by clinicians [20]. We collected monthly counts of HIV diagnoses and HCV cases reported to the MDPH from 2005 to 2018 for Lawrence and Lowell. These case counts were restricted to those with reported residence in the 2 cities. For the HIV cases, this means that these only represent a subset of the infections linked to the 2015–2018 HIV outbreak among PWID, some of which were reported with residence outside of these municipal jurisdictions. Finally, we collected monthly aggregate counts of fatal overdoses in Lawrence and Lowell utilizing the location of death for 2005–2018 from the Office of Population Health at the MDPH. Misclassification is always a risk. Here, we used the geographical areas of Lowell and Lawrence from the reported residence of the case (HIV and HCV), and these individuals may have been diagnosed at a different treatment facility. Likewise, although some fatalities may have occurred among individuals from other cities and towns, Lawrence and Lowell have longstanding local opioid-related problems, with opioid overdose deaths increasing rapidly in the 2 cities [5, 6].

Patient Consent Statement

Data collection occurred after ethical review of the study by the Massachusetts Department of Health and Yale University Institutional Review Boards (IRB) and the Partners Human Research Committee, all of which declared the study exempt from IRB review. Because the study used deidentified administrative data in aggregate, it did not include factors necessitating patient consent.

Statistical Analysis

To quantify the associations between the considered infectious diseases and other health outcomes (ie, overdose) preceding HIV case counts in Lawrence and Lowell, MA, we used a distributed lag Poisson regression model in the Bayesian setting. We modeled the expected number of HIV cases in a given month as a function of the counts of one of these other variables (see Figure 1 for the complete set) at several different lagged time periods (current-month, 1 month prior, 2 months prior, etc) simultaneously to determine whether and in which lagged period the variable is associated with HIV infection. Each component of the model was selected based on features of the data (eg, likelihood choice, correlation structure) and our primary study objectives (eg, inclusion of lagged predictors).

Figure 1.

Time series plots for human immunodeficiency virus (HIV) cases and each infectious disease and other health outcome variables. The first month corresponds to January 2005 and the last month corresponds to December 2018. Red dashed lines indicate the period of HIV outbreak in the area. Echo, echocardiogram; HCV, hepatitis C virus; IE, infective endocarditis; I&D, incision and drain procedures; OD, overdose; SSTI, skin, soft tissue, and musculoskeletal infections.

For each month of HIV data, we created lagged predictors of case counts for the selected variable and used them as predictors in the Poisson regression model. We specified a random walk model for the lagged regression parameters due to the large number of lags included in the model and the potentially high correlation between the lagged counts over time. The use of the random walk model led to more stable results during model fitting in this setting (ie, regularization), with similar smoothing methods used for distributed lag parameter estimation in previous work [21]. The model also included time series random effects to account for unexplained variability and correlation in the monthly HIV case counts across time. Prior distributions for the introduced model parameters are selected to be weakly informative to allow the data to drive the inference rather than our prior beliefs, and they represent standard choices made in past analyses using distributed lag hierarchical Bayesian regression models [22, 23].

Using output from the fitted model, we investigated the timing of the association by analyzing individual lagged regression parameters. In a sensitivity analysis, we also included results from a frequentist multivariable negative binomial regression model (same form as the previously described Poisson regression model without any random effects) in the Supplementary Web Appendix for each predictor. To assess model adequacy, we compute the Bayesian P value for each fitted model using an overall goodness-of-fit measure. The Bayesian P value compares properties of the observed data with data predicted from the fitted model, where values near zero or 1 suggest that the proposed model fails to adequately describe the observed data in some way [24]. See Supplementary Web Appendix for full details on the statistical model and Bayesian P value.

RESULTS

Figure 1 shows the time series plots for HIV cases and each infectious disease and other health outcome variables. The first month plotted in the figure corresponds to January 2005 and the last month corresponds to December 2018. Each plot highlights the magnitude and variability in each variable over time, but these plots are only descriptive in this setting. We used the distributed lag model to formally quantify the associations between these variables and HIV infection. We have bracketed the years associated with the HIV outbreak (2015–2018) with red lines for informational purposes only in the time series of HIV cases in Figure 1.

For the distributed lag modeling, we included monthly variables up to 12 months in the past to explain trends in current-month HIV infections. Figure 2 shows the estimates (posterior means) and 90% credible intervals (0.05 and 0.95 posterior quantiles) for the individual distributed lag regression parameters, presented on the relative risk scale. The interpretation of the relative risk in Figure 2 for one of the lagged parameters is based on a 1-case increase in the specific variable during the selected lag period. When these individual estimated effects are viewed together, the plot can inform about general trends in associations that vary across time.

Figure 2.

Results for the individual distributed lag regression parameters on the relative risk scale. Red indicates that the 90% lower bound of the credible interval is >1, whereas blue indicates that the upper limit is <1. Echo, echocardiogram; HCV, hepatitis C virus; IE,: infective endocarditis; I&D, incision and drain procedures; OD, overdose; SSTI, skin, soft-tissue, and musculoskeletal infections.

Based on Figure 2, current-month HIV counts are associated with overdoses cases approximately 8 months prior. An additional overdose occurring 8 months in the past is associated with a 3% increase (estimate, 1.03; 90% credible interval, 1.00–1.06) in the expected number of current-month HIV cases. In addition, current-month HIV counts are also associated with cases of infective endocarditis 10 months prior (estimate, 1.02; 90% credible interval, 1.00–1.05). Finally, current-month HIV counts are associated with skin and soft tissue infections (estimate, 1.00; 90% credible interval, 1.00–1.01), and incision and drain procedures associated with these infections (estimate, 1.02; 90% credible interval, 1.00–1.03) 12 months prior. Hepatitis C virus and echocardiograms show no individual lag associations with current-month HIV case counts.

Supplementary Web Figure 1 displays model comparison sensitivity analysis results where the first column for an outcome represents results from the hierarchical Bayesian distributed lag Poisson regression model and the second column is from the multivariable negative binomial regression model. Although the results from both models are similar, the distributed lag model yields estimates generally pulled towards the null with reduced uncertainty, allowing for a clearer view of association trends across time. The Bayesian P value results suggest that the proposed model is adequate for each of the conditions, with values ranging from 0.43 to 0.48.

Supplemental Web Figure 2 shows the results for the additional sensitivity analysis, in which infections and procedures were derived by zip code of residence of patients rather than the hospital from which they received their care. Most of the lag associations are similar to the previously described relationships, although the point estimates and credible intervals differ. In the zip code-based analysis, current-month HIV counts are associated incision and drain procedures with an additional lag 3 months prior. Current-month HIV counts in the zip code-based analysis are also associated with echocardiograms 8 months with a potentially spurious inverse correlation between current-month echocardiograms and current-month HIV counts.

DISCUSSION

The present study shows an association between overdose, infective endocarditis, skin and soft tissue infections, and incision and drain procedures with HIV cases in Lawrence and Lowell, MA from 2005 to 2018. We want to clarify that we are not suggesting that these health events were predictive of the HIV outbreak among PWID in these cities in 2015–2018. What we are suggesting is that over the entire 13-year period from 2005 to 2018, these health events maintained a lagged association with current-month HIV counts, trailing them in these relationships out to 1 year previously. In this study, we are extending the work of the critical case investigations by Cranston et al [5] and Alpren et al [6] further back in time to assess the dynamics of infectious diseases and other outcomes associated with injecting drug use in Lawrence and Lowell starting in 2005. We were particularly interested in examining the lagged patterns between current-month HIV case counts and all the included conditions and procedures (see Figure 1 for the complete set). Three major patterns emerged. Human immunodeficiency virus case counts in a given month are associated with overdoses 8 months prior, infective endocarditis 10 months prior, and skin and soft tissue infections and accompanying incision and drain procedures, occurring some 12 months earlier. The concurrent, overlapping timing of some of these conditions and procedures, in addition to the clinical plausibility of their clustering in this way, lends additional credibility to our findings.

Injection-related skin and soft tissue infections are the most common reason for hospitalization of PWID [25]. There are different correlates of these infections among PWID, but several are relevant in this context, including increased injection frequency, increased sharing of syringes and shorter injection histories, and mixing of drugs (eg, heroin and cocaine) [26]. The first 2 factors (ie, increased injection frequency and sharing) would be important for increased HIV transmission risk. The last 2 (ie, shorter injection histories, and mixing of drugs) would mirror the epidemiology of injection drug use in MA and the age range of the cohort in Lawrence and Lowell, where mixing of drugs is common and the age range of those infected in the outbreak were largely under 40 years of age [5, 27].

The association of cases of infective endocarditis 10 months prior with current-month HIV counts is consistent with reports of the rise of this condition among people who use drugs in Massachusetts [28]. This lagged association suggests that this condition may be an important marker of contemporary HIV risk as other studies have suggested as well [8]. The fact that skin and soft tissue infections and associated incision and drain procedures appear in lagged association 12 months before current-month HIV counts, but before infective endocarditis cases, is consistent with the role of repeated subcutaneous or intramuscular injection in allowing Staphylococcus aureus, streptococci, and other Gram-negative bacilli access to the bloodstream and subsequent infection of the heart valves and endocardium [29-31].

Overdose is a direct consequence of substance use, and its presence in lagged association with current-month HIV counts makes a direct link between substance use and HIV infection in our study. We were surprised to find that HCV infection does not show an association with current-month HIV case counts in our study, which examines lags out to 1 year. This may be due to even earlier penetration of HCV into drug-using communities in Lawrence and Lowell, which is consistent with the rapid and broad transmission of this virus among PWID even after a few years of injecting [32, 33]. In fact, HCV has been so extensively distributed among PWID for over a decade in the state, particularly in Lawrence and Lowell, that its ubiquity may reduce its utility as a proximal marker of current-month HIV counts, even as it is a robust indicator of injecting drug use [34, 35]. Although we could have extended our study to events further back in time, our choice to restrict this analysis to lagged associations within a year of current HIV counts was based on 2 considerations: first, to reduce the chance for type 1 error, false-positive signals from associations spuriously seen in the more distant past; and second, to strengthen the inferences made in our study by enlarging the sample size of health events under consideration through opening up additional months of cases for examination by keeping the lags to 1 year.

Our additional sensitivity analysis, in which lagged associations were investigated using hospital discharge data based on patient zip codes indicating residence in Lawrence and Lowell, Massachusetts (and not restricted to the 6 local hospitals), further support our primary analysis. The same 3 major patterns emerge in the zip code-based analysis, with overdose, infective endocarditis, skin and soft tissue infections, and accompanying incision and drain procedures matching the timing of associations with current-month HIV counts in the hospital-based one. In the zip code analysis, an additional lag emerges with incision and drain procedures 3 months prior associated with current-month HIV counts. Given the common nature of skin and soft tissue infections among PWID, it is not surprising to see this additional lagged association appear.

The notable difference between the zip code- and hospital-based analyses is with the emergence of lagged associations of echocardiograms with current-month HIV counts. The first lag, in which echocardiograms 8 months prior is consistent with the lagged association with current-month HIV counts seen with infective endocarditis 10 months prior. In fact, if patients from Lawrence and Lowell were referred to facilities that could offer the specialty care needed for these complex infections, which can involve prolonged courses of antibiotics and valve replacement surgery, these additional diagnostic procedures could be expected in our analysis [36]. In fact, at least one of the major hospitals caring for patients with infective endocarditis in the region conducts 30 000 echocardiograms a year, and many are for this condition [37]. The association between echocardiograms in the current month and current-month HIV counts is harder to explain and may represent a type 1 error. In future work on the dynamics of infectious diseases among PWID in Massachusetts and other states, we plan to examine the relationship among hospitals, zip codes, and these diagnoses in more detail.

Although previous research has described the presence of the same infectious diseases among PWID that we have examined here, no studies to date have described the temporal patterns of emergence of these infections in a community of PWID [7, 13]. Furthermore, although other investigators have suggested that conditions such as infective endocarditis and HCV infection should alert public health officials to the presence of a growing epidemic of injection drug use in a community and the risk for other blood-borne infections such as HIV, none have set out to quantify the trajectory of this risk [8, 38]. Although in 2016 the CDC identified 220 counties that were especially vulnerable to HIV and HCV outbreaks associated with injection drug use, the cities of Lawrence and Lowell, MA were not in the predominantly rural counties identified by CDC [39]. The CDC was trying to model the prevalence of PWID at the county level, and the variables included were drug overdose deaths, prescription opioid sales, per capita income, white, non-Hispanic race/ethnicity, unemployment, and buprenorphine prescribing potential by waiver [39]. Our research suggests that new models of HIV risk among communities of PWID should include specific other infections associated with injection drug use. The data used for this study were provided by the MDPH, and similar data are likely to be available to other state departments of health for these kinds of analyses. In fact, many states have their own All-Payer Claims Databases (APCDs) and the State Inpatient Databases (SIDs), which are part of the family of databases and software tools developed for the Healthcare Cost and Utilization Project (HCUP) under the aegis of the Agency for Healthcare Research and Quality (AHRQ) [40, 41]. Together, the State-specific SIDs encompass more than 97% of all US hospital discharges from 49 states and include principal and secondary diagnoses and procedures. Thus, both APCDs and state-specific SIDs represent existing sources for data on these other conditions with association with current-month HIV counts in our study.

Our study presents lagged associations of several infections, conditions, and therapeutic and diagnostic procedures with a given month’s HIV case counts. The major limitation of our study is that it is an examination of a single geographical area over time. We also realize that many of these infections and patterns may be representative of health risks common to PWID even in the absence of outbreaks of HIV infection. We also cannot attribute causality in this study to these patterns of infectious disease among PWID. To generalize these associations and make causal claims for their impact on HIV outbreaks, and to potentially use these conditions as ways to intervene earlier to prevent them, we hope to demonstrate that this lagged structure of infections (or a consistent association with others), with similar, although not necessarily identical, temporal relationships, is recapitulated in other settings with HIV outbreaks among PWID and, conversely, is not present in places with high prevalence of injection drug use, but without increases in cases of HIV among PWID. Because the epidemiology of substance use and its complications differ by place and time, we are currently extending this analysis to jurisdictions that have experienced similar outbreaks of HIV infection among PWID in Ohio, Florida, West Virginia, Kentucky, Washington, as well as in places in these states in which injecting drug use is common but no outbreaks of HIV have occurred to date. In addition, these new analyses will provide the opportunity to examine the impact of different approaches to the management of the complications of injection drug use in these settings and their influence on the patterns of health events we have described here. For instance, as we have noted Lawrence and Lowell, Massachusetts did not have a syringe exchange program before the 2015–2018 outbreak, whereas cities such as Seattle have had ongoing syringe exchange services, which could influence the patterns and dynamics of infectious disease transmission in these places.

Finally, the consistency of the results of the distributed lag model and the multivariable negative binomial regression model presented in the Supplementary Web Appendix gives additional strength to our conclusions. The choice of the distributed lag model with random walk distribution for the regression parameters for our primary analysis, which shrinks estimates towards the null while reducing uncertainty, likely results in fewer spurious associations being identified than a completely unstructured model of the same form [42]. In total, we observed 4 statistically significant associations using our primary analysis method versus 7 with the comparable multivariable negative binomial regression model.

CONCLUSIONS

The opioid epidemic in the United States, overdoses, and the infectious complications of injection drug use are an American health crisis [43]. Although expanding access to community-based harm reduction programs and services and improving the care continuum for individuals with infectious diseases and substance use disorders are critical, a better understanding of the epidemiology of all of these related health events among PWID may provide a way to further mitigate morbidity and mortality among this underserved population. Furthermore, if, upon additional research, these patterns of complications of injection drug use can serve as an early warning system for communities and decision makers alike, this may further strengthen the local rationale for the expansion of programs to serve PWID in the United States.

Supplementary Data

Supplementary materials are available at Open Forum Infectious Diseases online. Consisting of data provided by the authors to benefit the reader, the posted materials are not copyedited and are the sole responsibility of the authors, so questions or comments should be addressed to the corresponding author.

Click here for additional data file.