Form expert group

A group of MSK experts was recruited through a call for action in the Lancet18 and through identification of research leaders in the literature and across professional networks. The group consisted of rheumatologists, epidemiologists, nurses, physiotherapists, social scientists, a health linguist, patients, statisticians and BoD methodologists.

Establish case definitions and heath states

For each condition, a case definition was chosen on the basis of international recommendations19–23 and adequate alignment with the epidemiological literature to ensure there were sufficient data available to derive the burden estimates.

The main disabling sequelae for each condition were identified through a series of systematic reviews reported elsewhere.24–31 In GBD methodology, the term sequela encompasses the traditional clinical meaning and broader meanings, including severity and disabling consequences of a disease.8 Partitioning a condition into a series of sequelae enables disability variations within a condition to be taken into account through the application of different DWs to the different sequelae.

It was recognised that it was not possible to define every possible set of sequelae and, consequently, a parsimonious set was generated. These sequelae were chosen and defined according to the condition's natural history and the main functional states associated with the condition, as well as the availability of sufficient epidemiological data to estimate prevalence.8

The final list of MSK health states were acute LBP without leg pain; acute LBP with leg pain; chronic LBP without leg pain; chronic LBP with leg pain; acute mild NP; acute severe NP; chronic mild NP; chronic severe NP; mild MSK problems in the legs; mild MSK problems in the arms/mild OA; moderate OA; severe MSK problems in the legs/severe OA; moderate MSK problems in the arms/mild RA; moderate generalised MSK problems/moderate RA; severe generalised MSK problems/severe RA/chronic polyarticular gout; and acute gout.32 The specific health states of these sequelae are described in lay terms in table 2.

In preparation for formulating lay descriptions, the health state of each sequela was described according to a specific set of health domains, which were provided by the GBD Core Team (table 1). Fatigue and manual dexterity were also included by the MSK Expert Group as they are relevant for many MSK conditions. The health domains refer to body functions and structures (eg, lifting), as well as more complex human operations (eg, mobility), but they are not as inclusive as the International Classification of Functioning, Disability and Health. They do not refer to broader aspects of life such as participation, well-being, carer burden or economic impact.

The health state descriptions were designed to reflect the average case for the particular health state in the general population (table 2). This involved several iterations between the MSK Expert Group and the Core Team. Several rules were imposed to ensure consistency across the conditions, including avoidance of the condition's specific name (LBP and NP, however, were acceptable as symptoms); use of simple language; avoidance of emotive language; use of minimum words to convey health state; and inclusion of only the most significant aspects of the health state that would apply to the majority of people living with that condition.

While the process for deriving the health states was a substantial improvement on previous GBD studies, there were some limitations. The health states do not include aspects of life such as participation, well-being, carer burden and economic impact. While this was a deliberate decision, it is important that BoD estimates are supplemented with this information to consider the full impact of a condition in a population.

Establishment of DWs

DWs were derived for each health state in the GBD DW Measurement Survey 2010 (table 2). The Survey has been reported elsewhere.32 The DWs reflect the level of severity of each health state on a continuum from zero (equivalent to full health) to one (equivalent to death). While having short health state descriptions was important to avoid losing the attention of survey participants, a disadvantage was that some descriptions differed only slightly from one another, and, from the clinical perspective, they lacked the depth necessary to capture the full extent of disability and differences between health states. Moreover, while efforts were made on an informal basis to obtain patient input, a more structured approach to patient input, including from a diverse range of cultures from around the world, may have improved the accuracy of the descriptions.

Systematic reviews

A series of systematic reviews were conducted to determine the prevalence, incidence, remission, duration and mortality risk of each of the MSK conditions, as well as population levels of BMD. These have been reported previously.24–31 In brief, for each MSK condition, studies that investigated each of these parameters were identified for inclusion unless they (1) were clearly not representative of the general population (eg, clinic patients, pregnant women, miners); (2) were not population-based (eg, hospital or clinic-based studies); (3) provided no prevalence or incidence data (eg, a commentary piece or risk factor analysis); (4) were limited to a subset of sufferers (eg, for NP only considered those with whiplash); (5) had a sample size less than 150 or (6) were reviews. Searches were conducted in MEDLINE, Embase, CINAHL, CAB abstracts, WHOLIS and SIGLE databases to identify studies published or performed from 1980 to 2009. There were no age, gender or language restrictions. Reference lists of included studies were assessed for additional relevant studies. The potential for publication bias was minimised by conducting extensive reviews, including searches for unpublished data.

Bayesian meta-regression

There were a number of limitations to the prevalence data gathered, with variability between studies regarding case definition, age groupings and prevalence period used. In addition, data were often missing for specific age groups, regions and years of interest.

A Bayesian meta-regression tool, DisMod-MR, was developed by the Core Team to deal with these challenges.36 In brief, DisMod-MR is a tool that helps to pool heterogeneous data presented in different age groups; to adjust data for methodological differences; to check data on incidence, prevalence, duration, remission and mortality risk for internal consistency; and to predict values for countries and regions with little or no data using disease-relevant country characteristics (such as average BMI in OA) and random effects for country, region and super-region. DisMod-MR produced a full set of age/sex/region/year-specific estimates for prevalence, which were used to calculate YLD. For the calculation of TYLDs, ICD codes were mapped to the MSK conditions as shown in table 3.

Severity distributions

The literature search provided enough information on the distribution between mild, moderate and severe disease among prevalent cases of OA and RA that was consistent with the definitions used in the GBD DW surveys. Western Ontario and McMaster Universities Osteoarthritis Index and Health Assessment Questionnaire were the severity measures used for OA and RA, respectively. Random effects meta-analytical methods were used to pool the proportions mild, moderate and severe. This included running a thousand simulations of these proportions assuming β distributions. In each draw of the simulation, the sum of the three proportions to one was scaled and the proportions were multiplied with draws from the results of the DW surveys. This gave 1000 draws of a combined DW for the disorder that were then multiplied with the values of 1000 draws of prevalence to estimate YLDs and uncertainty ranges. For gout, the Expert Group used information from the literature on the frequency and duration of attacks in acute gout and the proportion of cases with gout that have chronic polyarticular gout.

There was no equivalent body of evidence on the range of severity for LBP, NP and the category of other MSK disorders. Instead, the US Medical Expenditure Panel Survey (MEPS) 2000–2009 was used.37 This is a data collection system with new panels being created annually, and each panel is then followed up over a 2-year period to collect information on health service use and related costs. For the purpose of GBD, MEPS contains a large amount of diagnostic information covering 156 of the GBD disease and injury categories, as well as two time points for health status information using the Short Form-12 (SF-12). Conditions are recorded in MEPS as verbatim text and coded to ICD-9CM three-digit codes by professional medical coders. They are recorded for the following reasons: (1) it was reported as a reason for a medical event, that is, inpatient or primary healthcare contact, or the purchase of a drug; (2) it was reported as the reason for one or more disability days; or (3) it was reported as ‘bothering’ the person during the reference period. The first of these reasons is by far the most common source of diagnostic information.

The GBD Core Team used a convenience sample to map SF-12 values to the 0–1 values of the GBD 2010 DWs by asking respondents to fill in a SF-12 questionnaire for 62 lay descriptions across the severity spectrum.

Next, regressions were run on these DWs for each diagnosis with dummies for all other diagnoses to determine the proportion of the combined amount of disability across comorbid conditions that could be attributed to each disorder. The distribution of the ‘comorbidity-adjusted’ DWs for each condition was grouped into categories of no disability and the specified health states by assuming boundaries between the midpoints between the DWs for these health states of increasing severity levels.

The lessons learned from this approach are that comorbidity matters a lot when defining severity and that for most diseases there is a category of no disability, that is, people who have a diagnosis but are not reporting health decrements on SF-12 that can be assigned to the condition of interest. The drawback is, of course, that the information only comes from one country and is, therefore, not representative of the world. Additionally, by applying the same severity distribution across all countries and over time, it implies an assumption that healthcare access does not alter severity of disease.

Final YLD estimates

After estimating YLD for each condition by multiplying prevalence and DW, a correction for comorbidity was applied using simulation methods. This was done because simple addition of YLDs for all conditions would assume that disability is additive if a person has comorbid health states. This could lead to a person with a number of more severe health states having a cumulative DW exceeding one, and thus equating to health loss worse than ‘being dead’. Assuming a multiplicative function between DWs for comorbid health states assures that a combined DW can never be greater than 1. Therefore, hypothetical populations for each country, year, age and sex were created picking draws from the prevalence distributions for all GBD sequelae. For members of these hypothetical populations who were prevalent cases of two or more sequelae, a multiplicative function was used to combine their DWs into a combined DW:

The DW for each condition present in a respondent was then adjusted downwards by the ratio of the combined DW and the original DW for the condition. The average adjustment of DWs for each condition in each age, sex, country and year stratum provided the factor by which to multiply YLDs in each stratum to get comorbidity-adjusted YLD estimates.

Mortality estimates

The GBD measure of mortality is in YLLs, computed as the product of the number of cause-specific deaths multiplied by a loss function dependent on the age at which the death took place.1 Mortality was explicitly modelled for the overall MSK disease category as well as RA and an ‘other MSK diseases’ category (which includes, eg, systemic lupus erythematosus, other autoimmune diseases and osteomyelitis), from 1980 to 2010 for 20 age groups, both sexes, and 187 countries. Cause of death estimates was developed using a comprehensive database of vital registration, verbal autopsy, surveillance and other sources. Ultimately, 2595 country-years of data from 124 countries were used for the estimation of MSK disease mortality.

Total deaths attributable to MSK diseases, deaths attributable to RA and deaths attributable to other MSK conditions were modelled using the Cause of Death Ensemble Modeling (CODEm) strategy, estimating deaths for males and females separately. CODEm is based on five general principles: identifying all available data; maximising the comparability and quality of the data set; developing a diverse set of plausible models; assessing the predictive validity of each plausible individual model and of ensemble models; and choosing the model or ensemble model with the best performance in out-of-sample predictive analysis. Possible models were identified based on statistical methods (mixed effects linear models or spatio-temporal Gaussian process regression applied to the log of death rates or the logit of cause fractions), and combining these with a range of plausible combinations of covariates. Covariates for all MSK disease models included a transformed measure of average income per capita, education level and an aggregate measure of health system access, which takes into account indicators of health service provision such as hospital beds per capita, in-facility deliveries and vaccination coverage rates. The RA model also included alcohol consumption, smoking and obesity because these have been associated with increased risk of RA, vegetable consumption because it has been inversely associated with risk of RA and cholesterol because elevated levels have been observed in RA populations. It was noted that the choice of covariates was based on how well any combination predicts the available data rather than having to pass a strict test of a causal relationship between the factor and the disease outcome.

A large number of possible models combining various statistical techniques and varying combinations of covariates were tested on how well they predicted 30% of the data that were deliberately left out of a model using statistical out-of-sample predictive validity measures. The model or combination of models that best predicted the data was eventually chosen.

Uncertainty in cause of death model predictions was captured using standard simulation methods by taking 1000 draws for each age, sex, country, year and cause. The subcause estimates for RA and the residual category were scaled to the total number of MSK deaths, taking into account the levels of uncertainty associated with each estimate. To create consistency between the sum of cause-specific mortality and all-cause mortality, the models of MSK diseases and all other causes included in the GBD 2010 were rescaled according to the uncertainty around the cause-specific rate.

The ICD category ‘Symptoms, signs and ill-defined conditions’ was not listed as one of the major causes in the GBD classification system. Deaths assigned to this category as well as some other codes used for ill-defined conditions, referred to as garbage codes (GCs), were reassigned to specific causes of death in the GBD classification scheme. This is an important step to allow non-biased comparisons of cause of death patterns across countries and regions. For MSK conditions, the percentage of each group of GCs that have been gone to deaths attributable to all MSK diseases, deaths attributable to RA and deaths attributable to other MSK conditions are shown in table 4. Using empirical algorithms, death due to the MSK conditions increased by 21.1% in total after the redistribution from these GCs (table 5).

Final disease burden estimates

The resulting YLDs were added to YLL to derive the overall burden for the particular condition. The uncertainty interval around each estimate of interest was calculated bounded by the 2.5th and 97.5th centile values of the 1000 iterations of all calculations systematically carrying forward uncertainty around data inputs and data manipulations. Further detail on how uncertainty was calculated can be found elsewhere.7 This is the first time that a comprehensive uncertainty analysis has been carried out as part of GBD estimation, a major advance compared with previous studies. However, there may still be unquantified measurement bias not captured in these uncertainty ranges.

CRA of risk factors

In order to estimate the burden of MSK disease arising from specific risk factors, population attributable fractions (PAFs) were determined based on estimates of exposure prevalence and risk arising from that exposure. In brief, the burden of a disease, for example, a hip fracture, arising from diminished BMD was compared with the hypothetically lower burden that would have occurred if the risk factor distribution had been at some optimal level. This optimal level is the so-called counterfactual exposure distribution that is associated with the lowest population risk and is also referred to as the theoretical-minimum-risk exposure distribution (TMRED). The risk arising from the exposure is based on effect sizes or relative risks (RR) per unit of exposure. Finally, PAFs were calculated for each age, sex, year, cause and country by comparing the current exposure distribution with the counterfactual exposure distribution, accounting for the RR at each exposure level.

For BMD, exposure distributions were modelled from BMD values measured by dual-X-ray-absorptiometry at the femoral neck. RRs for hip fractures and non-hip fractures were based on a meta-analysis that used individual-level data of 12 population-based studies from Western Europe, USA, Canada, Japan and Australia.38 Age-specific and gender-specific TMREDs were calculated as the 90th centile of standardised BMD values from updated NHANES III.39 PAFs of BMD were calculated for hip and non-hip fractures, and these PAFs were multiplied by estimates of the burden of fractures from falls to calculate the burden attributable to low BMD.

For occupation as a risk factor for LBP, occupation was used as a proxy for exposure to occupational ergonomic factors associated with back pain. The economically active population was categorised into seven occupational groups based on data from the International Labor Organization. These occupation groups were consistent with the exposure data in the studies that provided risk information. Clerical and related workers were the reference group (assumed to be the TMRED) in keeping with findings from the literature, and a meta-analysis was conducted to estimate the appropriate RRs for other occupation groups. PAFs and burden estimates were determined in the same way as described for low BMD. Details on the methodology for BMD and occupation are presented in a separate article published in Annals of Rheumatic Diseases.

Adiposity or excess body fat is a risk factor for OA and LBP, and this exposure is measured as BMI in GBD 2010. The exposure distribution for BMI was based on previously published work using measured height and weight from health examination surveys and epidemiological studies.40 Based on analyses of pooled cohort studies, a uniform distribution between 21 and 23 kg/m² and a SD of 1 kg/m² was used. New meta-analyses were performed to calculate RRs for LBP and knee OA per unit of BMI. RRs for hip OA were based on a 2010 meta-analysis of eight cohort studies.41 Using these data, PAFs were calculated for LBP, knee OA and hip OA. Note that a PAF can be interpreted as follows: a PAF of 15% for BMI compared with the TMRED means that 15% of mortality (or burden) in the population can be attributed to increased BMI.