In 2022, England reported its highest number of drug-related deaths on record. 1 Opioids were implicated in almost half of all drug-related deaths, and opioid use disorder (OUD) was a problem for half of all adults accessing community addiction services. 2 Over the past 8 years, between 1 and 2% of all adults accessing community addiction services with OUD have died each year while receiving treatment. 2 Professionals working in community addiction services have key roles in both the delivery of evidence-based treatment and the provision of prognostic information to individuals with OUD. However, despite a good understanding that people with OUD have, on average, standardised mortality rates up to ten times higher than those of the general population, Reference Larney, Tran, Leung, Santo, Santomauro and Hickman3,Reference Bahji, Cheng, Gray and Stuart4 uncertainty regarding individual prognosis and mortality risk presents challenges to addiction services in provision of individuals with accurate personalised risk information, prioritisation of finite resources and appropriate targeting of interventions.
Expansion in the use of clinical informatics and precision medicine has revolutionised the care provided in many healthcare sectors; Reference Castaneda, Nalley, Mannion, Bhattacharyya, Blake and Pecora5 however, development and validation of prognostic risk models in populations of people with OUD has been relatively limited. This is despite multiple systematic reviews having examined individual prognostic risk factors for mortality among people with OUD Reference Larney, Tran, Leung, Santo, Santomauro and Hickman3,Reference Bahji, Cheng, Gray and Stuart4,Reference Sordo, Barrio, Bravo, Indave, Degenhardt and Wiessing6–Reference Lyons, Yule, Schiff, Bagley and Wilens8 and several studies having recently developed models in populations routinely prescribed opioids (e.g. to examine the risk of developing OUD or the risk of opioid overdose). Reference Tseregounis and Henry9,Reference Garbin, Marques and Marques10 To our knowledge, no models have been developed to examine mortality risk in people with a diagnosis of OUD presenting to community addiction services. Models that examine both all-cause and drug-related mortality risks could provide useful information and assistance to both individuals and professionals in making collaborative treatment decisions at the clinically important point of entering addiction treatment.
Potential explanations for the relative paucity of prognostic modelling studies in this area include the required sample size and number of events, and a lack of centralised data repositories that include linked accurate prognostic and outcome information from healthcare and administrative agencies. England is unusual, having recently established and validated a 10-year national data linkage between all individuals presenting to community addiction services and their hospital admission and death records (n > 900 000). Reference Roberts, Doidge, Harron, Hotopf, Knight and White11 This is coupled with the fact that all people in England, regardless of overseas visitor or immigration status, are able to access community addiction services free of charge at the point of delivery; hence, there is a relative absence of a privately funded treatment system. 12 The availability and coverage of this national linked data-set thus provides a rare opportunity to develop and validate adequately powered prognostic models within this population.
Objectives
In this study, we aimed to develop and validate two models, one to estimate 6-month all-cause mortality risk and the other to estimate 6-month drug-related mortality risk for people with OUD on the basis of prognostic indicators routinely recorded during initial assessment at community addiction services in England.
Method
The complete protocol has been previously published, Reference Roberts, Strang, Horgan and Eastwood13 and the study was designed and reported in accordance with the Transparent Reporting of multivariable prediction models for Individual Prognosis Or Diagnosis that use regression or machine learning methods (TRIPOD + AI) statement. Reference Collins, Moons, Dhiman, Riley, Beam and Van Calster14 The completed TRIPOD + AI checklist is provided as supplementary Table 1 available at https://doi.org/10.1192/bjp.2025.10313. The work benefited throughout from input from the South London and the Maudsley Biomedical Research Centre Data Linkage Service User and Carer Advisory Group, which includes experts with lived experience of OUD. Reference Jewell, Pritchard, Barrett, Green, Markham and McKenzie15
Setting
In this study, we used a national English data-set containing linked individual records from three sources: (a) the National Drug Treatment Monitoring System (NDTMS), (b) Hospital Episode Statistics (HES) and (c) the Office of National Statistics (ONS). NDTMS is a centralised database collated and maintained by the Department of Health and Social Care (DHSC) that receives monthly input from all adult statutory community addiction services in England 16 and contains individual-level data on sociodemographic characteristics (date of birth, sex, housing status, etc.), what substances the individual is using problematically, any interventions received and measures of treatment success. HES is a centralised database collated and maintained by the National Health Service (NHS) that collects all information pertaining to NHS in-patient hospital admission in England 17 and covers all NHS in-patient admissions, including any admissions to private or third-sector hospitals subsequently reimbursed by the NHS; it is estimated to contain >99% of all in-patient hospital activity in England. An in-patient hospital admission includes any secondary-care-based activity requiring a hospital bed; this includes day cases and both planned and emergency admissions in physical and mental health settings. HES does not cover accident and emergency (emergency department) attendances, nor out-patient bookings, as these data are held in separate databases. The ONS is a centralised database that contains official death certification records for those individuals that have died. The overall structure of the linked NDTMS-HES-ONS data is clustered with individuals attending one of 150 uniquely commissioned statutory community addiction services across each local authority area in England.
Approval to conduct the linkage analyses was granted under regulation 3 of the Health Service (Control of Patient Information) Regulations 2002, following review by the Caldicott Advisory Panel (ref: CAP-2019-06) and the DHSC Office of Data Protection. NDTMS data were available from 1 April 2013 to 1 April 2022 and contained information for 236 064 unique individuals aged 18 years or over who had attended community addiction services for treatment of problematic opioid use at least once within that timeframe. A standardised clinical history supported by urine drug screen investigation is used within community addiction services to confirm and record a diagnosis of OUD. Linked HES and ONS data are available for individuals; these data include details of any subsequent death records and any individual hospital admissions since HES database inception in 1997. Reference Roberts, Doidge, Harron, Hotopf, Knight and White11 The linked database can only be accessed by DHSC staff working on the project, and all records are stored for a minimum of 5 years after study completion. The studied timeframe deviated from the published protocol, as the end date of the study window was previously specified as 1 April 2023. This deviation was necessary because accurate death record outcome data were not yet available for all individuals presenting to addiction services after 1 April 2022.
Candidate predictor variables
The prognostic indicators for consideration in the multivariable model were initially identified from a systematic search of review articles and their included studies which examined demographic and clinical features associated with increased mortality for people with OUD. Reference Larney, Tran, Leung, Santo, Santomauro and Hickman3,Reference Bahji, Cheng, Gray and Stuart4,Reference Sordo, Barrio, Bravo, Indave, Degenhardt and Wiessing6–Reference Lyons, Yule, Schiff, Bagley and Wilens8,Reference McAuley, Fraser, Glancy, Yeung, Jones and Vickerman18 All candidate predictor variables that were significantly associated with all-cause mortality were extracted, and those which were available within the linked NDTMS-HES-ONS records were discussed with patients and clinicians over a series of three consensus meetings. Given that the aim was to create a model that could be readily incorporated into routine clinical care within time-pressured services, a parsimonious approach was taken to select prognostic indicators, with clinician and patient involvement suggesting that ideally no more than ten variables should be included in the final model. An agreed consensus set of prognostic indicator variables was subsequently extracted from NDTMS-HES-ONS records retrospectively from the time of baseline assessment (i.e. the point at which the patient initially presented for treatment) at the community addiction service, designated time zero (t 0). Following discussion, the protocol initially identified 12 candidate prognostic indicator variables; however, on the basis of additional patient and clinician input during model development, a supplementary candidate variable, the binary of whether a person had ever previously been in addiction treatment, was added as a candidate predictor. Descriptions and structures of the candidate variable can be found in Table 1. 19,20
Table 1 Candidate predictor variables
NDTMS-HES-ONS, National Drug Treatment Monitoring System, Hospital Episode Statistics and Office of National Statistics.
a. Available from 2020 onwards.
b. Problematic use as deemed by the assessing clinician.
c. Variation from protocol which specified that the acute and mental health in-patient hospital admissions be within the past 6 months; following continued patient and clinician input during model development, this was revised to lifetime admission to hospital, as this was thought to be an easier question to ask individuals with opiod use disorder, and the strict timeframe could have led to implementation issues with checking records.
d. Variation from protocol: this predictor was added following patient and clinician input during model development.
Outcome measures
The binary outcomes of all-cause and drug-related mortality were assessed prospectively for each individual at any point up to 6 months after t 0. This timepoint was chosen following clinician and patient feedback, as it was thought to reflect a time horizon of sufficient duration to potentially encourage risk factor modification. For drug-related death, we followed the definition used by the ONS when reporting official national statistics for deaths related to drug-poisoning. The death certificate ICD-10 codes for drug-related deaths can be found in supplementary Table 2; these include codes for deaths due to mental and behavioural disorders caused by drug use (ICD-10: F11–F16, F18–F19), assault (ICD-10: X85), and poisoning of accidental, intentional and undetermined intent (ICD-10: X40–X44, X60–X64 and Y10–Y14).
Sample size
The minimum required sample size for time-to-event model development was based on estimated event rates of prediction model outcomes. Reference Riley, Snell, Ensor, Burke, Harrell and Moons21 Given that the drug-related death event rate was by definition smaller than the all-cause death rate and thus required a larger sample size, this outcome was chosen for sample size calculation. Estimation was performed using the ‘pmsampsize’ command, and, owing to the absence of any reported Cox–Snell R-squared values from previously developed models, we aimed to develop a model with a minimal anticipated Harrel’s C statistic (a measure of discrimination that is similar to the area under a receiver operating characteristic curve but takes account of the censored nature of the data) of 0.70, allowing a maximum shrinkage of 10% to minimise potential overfitting. Reference Ensor, Martin and Riley22 A maximum total of 12 candidate predictors was originally planned, with an estimated event rate based on a previous cohort study which reported 0.0134 drug-related deaths per person year. Reference McAuley, Fraser, Glancy, Yeung, Jones and Vickerman18 The estimated minimum required sample size was 2487 participants and 51 events.
Missing data
Complete outcome data were available for all individuals, and complete candidate predictor information was available for eight variables. Of the five candidate predictors with missing data (injecting status, HIV positivity, hepatitis C RNA positivity, prison referral and accommodation need) the fraction of missing information and its assumed missingness mechanism were assessed for each variable, and all were deemed to at least reasonably fulfil the missing at random assumption (i.e. the probability of a value’s being missing in one variable was not deemed to be related to the probability of missing data in another variable). Missing data were addressed using multiple imputation by chained equations, with the number of imputations set to m = 50 on the basis of the highest fraction of missing information. Reference White, Royston and Wood23 Rubin’s rules were used to combine the results across imputed data-sets. Reference Rubin and Schenker24
Statistical analysis
Multivariable Cox regression models were developed using backwards elimination with the level of alpha for variable exclusion set to 0.157, as recommended on the basis of the Akaike information criterion. Reference Lin and Wei25,Reference Heinze, Wallisch and Dunkler26 Non-linearity of continuous variables was addressed by a multivariable fractional polynomial approach, an established technique for transforming non-linear continuous variables when developing a backwards elimination model. Reference Royston, Moons, Altman and Vergouwe27 Model discrimination was assessed through calculation of Harrel’s C statistic and D statistic (a measure of discrimination, with higher values indicating better discrimination), and we examined calibration curve slopes and the ratio of observed to predicted event rates. Reference Royston, Moons, Altman and Vergouwe27,Reference Riley, Ensor, Snell, Debray, Altman and Moons28 Validation was undertaken using bootstrapping resampling methods, which account for bias due to overfitting more accurately than split sample cross-validation approaches. The model development process was repeated in 1000 bootstrap samples to allow calculation of optimism-adjusted discrimination and calibration measures. Reference Lee and Mark29 Performance was also evaluated by calculation of Harrell’s C statistics for each cluster (i.e. each of the 150 individual statutory community addiction services), and the results were combined using random effects meta-analysis, weighted by the number of events per service. Between-cluster heterogeneity was assessed using the I2 statistic. Reference Riley, Ensor, Snell, Debray, Altman and Moons28 All analyses were conducted in Stata version 18.0 (StataCorp, College Station, Texas, USA).
Results
Data were available for 236 064 people with OUD. There were 2427 deaths due to any cause and 1289 deaths due to a drug-related cause within 6 months of the individual’s most recent presentation to community addiction services in England. Baseline characteristics of the whole sample and those dying due to any or drug-related causes within 6 months are available in Table 2.
Table 2 Baseline characteristics of people with opioid use disorder presenting to community addiction services in England between 1 April 2013 and 1 April 2022
a. Proportion of missing values imputed as variable values can be found in the online supplementary material in Table S4 in the online supplementary material.
Development
Table 3 shows the optimism-adjusted hazard ratios (aHRs) for both final models. In both models, the final model and all variables met the assumption of proportional hazards.
Table 3 Optimism-adjusted hazard ratios (95% CIs) for 6-month all-cause and drug-related mortality in individuals with opioid use disorder presenting to community addiction services in England
a. Modelled as a cubic function in all-cause mortality model.
b. Continuous covariates set at their means and categorical or binary variables set at their reference values.
All-cause mortality
All variables were preserved in the all-cause mortality model except HIV positivity and polysubstance use, leading to a final model event per variable (EPV) rate of 2427/13 = 186.69. Supplementary Fig. 1 shows a graphical representation of the aHRs for the fractional polynomial terms for age; this was treated as a cubic function in the final model. For the variables included in the final model, the optimism-adjusted increase in individual risk of all-cause mortality ranged between a 13% increase (95% CI 2–26%) for having a non-urgent housing problem to a 123% increase (95% CI 34–271%) for having a previous involuntary mental health admission.
Drug-related mortality
All variables were preserved in the drug-related mortality model except HIV positivity, polysubstance use and age, for which fractional polynomial terms did not reach the prespecified significance for inclusion at either the first or the second degree. This led to a final model EPV rate of 1289/12 = 107.42. For the variables included in the final model, the optimism-adjusted increase in individual risk of drug-related mortality ranged between a 12% increase (95% CI 0–27%) for having a previous acute hospital admission to a 187% increase (95% CI 144–238%) for being a person who currently injects.
Supplementary Table 3 shows the complete case analysis (i.e. the results based only on people with complete candidate predictor variable data); the aHRs showed broadly similar trends in both final models using the multiply imputed data.
Validation
Table 4 shows the optimism-adjusted performance statistics for both final models.
Table 4 Mean (95% CI) performance of 6-month all-cause and drug-related mortality prediction models in people with opioid use disorder presenting to community addiction services in England
Discrimination
All-cause mortality
The optimism-adjusted final model explained 20% of the variation in time to all-cause mortality (R 2 D), the D statistic was 1.01 and Harrell’s C statistic was 0.73. Supplementary Fig. 2 shows a forest plot of Harrel’s C statistics across individual community addiction services in England; services with fewer deaths had wider variation in Harrel’s C statistic than services with more events. Four services had no deaths over the studied timeframe and were not included in the meta-analysis. The summary-pooled Harrel’s C statistic was 0.77 (95% CI 0.75–0.79), ranging from 0.66 (95% CI 0.57–0.74) to 0.99 (95% CI 0.98–1.00). The I 2 value (i.e. the percentage of total variation in Harrel’s C statistics explained by between service heterogeneity) was 99.5% (95% CI 96.0–99.8%).
Drug-related mortality
The optimism-adjusted final model explained 21% of the variation in time to drug-related mortality (R 2 D), the D statistic was 1.07 and Harrell’s C statistic was 0.74. Supplementary Fig. 3 shows a forest plot of Harrel’s C statistics. Six services had no drug-related deaths over the studied timeframe and were not included in the meta-analysis. The summary-pooled Harrel’s C statistic was 0.81 (95% CI 0.80–0.82), ranging from 0.63 (95% CI 0.47–0.79) to 0.99 (95% CI 0.99–1.00). The I 2 value was 99.9% (95% CI 99.7–100.0%).
Calibration
All-cause mortality
The optimism-adjusted calibration slope was 1.01 (0.95–1.08). The observed all-cause mortality risk at 6 months was 1.03% (0.99%–1.07%), compared with a mean individual predicted risk of 5.02% (5.00–5.04%).
Drug-related mortality
The optimism-adjusted calibration slope was 1.01 (0.94–1.10). The observed drug-related mortality risk at 6 months was 0.55% (0.52–0.58%), compared with a mean individual predicted risk of 2.19% (2.18–2.20%).
Table 5 depicts clinical examples of 6-month all-cause and drug-related mortality risk for individuals with OUD presenting to community addiction services in England. A web calculator for individualised risk is available at https://connect.calcapp.net/?app=4pekem.
Table 5 Clinical examples of all-cause and drug-related 6-month mortality risk for people with opioid use disorder presenting to community addiction services in England
Discussion
We have developed and internally validated Roberts’ OUD mortality risk using two multivariable prognostic models to estimate 6-month all-cause and drug-related mortality risk for people with OUD presenting for baseline assessment at community addiction services in England. To our knowledge, no previous models have been developed to examine these outcomes in the studied population. Therefore, our models may provide clinically useful information and assistance to both patients and professionals when making treatment and care decisions.
Both models performed well in terms of discriminatory ability, with optimism-adjusted Harrell’s C above 0.7 and D statistics above 1.0. Only around one-fifth of the variation in time to death was explained by each model. Although addition of further predictors has the potential to increase this statistic, we were mindful throughout model development and consultation with clinicians and patients of the need for balance between model parsimony and performance. Models that have excessive numbers of parameters or use parameters that may require complex questioning or interpretation are likely to face implementation problems and thus lack clinical utility within time-pressured services, particularly if risk scoring is completed by professionals with a range of levels of experience and training. Reference Markowetz30 As there was substantial heterogeneity in discriminative ability among services across England, the pooled Harrel’s C statistics should be interpreted with caution owing to the large number of services with a low EPV rate; the pooled estimate is likely to have been an overestimate of discriminative performance. However, it was reassuring that all individual service Harrel’s C estimates remained >0.65 for both models. Patient and clinician consultation before model development suggested that a ten-variable model would be optimal. Both final models contain 11 variables or fewer and are thus roughly in accord with this consensus. Reference Markowetz30 Mean predicted mortality risks were higher in both models than the observed risk; this is probably indicative of the risk reduction associated with being engaged in addiction treatment, in particular, for individuals receiving opiate agonist therapy (OAT). Indeed, recent estimates demonstrate a three and a half times higher drug-related death rate for those not in receipt of OAT compared with those on OAT after adjustment for confounders, Reference McAuley, Fraser, Glancy, Yeung, Jones and Vickerman18 and a fully adjusted post hoc model, which we developed to assess the impact of the binary variable of whether individuals had remained engaged in addiction treatment at 6 months from baseline assessment, demonstrated a significantly decreased risk of both all-cause and drug-related mortality at 6 months if individuals remained in treatment. As these models are intended to calculate risk at baseline assessment in community addiction services – a point at which, by definition, individuals are currently not in receipt of any treatment, including OAT – they allow individualised risk calculation at a clinically important point in time and allow individuals and professionals to contemplate risks without subsequent treatment. Provision of individualised risk information at this juncture may thus promote internally generated behaviour change to address modifiable risk factors for individuals with OUD and increase the likelihood of their remaining engaged in treatment, Reference French, Cameron, Benton, Deaton and Harvie31 as well as assisting professionals in the prioritisation of finite resources such as prescriber availability and ensuring that individuals with elevated risks are actively supported to remain in treatment services. Reference Tai, Kim, Schmeckenbecher, Kitchin, Wang and Kazemi32 There may be concern that providing individualised mortality risk could result in emotional distress for both patients and professionals; however, patient and clinician involvement suggested that this information would be welcomed, compared with the current clinical reality of stating that there is increased risk but with little personalised quantification. Indeed, studies in other settings have demonstrated perceived utility of provision of this type of individualised health information. Reference Nathan and Bhandari33 It is notable that no fractional polynomial of any degree resulted in the inclusion of age in the drug-related mortality model; this is in accordance with previous research that has challenged the ‘ageing cohort’ theory, indicating instead that the observed increase in drug-related deaths does not appear to be driven by age. Reference Lewer, Brothers, Van Hest, Hickman, Holland and Padmanathan34
This study had several strengths, including the comprehensive, non-selective and national nature of the data-set and the substantial involvement of clinicians and patients from the outset and throughout model development, validation and interpretation. The pre-publication of the protocol, alongside documentation of any deviation, and reporting of the study in accordance with the TRIPOD + AI statement also are notable strengths. Reference Collins, Reitsma, Altman and Moons35 However, there were also several potential limitations. All prognostic indicator variables were collected retrospectively from an administrative data-set, the underlying data for which was supplied by addiction treatment, hospital and death registration services. There was therefore potential information bias and risk of lack of availability of some predictor variables, as noted particularly with HIV and Hepatitis C, if submitted documentation was incomplete. Although relying on routinely documented clinical information as the source of prognostic information has limitations, this approach has been used frequently and reflects how the model would be used in clinical practice; for instance, some information may not be available to professionals or patients at the time of baseline assessment. Reference Hippisley-Cox, Coupland, Vinogradova, Robson, May and Brindle36 The model will benefit from independent validation in other samples and subsequent examination of its utility in clinical practice and its acceptability among professional and patient groups; potentially suitable data-sets have been identified in both Wales and Australia. Reference Larney, Jones, Fiellin, Nielsen, Hickman and Dobbins37,Reference Trefan, Akbari, Paranjothy, Farewell, Gartner and Fone38 Given this, we have ensured reporting of model baseline survivor and linear predictor values. Continued co-production through independent validation and implementation with clinicians and patients will remain a key requirement.
Despite significant expansion and understanding in the use of machine learning methods to develop prognostic models across healthcare sectors, initial patient and clinician feedback has demonstrated reticence with respect to employing these in the context of mortality prediction in OUD. The perception of a ‘black box’ or lack of transparent understanding of what prediction outcome scores are based on, as well as the relative infancy of clinical informatics within the OUD field, led to concerns about clinical utility and implementation within community addiction services. Reference Rudin39 Clinicians were comfortable with clinical risk tools developed using classical statistical methods and their corollaries used in other areas of healthcare Reference Hippisley-Cox, Coupland, Vinogradova, Robson, May and Brindle36 and welcomed their potential expansion within addiction settings. However, there was concern among patients that results from machine learning methods would not be believed and that explanation of algorithms could create difficulties in conveying the predictive information to individuals accessing services. As such, traditional statistical methods were chosen to develop this model.
Standardised all-cause and drug-related mortality rates are significantly elevated among people with OUD, and, despite a significant body of literature describing individual prognostic risk factors, often clinical judgement alone is used to consider individual prognosis and the prioritisation of treatment interventions in addiction treatment services. Whereas other areas of medicine routinely incorporate risk tools into care to assist clinical decision-making, Reference Hippisley-Cox, Coupland, Vinogradova, Robson, May and Brindle36,Reference Olesen, Torp-Pedersen, Hansen and Lip40 clinical informatics within the addiction field has been somewhat slower to progress. Given the significant elevated mortality risks within this population, the development and validation of individualised prediction models that demonstrate good optimism-adjusted discrimination and calibration is timely, warranted and urgent. This is the first stage of model assessment. Independent validation and demonstration of the clinical utility of both models are necessary next steps to ensure buy-in from professionals, policy makers and patients and for the models to be valued and successfully implemented.
Supplementary material
The supplementary material is available online at https://doi.org/10.1192/bjp.2025.10313
Data availability
The data-sets generated and analysed during the current study are not publicly available as they contain sensitive patient-identifiable data. Although access to the linked data-set is only available within the DHSC subject to approval, aggregated extracts of NDTMS data are publicly available at www.ndtms.net, and extracts of HES-ONS data are available through the Data Access Request Service at NHS England. The full published protocol and statistical analysis plan are also available. Reference Roberts, Strang, Horgan and Eastwood13
Author contributions
All authors met all criteria for authorship according to the ICMJE recommendations. Contribution of article co-authors as per the Contributor Roles Taxonomy (CRediT) author statement were as follows. E.R.: conceptualisation, methodology, formal analysis, investigation, data curation, writing – original draft, writing – review and editing, visualisation, and project administration; J.S.: writing – review and editing, and supervision; E.T.: writing – review and editing; J.C.: methodology, project administration, and writing – review and editing; T.L.: methodology, project administration, and writing – review and editing; C.A.: methodology, project administration, and writing – review and editing; B.E.: conceptualisation, methodology, data curation, writing – review and editing, and supervision.
Declaration of interest
All authors have completed the ICJME Unified Competing Interest form (available on request from the corresponding author). E.R. is a member of the British Journal of Psychiatry editorial board; he did not take part in the review or decision-making process for this paper.
Funding
This paper represents independent research funded by the corresponding author’s National Institute of Health and Care Research Advanced Fellowship (NIHR302215). The funders had no contribution to the study design; the collection, analysis and interpretation of data; the writing of the report; or the decision to submit the article for publication. All authors were independent from funders, had full access to all of the data (including statistical reports and tables) in the study and take responsibility for the integrity of the data and the accuracy of the data analysis. The views expressed are those of the authors and not necessarily those of the NHS, the NIHR or the DHSC.
Open access
eLetters
No eLetters have been published for this article.