Participants
A convenience sample of 41 adults with ABI were recruited at the time point of discharge from a specialist brain injury inpatient rehabilitation unit, at an Australian tertiary hospital. Participants were recruited between March 2017 and March 2018, as part of a larger project [24]. Forty-nine individuals involved in the larger study did not complete the surveys, and therefore, are not included in the current analysis. There was a lower proportion of people with severe traumatic injuries in the non-responder group (33%) than in the group who completed the surveys (61%; see Supplement 1). Participants were similar on all other sociodemographic, injury and discharge variables (Supplement 1). Participants were included if they met the following criteria: a newly diagnosed ABI, capacity to provide informed consent (or consent by a substitute decision maker on behalf of the participant) and communication skills to complete a telephone survey (or with the assistance from their substitute decision maker). Ethical clearance for the project was obtained from the Princess Alexandra Hospital Human Research Ethics Committee (HREC/16/QPAH/684; SSA/16/QPAH/685) and the Griffith University Human Research Ethics Committee (2016/915). All methods were carried out in accordance with relevant guidelines and regulations. Participants, or their substitute decision makers, provided informed consent.
Study design
Using a prospective observational cohort design, participants were followed for 6-months, after discharge from specialist inpatient rehabilitation. Health service use was continuously recorded for the 6-month period, with data retrieved via data linkage methods [24]. Telephone surveys were conducted at 6-months, to collect data on a range of outcome measures [24]. All surveys were administered by the same research assistant.
Outcomes measures
The EuroQol-5D (EQ-5D) was used to capture health-related quality of life [25]. The EQ-5D measures health in five dimensions: mobility, self-care, usual activity, pain/discomfort and anxiety/depression. Participants ascribe a rating to each dimension, using a scale ranging from 1 ‘no’ to 5 ‘extreme problems or unable’. [25] The five responses can be converted to a single utility value, by summing the weights associated with each dimension rating, and then subtracting the summed value from one [26]. For example, dimension ratings of 5–5–5-5-5—the worst ratings on each dimension—correspond to weights of 0.274, 0.203, 0.184, 0.335 and 0.289, respectively [26]. The utility score is calculated as: 1–1.285, giving a utility score of − 0.285. Utility values ranged from − 0.285 ‘worst health’ (i.e., ratings of 5–5–5-5-5) to 1 ‘best health’ (i.e., ratings of 1–1–1-1-1) [26].
Psychological wellbeing was measured using the Depression, Anxiety and Stress Scale short-form (DASS-21) [27, 28]. The DASS-21 contains three, 7-item subscales: depression, anxiety and stress [27, 28]. Participants rate the extent to which each item (statement) applied to them in the past week, using a 4-point scale ranging from 0 ‘never’ to 3 ‘almost always’ [28]. Subscale responses are summed and doubled [27]. The internal consistency of the depression (Cronbach’s alpha [α] = .84), anxiety (α = .76), and stress (α = .89) subscales for our sample was acceptable or greater.
Global function (i.e., functional abilities, psychological adjustment and community reintegration) was evaluated using the Mayo-Portland Adaptability Inventory (MPAI-4) [29]. The 29-item MPAI-4 comprises three subscales: ability, adjustment, and participation—which measure interference with daily or valued activities in key areas [29]. Each item is rated on a 5-point scale, ranging from 0 ‘no problem’ to 4 ‘severe problem’ [29]. Subscale responses were summed, with higher scores indicating greater problems. The internal consistency of the abilities (α = .81), adjustment (α = .82) and participation (α = .86) subscales for our sample was acceptable or greater.
The Sydney Psychosocial Reintegration Scale (SPRS-2) captured the change in participants’ participation after injury, in three domains—occupational activity, interpersonal relationships, and independent living skills [30]. Domains comprise 4-items, with each item rated on a 5-point scale ranging from 4 ‘not at all’ to 0 ‘extreme’. Domain responses were summed, with higher scores indicating less disturbance compared to pre-injury [30]. The internal consistency of occupational activity (α = .80), interpersonal relationships (α = .79), and independent living skills (α = .81) was acceptable or greater.
Demographic variables and injury characteristics
Participants’ demographic (i.e., age, gender, indigenous status, employment status at the time of injury), injury (i.e., length of stay, traumatic or non-traumatic, severity) and personal context (i.e., marital status, place of residence at discharge) characteristics, and funding support (i.e., whether in receipt of funding through a national injury insurance scheme, or other government funded support) were retrieved from electronic hospital records. The severity of traumatic injuries was classified according standard guidelines [31], as follows: mild traumatic injury (post traumatic amnesia < 24 h; and/or Glasgow Coma Scale 13–15), moderate traumatic injury (post traumatic amnesia > 1 and < 7 days; and/or Glasgow Coma Scale 9–12) or severe traumatic injury (post traumatic amnesia > 7 days; and/or Glasgow Coma Scale 3–8) [31]. When only Glasgow Coma Scale scores were available (n = 6) the best score in the first 24 h was used.
Participants completed the Functional Independence Measure at discharge [32]. The 18-item Functional Independence Measure comprises a motor and cognitive subscale. The motor subscale includes 13-items that evaluate functions of self-care, bowel and bladder management and mobility. The cognitive subscale comprises 5-items related to communication and social cognition. Each item is rated on a 7-point scale, ranging from 1 ‘total assistance’ to 7 ‘complete independence’. Items across the two subscales are summed, with higher scores indicating greater independence, and therefore, less physical or cognitive disability [32]. A total score is also calculated. We used the motor subscale in our modelling, as an indicator of a person’s level of physical disability [32]. Cognitive and total scores were used for descriptive purposes only.
The index of relative socio-economic advantage and disadvantage (IRSAD) value of participants’ place of residence at discharge was also included [33]. An area with relatively more disadvantage than advantage could be expected to have poorer healthcare resourcing [21]. IRSAD values were constructed from 2016 Census data, collected by the Australian Bureau of Statistics [33]. An IRSAD value summarizes the relative advantage or disadvantage for a given area (i.e., postcode). Higher values (state deciles) indicate an area with a relatively high incidence of advantage and a relatively low incidence of disadvantage [33].
Health service use
Participant’s health service use was determined from sources retrieved through the Queensland Health Data Linkage Group. Data sources were: (i) The Queensland Health Non-Admitted Patient Data Collection; (ii) The Queensland Hospital Admitted Patient Data Collection; and (iii) The Emergency Data Collection. The The Queensland Health Non-Admitted Patient Data Collection captured clinic-based appointments at a public hospital or local health network service, for medical specialist, nursing and allied health services [34]. The The Queensland Hospital Admitted Patient Data Collection recorded acute medical service use—i.e., admitted services, including day-admissions—at a hospital service, at public hospitals, licensed private hospitals and private day surgeries [34]. The Emergency Data Collection recorded re-hospitalizations via data from all Queensland Health emergency departments [34]. A re-hospitalization was defined as at least one overnight hospital stay. The type of service used reflects the predominant nature of a service event.
Whether participants accessed a specialist ABI transitional rehabilitation service provided as part of the ABI rehabilitation continuum of care at the study site, after leaving hospital, was also recorded. If attended, transitional rehabilitation involved 10-weeks of intensive rehabilitation in the community, with multi-disciplinary care provided 2–4 times per week. Therapy was delivered on an individual or group basis, and telehealth was available where suitable. Therapies were goal directed, with a focus on return to meaningful activity (e.g., work, study). Family members and/or carers were encouraged to be involved where possible.
Unmet need for services and service obstacles
Unmet need for services was collected using the statement ‘Have you needed health care [in the past 6-months] but did not get it’, to which a yes/no response was recorded [35]. The Service Obstacles Scale captured information about participants’ difficulties in relation to accessing services, knowledge of and availability of resources, and satisfaction with the quality of care [36]. Comprising three subscales: transport (1-item), finances (1-item) and treatment (4-items), Service Obstacles Scale items are rated on a 7-point scale, where responses range from strongly disagree ‘1’ to strongly agree ‘7’ [36]. Responses on the transport and finances subscales were re-categorized to ‘disagree’ 1–3 and ‘agree’ 5–7 for the analysis, due to limited ratings in some categories. The treatment items were summed (range: 4–28), with higher scores indicating less satisfaction with the quality of care received.
Data analysis
All analyses were performed in R [37]. Missing data were visually inspected (Supplement 2) [38]. Data were assumed missing-at-random, with missing values imputed using predicted mean matching, via the mice package [39]. The results were averaged over five imputed datasets.
There were 10 outcome variables: EQ-5D utility scores, and the subscales from the DASS-21 (depression, anxiety and stress), MPAI-4 (ability, adjustment and participation) and SPRS-2 (occupational activity, interpersonal relationships and independent living skills). All models included the predictor variables: outpatient medical service use, outpatient nursing service use, outpatient allied health service use; medical acute service use; re-hospitalization; transitional rehabilitation service use; age; gender; place of residence at discharge; marital status; employment status at the time of injury; length of stay; injury type; comorbidities; funding support; IRSAD value; Functional Independence Measure motor subscale score; unmet need for services; and Service Obstacles Scale transport, finances and treatment scores.
Our primary interest was to investigate the relationship between health service use variables and EQ-5D utility scores, DASS-21 depression, anxiety and stress scores, MPAI-4 ability, adjustment and participation scores, and SPRS-2 occupational activity, interpersonal relationships and independent living skills scores, with a specific focus of the influence of the covariates, including their interactive effects, on the relationship between service use and these outcome variables. To model such a relationship, we used Bayesian Additive Regression Trees (BART), to allow for interactions between service usage variables and the remaining covariates [40]. BART models allow higher-order interactions to be included if substantiated from the data, resulting in a flexible, but parsimonious, regression model [40]. BART models are also useful in this context, as there are many covariates to consider.
To assess the primary relationship between service use and EQ-5D utility scores, and DASS-21, MPAI-4 and SPRS-2 subscale scores, a heterogenous treatment effects approach was taken [41]. Using a heterogenous treatment effects approach, we are able to investigate whether the average response to a particular treatment is mitigated or exacerbated by certain characteristics in subgroups of the subjects, hence the effect is heterogenous across the dataset or population of interest. These ‘treatments’ were represented by the discrete categorization of service use variables, if necessary, at the first quartile. Discrete categories were: outpatient medical (< 2 and ≥ 2 contacts), outpatient nursing (0 and ≥ 1 contact), outpatient allied health (< 2 and ≥ 2 contacts), medical acute (0 and ≥ 1 contact), re-hospitalization (0 and ≥ 1 contact), transitional rehabilitation (0 ‘did not attend’, 1 ‘attended’). Treatment effects are calculated from the model by considering the difference between each participants’ estimated response and their counterfactual response––i.e., their estimated response having received the alternate ‘treatment’. An individuals’ treatment effect is represented by ‘ τi ’ in Eq. 1. The values ‘ yi(1) ’ and ‘ yi(0) ’ are the responses to treatment and non-treatment respectively. The response ‘ yi(z) ’ will a counterfactual response when ‘ z ’ is not the value observed in the data for subject ‘ i ’.
Treatment effect equation.
$${\tau}_i={y}_i(1)-{y}_i(0)$$
(1)
In the Bayesian framework, parameters are treated as random variables with unknown values [42]. The uncertainty of a parameters value is described by the posterior probability distribution which updates the prior distribution with the likelihood once data has been observed (the posterior is proportional to the product of the likelihood and prior distribution). The prior is a statistical distribution that captures the uncertainty in a population parameter before data collection. Using Bayesian modelling methods, we estimated posterior distributions of these treatment effects at the individual level, which can be aggregated to estimate average treatment effects, or aggregated for certain subgroups to estimate conditional average treatment effects. BART estimates the responses ‘ yi(z) ’ which depend on many covariates, and the treatment variable. The use of BART models made it is possible to assess whether there was heterogeneity in these effects, as the interactions allow particular subgroups of participants to have varying responses to treatment (i.e., health service use). Where no heterogeneity was observed, effects were averaged over all 41 participants.
In interpreting the results, we focused on marginal effects where the 66% credible interval did not include zero. This standard is one suited to the exploratory nature of the paper, but to emphasise the relatively low threshold, we use the terminology ‘weak evidence’ when referring to the results. Our results are likely to be specific to the population for which data was available but provide direction for future research.