A multi-state model to estimate incidence of heroin use
© Sanchez-Niubo et al.; licensee BioMed Central Ltd. 2013
Received: 21 June 2012
Accepted: 11 January 2013
Published: 14 January 2013
Existing incidence estimates of heroin use are usually based on one information source. This study aims to incorporate more sources to estimate heroin use incidence trends in Spain between 1971 and 2005.
A multi-state model was constructed, whereby the initial state “heroin consumer” is followed by transition to either “admitted to first treatment” or to “left heroin use” (i.e. permanent cessation or death). Heroin use incidence and probabilities of entering first treatment ever were estimated following a back-calculation approach.
The highest heroin use incidence rates in Spain, around 1.5 per 1,000 inhabitants aged 10–44, occurred between 1985 and 1990; subdividing by route of administration reveals higher incidences of injection between 1980 and 1985 (a mean of 0.62 per 1.000) and a peak for non-injectors in 1990 (0.867 per 1,000).
A simple conceptual model for heroin users’ trajectories related to treatment admission, provided a broader view of the historical trend of heroin use incidence in Spain.
KeywordsBack-calculation Epidemiology Heroin Incidence Multi-state model
Recently, there has been increasing interest in ascertaining illegal drug use incidence for planning and evaluating prevention strategies . In the case of heroin use, as survey data is not effective [2, 3], incidence has been estimated from users eventually showing up mostly in treatment registers [1, 4–9]. This incidence can be referred to as “problematic use incidence” and its trend could provide a satisfactory overview of problematic heroin use, assuming a constant proportion over total incidence.
The present study is an attempt to take the definition of problematic use incidence further by incorporating a proportion of individuals that used heroin who never show up in the main source under study (usually treatment). This proportion is based on other kinds of information, aggregated sources, estimates, or assumptions, such as mortality and cessation rates.
The idea is to study the unobserved entry (or immigration) of people to the state of “consumer”, based on a later first entry to treatment. Since heroin users may exit the state of consuming heroin before entering treatment, whether due to death or permanent cessation of their consumption, a “left heroin use” state is added to the model. The situation is thus represented by a set of mutually related states, in a so-called multi-state model with immigration.
This approach has parallels in studies in the HIV field, where a multi-state model was presented describing progression of HIV disease from infection to AIDS in several stages [10, 11] or Rossi’s dynamic “mover-stayer” model as a theoretical approach applied to simulate a complete drug user “career” [12, 13].
We shall use a back-calculation type approach to estimating the incidence. This is similar in spirit to de Angelis et al. . However, following the presentation of Aalen et al. , we have found it useful to display the approach as a simple multi-state model.
To our knowledge, previous approaches had to assume that treatment availability was stable over time [1, 4–9]. As treatment data is usually only available for a limited period of time, a lag time distribution between heroin use onset and first treatment ever has been employed to avoid underestimating incidences. However, this assumption may be too strong if important changes in treatment availability had occurred. Therefore, the back-calculation approach used in this study avoids this assumption.
In Spain, heroin use generated important health problems in the eighties and early nineties , thereafter all indicators (mortality, treatment admissions, hospital emergencies, surveys, etc.) showed a decreasing trend until 2006 . Efforts to calculate incidence of heroin use have been and still are considerable to understand the overall trend and assess its consequences . To date only one previous study has estimated heroin use incidence in Spain. However, due to important changes of treatment availability, the reliability of the estimation was questionable .
The multi-state model approach may be successfully used to estimate heroin use incidence with the available Spanish heroin users’ treatment data. We use the same treatment database as the previous study to assess differences in the incidence estimates depending on the approach used.
The objective of the present study was to estimate heroin use incidence in Spain through a multi-state model with immigration and assess differences with previous study’s estimates.
Incidence estimation: the multi-state model with immigration
We shall now formulate our back-calculation type model. We let t be the calendar time and assume that the number of new heroin users per unit of time is Poisson distributed. The Poisson assumption is standard and well justified on statistical grounds . Moreover, we let h t denote the expected number of people entering state 1 (heroin use) at t. The number p t denotes the probability that a given heroin user initiates their first treatment ever at time t, given that they were in state 1 at the previous time. We assume this probability independent of the era when their heroin use began. The number q t denotes the probability of an individual leaving heroin use at time t, given that they were in state 1 at the previous time. The cause could be death or other permanent cessation of heroin use. We want to estimate the parameters h t , the expected number of new heroin users at time t. The probability q t is assumed known.
Maximizing this likelihood yields estimates for h t and p t.
The Spanish Drug Observatory maintains a drug information system. Its indicator “treatment” is based on data from all treatment starts in public and publicly funded centres. The health coverage is universal and all kinds of treatment are considered. In the present study, we included 169,257 persons who entered first treatment for heroin use from 1991 to 2006 when aged 15 to 54 (mean 28 years), and who started their heroin use between 1971 and 2006 when aged 10 to 44 (mean 21 years). The database was split into two subsets: people who declared injection as the most frequent route for heroin use in the last 30 days before first treatment admission (29%) and people who declared using routes other than injection (68%). Route of administration was missing in 3% of persons.
Assumptions about parameters
As in all applications of back-calculation where a detailed history of individuals is not observed , we have to make some simplifying assumptions, presented below.
First treatment data for the period 1971–1990
In the estimation of p t , first treatment data was available for t between 1991 and 2006, thus restricting estimates to this period. For the preceding part of the study period (t in 1971–1990), we made an educated guess of p t based on general heroin use information in Spain. Based on the first appearances of admissions for heroin use in the emergency units in Spain in 1982 [17–19]; we assumed probabilities of entering treatment (p t ) as low as 0.01 between 1971 and 1981, as there were still no specific treatments available. Thereafter we assumed a linear increase to the value estimated for the parameter p t in 1991.
Mortality for heroin users
Mortality rates for heroin users were only available from two local cohort studies covering the period 1985 to 1999 [20, 21] in an area where injecting was the predominant route of administration . As we did not have better approximation, yearly rates from these studies (minimum 1.4% in 1985, maximum 6.6% in 1995) were extrapolated to the whole country for the corresponding year. For the period 1971 to 1984 a smooth increasing trend from a mortality rate of 1% to 1.4% was applied. Mortality rates from 2000 to 2006 decreased from a rate of 1.5% in 1999, to 1% in 2006. In the analysis by route of administration the same mortality rate was used for injectors, but a constant mortality rate of 1% for non-injectors since they have lower risk .
Degenhardt et al. reported a pooled crude mortality rate of 2.09 per 100 person-years and that mortality risk was increased among out-of-treatment heroin users . In a sensitivity analysis (see analysis section) we considered alternative mortality rates obtained by adding 0.01 to the yearly mortality rates in order to ensure a minimum rate of at least 2%. Note that, in the multi-state model we are imputing mortality rates before first treatment.
Owing to the impossibility of obtaining permanent cessation rates, we looked for lasting cessation rates from long-term cohort studies. As such studies are not available in Spain, we considered yearly cessation rates from a thorough review which reported a range of 0.02-0.04 . Our analyses considered these two extreme values.
For all heroin users and for injectors and non-injectors separately, we applied the aforementioned multi-state model with the Spanish treatment data and the assumed leaving rates to estimate the heroin use incidence (h t , t ranging from 1971 to 2006). As explained, the model also estimates the probability of entering first treatment (p t , t in range 1991 to 2006). We considered the yearly cessation rate of 0.04 and the non-modified mortality rate derived from the local cohort studies. Note that the probability of leaving drug use without having ever been registered for first treatment (q t ) is the sum of the cessation and mortality rates for each year from 1971 to 2006.
In equation 1 when i=j, it means that users began treatment in the same year as they started heroin use. This gives on the average about half a year of observation, and so we must weight μ ij by 0.5.
To assess the fit of the expected incidence values μ ij with their observed values Nij, we have drawn their curves stratified by year of heroin use onset (i).
As results can be dependent on assumptions, a sensitivity analysis was performed to evaluate the two chosen mortality and cessation rates obtaining four combinations of q t , that are reflected in four curves of estimated incidence rates. These combinations were: firstly and as a matter of choice, the available mortality rates and a yearly cessation rate of 0.04; secondly, the same mortality rates and a yearly cessation rate of 0.02; thirdly, the same mortality rates modified by adding 0.01 to the rate for each year and a yearly cessation rate of 0.04; and finally, the modified mortality rates and a yearly cessation rate of 0.02.
Statistical uncertainty was estimated using a bootstrap technique with 500 re-samples, where each re-sample was made up of two parts: 1) the treatment database was re-sampled with replacement and, 2) both our “best guesses” for p t in the period from 1971 to 1990, and the cessation rate for all years were sampled from gamma distributions. The shape and scale parameters were derived from the mean and standard deviation, taking the mean as the “best-guess” value, and the standard deviation was established as 0.01.
The expected number of new heroin users per year was obtained, and converted into rates per 1,000 inhabitants, based on Spanish population yearly census data for people aged 10-44 .
Incidence estimates from the previous study were retrieved to compare with the present estimates. Both the period of years covered and the census figures were the same.
The software used for the statistical analysis was R version 2.13.0 . The study was approved by the Ethics Committee of the Institution with number 2004/1828/I.
Estimates of incidence and probability of entering treatment
Regarding probabilities of entering treatment for the first time, from 1991 to 2000 estimates only varied slightly, with a maximum difference of only 0.02 (data not shown). From 2001 differences increased progressively reaching a final value of 0.2 in 2005 between the lowest and highest combinations of cessation and mortality rates (probability of 0.22 versus 0.42).
Confidence intervals for the various estimates overlap, except for the combination yielding the lowest estimates (results not shown).
Comparison of incidence estimates
Incidence estimates from a previous study  are shown in Figure 6. These estimates had an earlier peak around 1980 and, although they were lower than present ones in the 90’s, they overtook them in the last few years.
We have established a conceptually simple multi-state model to obtain estimates for the incidence rates of heroin use and applied it over a long period in Spain. The highest incidences were observed from 1980 to 1985 corresponding to injectors and a peak in 1990 to non-injectors.
In comparison with previous studies, our estimates are wider in scope since by including mortality and other permanent cessation into the multi-state model it is possible to account for almost all problematic heroin users after drug use onset.
The conceptual model employed in this study focuses on the first phase of a heroin user’s “career”: from heroin use initiation to treatment. Other more complete models based on a theory of compartmental epidemic models over drug user “career” have also been described . Adding more states into the model using the data available would, however, make estimating incidence too complex because heroin use cessation and relapse are frequent and difficult to follow up. Knowing first entry to treatment and first use, we only need to account for quitting heroin use before first treatment through complete recovery (cessation) or death.
As a Spanish anthropologist described, in Spain heroin use had its first phase in the years 1977 and 1978, when the first users became visible, being endemic in the second phase between 1979 and 1982, and reaching its zenith between 1983 and 1986 leading to the institutionalization of the problem . Therefore, trends in overall heroin use incidence obtained in the present study seem reasonable and consistent with previous knowledge about the Spanish heroin epidemic and the HIV-AIDS epidemic . Specifically, the inflection point observed in 1985, when the incidence rate of injected heroin use fell below the rate for non-injection, is consistent with the trend of decreasing HIV incidence among injectors in Spain. However, we observe that estimates from the previous study had higher incidence figures earlier than the present ones (Figure 6). They reflect the fact that the availability of treatment was assumed stable throughout the entire period, leading to high estimates too soon.
We found a decreasing trend in the incidence estimates for the last years observed, which is probably related to the decreasing trend observed in all indicators towards the end of the period studied, as mentioned in the introduction. However, estimates for these last years from the previous study became stable overtaking the estimates from the present one (Figure 6). This is due to the two studies employing different approaches. Equation 1 in the present study was formulated assuming that p t , the probability of entering first treatment, was independent of the era when a person's first heroin use began. Actually, this would be not entirely true if lag time between the drug use onset and first treatment followed a determined pattern, as previous studies assumed [1, 6]. However, if we observe Figure 5, the lag time distribution for the observed values N i,j-i and for the expected values μ i,j-i (j-i represents the lag time), for each year of heroin use onset from 1991 to 2004 all fitted well. So, to modify the equation 1 including the probability of entering first treatment conditional on the initiation of heroin use would be too complex and may not have great practical importance. Therefore, bias can be inherent in both the independence assumption and assuming a determined pattern of lag time, although the direction of such bias cannot be determined.
As observed by de Angelis , results are dependent on assumptions. However, the sensitivity analysis showed that the different incidence curves generated by varying the cessation and mortality rates had similar shapes (i.e. trends) although different levels, suggesting that the model’s estimates were stable. Moreover, we observed that the confidence intervals of incidence figures estimated using non-modified mortality rates and a cessation rate of 0.04 completely contain those of the other three estimates. Thus the chosen incidence estimates do not differ significantly from the other three estimates, which resulted from varying the rates involved.
Concerning the assumptions made about the model parameters, such as 1) mortality rates, 2) permanent cessation, both before first treatment and 3) heroin users that started their first treatment before the observation period, we need to consider possible limitations:
1) Mortality rates were extrapolated by applying to the whole country figures from the North-East of Spain where heroin use injection was more frequent than in the rest of the country . The extrapolation appears to be appropriate, since the period where the highest mortality rates are found for the two cohorts studied (1985 to 1999) coincides with the period when there were more HIV and drug injection related deaths in Spain . However, if the extrapolation is not appropriate it would lead to over-estimation of the total incidence of heroin consumption for the whole country. Note that adding an additional 0.01 to the yearly mortality rates, i.e. to account for the risk of dying when out of treatment being greater, would lead to even greater over-estimations of the incidence.
2) Using lasting cessation rates from long-term cohort studies would overestimate incidence as they include persons with long cessation periods who finally may relapse. On the other hand, the fact that experimental users were not included in studies estimating cessation rates would produce underestimates. Nevertheless, these experimenters are only of anecdotic value for policy interventions.
3) Although the exact dates and figures we have taken for first treatment probabilities prior to the observation period may not apply to the whole country, the values assumed seem plausible as we obtained an increasing sequence of probabilities from 1982 up to the first one estimated based on observed data in 1991, which happens to have a similar slope to that estimated in the observed period.
Besides model building, it is important to consider other limitations related to treatment data, both to its overall availability and to its accuracy. Treatment register data covers public and publicly funded centres, missing people using private treatment centres. This entails a small proportion especially once public substitution treatment centres were widely implemented all over the country following legislation in 1990 . In relation to treatment variables used, we acknowledge the possibility of error in the reported year of heroin use onset, in which we cannot discern any systematic trend except perhaps a certain propensity to round to years ending in 0 or 5.
Incidence trends by route of administration do not necessarily reflect the route used at the time of onset, as the variable was collected referring to the 30 days prior to first treatment. However, in a previous study involving heroin users, both in and out of treatment, and a mean length of use of 10 years, more than 50% did not change their initial route of heroin administration . Thus the study of incidence trends by route of administration in the period immediately previous to first treatment can provide an idea of the different patterns of heroin administration during the heroin epidemic in Spain . The higher probability of entering treatment among individuals declaring injection in the previous month may possibly be related to a change to a more harmful route of administration.
With a simple conceptual model of heroin users’ trajectories related to treatment demand, it has been possible to obtain approximations of heroin use incidence trends. Moreover, different assumptions made do not systematically skew the conclusions. However, enhancing accuracy of drug users’ trajectories and an updating of new treatment admissions will further contribute to better incidence estimates.
Human immunodeficiency virus
Acquired immunodeficiency syndrome
This work was supported by the Fondo de Investigaciones Sanitarias del Instituto Carlos III [PI041783]; Ministerio de Sanidad y Política Social [PNSD 2007-I044] and the Agència de Gestió d’Ajuts Universitaris i de Recerca [AGAUR 2009 SGR 718]. Albert Sánchez-Niubó is supported by ISCIII grant CA08/00214. The authors thank the European Monitoring Centre of Drugs and Drug Addiction (EMCDDA) for support given to projects on drug use incidence. They also are grateful to Dave Macfarlane for English revision, and finally, to clients and professionals from drug treatment centres, as well as professionals from regional and state drug administrations, who have collected and prepared the database.
- Scalia Tomba GP, Rossi C, Taylor C, Klempova D, Wiessing L: Guidelines for Estimating the Incidence of Problem Drug Use. 2008, Lisbon: EMCDDAGoogle Scholar
- Gfroerer J, Brodsky M: The incidence of illicit drug use in the United States, 1962–1989. Br J Addict. 1992, 87: 1345-1351. 10.1111/j.1360-0443.1992.tb02743.x.View ArticlePubMedGoogle Scholar
- Wu LT, Korper SP, Marsden ME, Lewis C, Bray RM: Use of Incidence and Prevalence in the Substance Use Literature: A Review. 2003, Rockville, MD: Substance Abuse and Mental Health Services Administration, Office of Applied StudiesGoogle Scholar
- Hickman M, Seaman S, de Angelis D: Estimating the relative incidence of heroin use: application of a method for adjusting observed reports of first visits to specialized drug treatment agencies. Am J Epidemiol. 2001, 153: 632-641. 10.1093/aje/153.7.632.View ArticlePubMedGoogle Scholar
- Sanchez-Niubo A, Domingo-Salvany A, Gomez-Melis G, Brugal MT, Scalia-Tomba G: Two methods to analyze trends in the incidence of heroin and cocaine use in Barcelona [Spanish]. Gac Sanit. 2007, 21: 397-403. 10.1157/13110444.View ArticlePubMedGoogle Scholar
- Sanchez-Niubo A, Fortiana J, Barrio G, Suelves JM, Correa JF, Domingo-Salvany A: Problematic heroin use incidence trends in Spain. Addiction. 2009, 104: 248-255. 10.1111/j.1360-0443.2008.02451.x.View ArticlePubMedGoogle Scholar
- Nordt C, Stohler R: Incidence of heroin use in Zurich, Switzerland: a treatment case register analysis. Lancet. 2006, 367: 1830-1834. 10.1016/S0140-6736(06)68804-1.View ArticlePubMedGoogle Scholar
- Nordt C, Stohler R: Estimating heroin epidemics with data of patients in methadone maintenance treatment, collected during a single treatment day. Addiction. 2008, 103: 591-7. 10.1111/j.1360-0443.2007.02055.x.View ArticlePubMedGoogle Scholar
- Hunt LG, Chambers CD: The Heroin Epidemics: A Study of Heroin Use in the United States, 1965–75. 1976, New York: SpectrumGoogle Scholar
- Aalen OO, Farewell VT, de Angelis D, Day NE, Gill ON: A Markov model for HIV disease progression including the effect of HIV diagnosis and treatment: application to AIDS prediction in England and Wales. Stat Med. 1997, 16: 2191-2210. 10.1002/(SICI)1097-0258(19971015)16:19<2191::AID-SIM645>3.0.CO;2-5.View ArticlePubMedGoogle Scholar
- Amundsen EJ, Aalen OO, Stigum H, Eskild A, Smith E, Arneborn M, et al: Back-calculation based on HIV and AIDS registers in Denmark, Norway and Sweden 1977–95 among homosexual men: estimation of absolute rates, incidence rates and prevalence of HIV. J Epidemiol Biostat. 2000, 5: 233-243.PubMedGoogle Scholar
- Rossi C: The role of dynamic modelling in drug abuse epidemiology. Bull Narc. 2002, 54: 33-44.Google Scholar
- Rossi C: Operational models for epidemics of problematic drug use: the Mover-Stayer approach to heterogeneity. Socioecon Plann Sci. 2004, 38: 73-90. 10.1016/S0038-0121(03)00029-6.View ArticleGoogle Scholar
- de Angelis D, Hickman M, Yang S: Estimating long-term trends in the incidence and prevalence of opiate use/injecting drug use and the number of former users: back-calculation methods and opiate overdose deaths. Am J Epidemiol. 2004, 160: 994-1004. 10.1093/aje/kwh306.View ArticlePubMedGoogle Scholar
- de la Fuente L, Brugal MT, Domingo-Salvany A, Bravo MJ, Neira-Leon M, Barrio G: More than thirty years of illicit drugs in Spain: a bitter story with some messages for the future [Spanish]. Rev Esp Salud Publica. 2006, 80: 505-520. 10.1590/S1135-57272006000500009.View ArticlePubMedGoogle Scholar
- Delegación del Gobierno para el Plan Nacional sobre Drogas. Observatorio Español sobre Drogas (OED): Informe 2009. Situación y tendencias de los problemas de drogas en España [Spanish]. 2009, Madrid: Ministerio de Sanidad y Política SocialGoogle Scholar
- Cami J, Alvarez F, Monteis J, Caus J, Menoyo E, De Torres S: Heroin as a new cause of toxicological emergencies [Spanish]. Med Clin (Barc). 1984, 82: 1-4.Google Scholar
- Gamella JF: Heroin in Spain (1977–1996). Assessment of a drug crisis [Spanish]. Claves de Razón Práctica. 1997, 72: 20-30.Google Scholar
- Castilla J, de la Fuente L: Trends in the number of human immunodeficiency virus infected persons and AIDS cases in Spain: 1980–1998 [Spanish]. Med Clin (Barc). 2000, 115: 85-9.View ArticleGoogle Scholar
- Orti RM, Domingo-Salvany A, Munoz A, Macfarlane D, Suelves JM, Anto JM: Mortality trends in a cohort of opiate addicts, Catalonia, Spain. Int J Epidemiol. 1996, 25: 545-553. 10.1093/ije/25.3.545.View ArticlePubMedGoogle Scholar
- Brugal MT, Domingo-Salvany A, Puig R, Barrio G, Garcia De Olalla P, De la Fuente L: Evaluating the impact of methadone maintenance programmes on mortality due to overdose and aids in a cohort of heroin users in Spain. Addiction. 2005, 100: 981-989. 10.1111/j.1360-0443.2005.01089.x.View ArticlePubMedGoogle Scholar
- de la Fuente L, Lardelli P, Barrio G, Vicente J, Luna JD: Declining prevalence of injection as main route of administration among heroin users treated in Spain, 1991–1993. Eur J Publ Health. 1997, 7: 421-426. 10.1093/eurpub/7.4.421.View ArticleGoogle Scholar
- Darke S, Degenhardt L, Mattick RP: Mortality amongst illicit drug users: epidemiology, causes, and intervention. 2007, Cambridge, United Kingdom: Cambridge University Press, 1Google Scholar
- Degenhardt L, Bucello C, Mathers B, Briegleb C, Ali H, Hickman M, et al: Mortality among regular or dependent users of heroin and other opioids: a systematic review and meta-analysis of cohort studies. Addiction. 2011, 106 (1): 32-51. 10.1111/j.1360-0443.2010.03140.x.View ArticlePubMedGoogle Scholar
- Amundsen EJ, Bretteville-Jensen AL, Kraus L: A Method to Estimate Total Entry to Hard Drug Use: The Case of Intravenous Drug Use in Norway. Eur Addict Res. 2011, 17 (3): 129-135. 10.1159/000324341.View ArticlePubMedGoogle Scholar
- Instituto Nacional de Estadística: 2010, Madrid: INE, Accessed 18 Jan 2010. Available from: http://www.ine.es,
- R Development Core Team. R: A Language and Environment for Statistical Computing. [2.13.0]. 2011, Vienna, Austria: R Foundation for Statistical ComputingGoogle Scholar
- Domingo-Salvany A, Perez K, Torrens M, Bravo MJ, Anto JM, Alonso J: Methadone treatment in Spain, 1994. Drug Alcohol Depend. 1999, 56: 61-66. 10.1016/S0376-8716(99)00011-3.View ArticlePubMedGoogle Scholar
- De la Fuente L, Barrio G, Royuela L, Bravo MJ, Spanish Group for the Study of the route of heroin administration: The transition from injecting to smoking heroin in three Spanish cities. Implications for control of the HIV epidemic. Addiction. 1997, 92: 1749-1763. 10.1111/j.1360-0443.1997.tb02895.x.View ArticlePubMedGoogle Scholar
- de la Fuente L, Barrio G, Bravo MJ, Royuela L: Heroin smoking by "chasing the dragon": its evolution in Spain. Addiction. 1998, 93: 444-446.PubMedGoogle Scholar
- The pre-publication history for this paper can be accessed here:http://www.biomedcentral.com/1471-2288/13/4/prepub
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.