The implications of competing risks and direct treatment disutility in cardiovascular disease and osteoporotic fracture: risk prediction and cost effectiveness analysis

Competing mortality risk and multimorbidity

Most risk prediction tools do not account for competing mortality risk. Competing mortality risks occur when someone is at risk of death from an unrelated condition to the one being studied (e.g. death from Parkinson’s disease in a study of ischaemic heart disease). Competing mortality risk is more common among people with multimorbidity than those with single conditions, and more common among older people than younger people. 24 Risk prediction models can be designed to account for competing mortality risk, but most do not. In the context of this study, QRISK325 and QFracture-201226 do not account for competing mortality risk. QRISK®-Lifetime (ClinRisk Ltd) does account for competing mortality risk,27 as does FRAX22,28 (but FRAX cannot be externally validated because the risk prediction algorithm is not publicly available).

Competing risk is a well-recognised problem in survival analysis. 29–33 In essence, conventional survival analysis and risk prediction tool development treats death from unrelated conditions as a censoring event equivalent to loss to follow-up for other reasons. The underlying assumption is that censored individuals have the same risk of the outcome being examined as those remaining in the study. Although this may be reasonable for ordinary loss to follow-up, this is clearly not true for people who have died. The consequence is that prediction models that do not account for competing risk in their development overestimate the risk of the outcome among those with competing risks. 29,31,32

Studies comparing conventional models with models that are adjusted for competing risk have shown that conventional models, on average, overpredict CVD risk,27,34–36 but the effect of competing risks has not been systematically examined in the context of fracture risk prediction (although the NICE bisphosphonate guidance23 acknowledges that QFracture-2012 and FRAX risk predictions are difficult to compare as a result). Of note is that conventional external validation of risk prediction tools such as QRISK3 and QFracture-2012 also does not account for competing mortality risk, so the excellent observed calibration in external validation is in the context of making the same assumptions about competing mortality risk as prediction tool derivation. 16,18,19

Therefore, the risk prediction modelling element of this project is concerned, first, with external validation of existing risk prediction tools recommended by NICE (QFracture-2012) or being considered by NICE (QRISK3), accounting for competing mortality risks and examining performance in important subgroups [objective 1, whereby examination of the subgroups of people with diabetes and chronic kidney disease (CKD) was funded by the costed extension to the original grant], and, second, with the derivation and internal validation of new risk prediction tools that account for competing mortality risk in derivation (objective 2). In addition, funded by the costed extension to the original grant, we have completed an external validation of the QRISK-Lifetime tool (objective 3).

Data sources

Data used in this study were taken from CPRD GOLD,66 which derives data from general practices using INPS Vision electronic health records and is distinct from the derivation data set, which is derived from practices using the EMIS system. Identical to QRISK3 derivation and internal validation, patients were eligible for inclusion if they:

• were permanently registered with a general practice, contributing up-to-standard data for at least 1 year and with consent to link GP data to hospital discharge (Hospital Episodes Statistics Admitted Patient Care) and mortality [Office for National Statistics (ONS) mortality registration] data

• were aged ≥ 25 years and < 85 years

• had no prior history of any CVD

• had no prior history of statin prescription.

Cohort entry was defined as the latest date of an individual’s date of registration plus 1 year, the individual’s 25th birthday or 1 January 2004. Cohort exit was defined as the earliest of:

• the first non-fatal or fatal cardiovascular event

• receipt of a statin prescription

• deregistration from their participating general practice

• date of last data collection from their participating general practice

• end-of-study follow-up on 31 March 2016.

Of note, the QRISK3 derivation and internal validation did not censor on statin prescription, but we chose to, as the primary purpose of the tool is to inform decisions on statin initiation.

Sample size

The sample size is fixed by the size of the CPRD GOLD data set. Therefore, no formal power calculation was carried out, as it could not alter study design and the available sample size was considered sufficient for the purpose. 54

Outcome definition

The outcome was the first CVD event experienced by an individual, defined as the earliest GP, hospital or mortality record of non-fatal coronary heart disease, ischaemic stroke or transient ischaemic attack (TIA). Outcomes were defined using Read codes (for GP data) and International Statistical Classification of Diseases and Related Health Problems, Tenth Revision (ICD-10) codes (for hospital discharge and mortality data). ICD-10 codes are those listed in the published QRISK3 derivation paper;25 however, there are no published Read codesets available. We, therefore, derived our own Read codeset, and this and the ICD-10 codes are listed in appendices to the published paper67 reporting our external validation of QRISK3.

Other variable definitions

At the time the analysis was done, there were no published Read codesets for other variables included in the QRISK3 model, which we defined using Read codes in GP data which we created for this study and values [e.g. systolic blood pressure (SBP), cholesterol] in GP data (listed in our published paper). 67 There were several data-handling variations compared with the original QRISK3 derivation:

• We chose a later cohort entry date (1 January 2004 vs. 1 January 1998).

• Where no cholesterol value was available at baseline, then QRISK3 allowed cholesterol values after the cohort entry date to be used provided that they were before any CVD event. In contrast, we used values from before the cohort entry date only to avoid using future information in prediction.

• CPRD makes only group Townsend deprivation scores available as vigintile (equal 20th) of Townsend score. Therefore, we estimated the median Townsend score of national vigintiles and used that in prediction.

In addition, we calculated a mCCI based on Read codes in the GP data, using a published codeset. 65 The mCCI was modified in that the original Charlson Comorbidity Index (CCI) includes several CVDs that are, by definition, excluded at baseline in this study because participants are CVD free at baseline. The mCCI was not used in prediction but was used to examine discrimination and calibration in subgroup analysis (categorised as 0, 1, 2 and ≥ 3), along with age group (categorised as 25–44 years, 45–64 years, 65–74 years and 75–84 years).

Missing data

Missing data handling and the proportion of each variable with missing data are shown in Appendix 1, Table 24. As with QRISK3 derivation, patients were excluded if the Townsend deprivation score was missing, patients with missing data on ethnicity were assumed to be white and patients with no record of a condition were assumed not to have the condition. For continuous variables [i.e. body mass index (BMI), total cholesterol : high-density lipoprotein ratio (TC : HDL), SBP, SBP variability] and for smoking status, multivariate imputation via chained equations68 was used to generate five imputed data sets. Analyses of these imputed data sets were combined using Rubin’s rules to account for the uncertainty association with imputation. 69

Analytical methods: external validation

The published QRISK3 2017 prediction model was implemented (under GNU Lesser General Public Licence v3), and the predicted 10-year risk of experiencing a CVD event was calculated for each patient without recalibration of baseline risk. Model performance was evaluated by examining discrimination and calibration.

Discrimination evaluates how well the risk score differentiates between patients who experience a CVD event (or more generally, the event of interest) during the study and patients who do not. We primarily examined discrimination using the truncated version of Harrell’s c-statistic to include only pairs where the earliest survival time was no later than 10 years after entry. A c-statistic of 0.5 indicates that the risk score performs no better than chance, whereas a c-statistic of 1 indicates perfect discrimination. Evaluating how good discrimination is for values between 0.5 and 1 is arbitrary and involves judgement. We considered c-statistic values of 0.5–0.599 poor, values of 0.6–0.699 moderate, values of 0.7–0.799 good and values > 0.8 excellent.

Two additional measures of discrimination were calculated. First, we calculated Royston and Sauerbrei’s D-statistic (based on the separation in event-free survival between patients with predicted risk scores above and below the median, where higher values indicate greater discrimination and a difference of ≥ 0.1 is suggested as indicating a meaningful difference in discrimination). 70 Second, we calculated a related R²-statistic designed for estimating explained variation in censored survival data. 71

Models may have good discrimination but imperfectly predict risk, for example by systematically overpredicting or underpredicting. Examination of calibration is, therefore, important, particularly where predicted risk is used to determine offers of treatment, as it is for primary prevention of CVD. 10 Calibration refers to how closely the predicted and observed probabilities agree at group level, and for this purpose participants were grouped into 10 equal-sized groups (deciles) of predicted risk. Calibration of the risk score predictions was assessed by plotting observed proportions against predicted probabilities. For both men and women separately, plots were generated for all patients and for prespecified subgroups of age, mCCI, diabetes type and CKD, based on summary statistics pooled across the imputed data sets. Subgroups were defined to ensure that there were enough events in each subgroup to ensure stable estimates of observed risk (and for diabetes and CKD, analysis was, therefore, in the whole subgroup without further stratification for age group or mCCI). CKD was defined in two ways: (1) only using Read codes,67 as per QRISK3 derivation25 and (2) using the same set of Read codes or the last recorded eGFR or eGFR based on last recorded serum creatinine, where an eGFR < 60 ml per minute defined CKD.

The following summary statistics and their standard errors (SEs) were obtained by decile of predicted risk score and for each imputed data set, in turn: non-parametric measures of observed risk or proportions of patients with a CVD event, the Kaplan–Meier estimator (the conventional measure ignoring competing risks), the Aalen–Johansen estimator (an extension to allow for competing events, non-CVD death in this case)16 and the mean predicted risk score. All models were fitted in R 4.0.0 (The R Foundation for Statistical Computing, Vienna, Austria) and Stata® 11.2 (StataCorp LP, College Station, TX, USA).

Analytical methods: competing risk model derivation and internal validation

Competing risk model derivation and internal validation was carried out in the same data set as QRISK3 external validation. For this purpose, participants were randomly allocated to distinct derivation and test data sets in a 2 : 1 ratio, with allocation balanced in terms of age and final event status. The derivation data set was used to derive CRISK, that is a Fine–Gray model to predict the 10-year risk of experiencing a CVD event, accounting for the competing risk of non-CVD death. Separate models were estimated for men and women. Reflecting the overall aim of the project, where we wished to explicitly compare prediction in models accounting for competing risk compared with ignoring competing risk, we included all of the same main effects (i.e. predictors) and age interactions as QRISK3, modified as follows. First, we accounted for non-CVD death as a second (competing) outcome using the Fine–Gray model, and we re-estimated fractional polynomial terms for continuous variables, including in QRISK3, selecting terms based on those performing best (as measured by the c-statistic) in balanced 10-fold cross-validation and showing consistency of model fit [i.e. Akaike information criterion (AIC)] across folds of the derivation data set (this model is called CRISK). Second, as QRISK3 predictors are focused on CVD events, we derived a further model [i.e. the competing mortality risk model with Charlson Comorbidity Index (CRISK-CCI)], which additionally included the mCCI score in the model (categorised as 0, 1, 2, ≥ 3), as CCI is a well-validated predictor of total mortality. 12 Fine–Gray models allow the cumulative incidence function or probability of a CVD event occurring over time to be directly predicted; however, the subdistribution hazard ratios (HRs) in the Fine–Gray models do not have a straightforward interpretation, as they describe the direction but not the magnitude of the effect of predictors on the cumulative incidence function. The use of fractional polynomials and the inclusion of complex interactions with age further complicate direct interpretation of model coefficients. Model coefficients are, therefore, not straightforwardly interpretable, but the derived model is provided in Appendix 1, Tables 25 and 26, to allow replication.

The performance of all three models (i.e. CRISK, CRISK-CCI and QRISK3) was evaluated in the independent validation data set by examining discrimination and calibration, as described above. R 4.0.0 was used for all analyses.

Type 1 diabetes

There were 6025 women with type 1 diabetes potentially eligible for inclusion, of whom 646 (10.7%) and 1627 (27.0%) were excluded because of prior CVD or prior statin prescribing for primary prevention, respectively. There were 8260 men with type 1 diabetes potentially eligible for inclusion, of whom 953 (11.5%) and 2464 (40.9%) were excluded because of prior CVD or prior statin prescribing for primary prevention, respectively. Therefore, 3752 women (62.3% of potentially eligible women with type 1 diabetes) and 4843 men (48.6% of potentially eligible men with type 1 diabetes) were included in analysis of type 1 diabetes.

During follow-up of the type 1 diabetes cohort, there were 108 CVD events in 13,098 person-years’ follow-up in women [8.25 (95% CI 6.83 to 9.94) events per 1000 person-years] and 172 CVD events in 15,824 person-years’ follow-up in men [10.90 (95% CI 9.40 to 12.60) events per 1000 person-years].

Discrimination in people with type 1 diabetes was excellent [Harrell’s c-statistic: women, 0.830 [95% CI 0.768 to 0.891; men, 0.853 (95% CI 0.803 to 0.902)] and explained variance in the model was 51.6% in women and 48.0% in men (see Appendix 1, Table 29).

Ignoring competing risks (see Appendix 1, Figure 28, parts a and c), calibration in women with type 1 diabetes was good (allowing for small number of events and, therefore, a relatively noisy plot), but there was some overprediction in men at higher predicted risk. Accounting for competing risks (see Appendix 1, Figure 28, parts b and d), there was overprediction in women at higher predicted risk and greater overprediction in men at higher predicted risk.

Type 2 diabetes

There were 53,284 women with type 2 diabetes potentially eligible for inclusion, of whom 12,068 (22.6%) and 24,194 (45.4%) were excluded because of prior CVD or prior statin prescribing for primary prevention, respectively. There were 68,236 men with type 2 diabetes potentially eligible for inclusion, of whom 19,777 (28.9%) and 27,382 (40.1%) were excluded because of prior statin prescribing for primary prevention, respectively. Therefore, 24,194 women (i.e. 31.9% of the potentially eligible women with type 2 diabetes) and 21,077 men (i.e. 30.9% of the potentially eligible men with type 2 diabetes) were, therefore, included in analysis of type 2 diabetes.

During follow-up of the type 2 diabetes cohort, there were 1167 CVD events in 44,678 person-years’ follow-up in women [26.12 (95% CI 24.68 to 27.64) events per 1000 person-years] and 1682 CVD events in 57,160 person-years’ follow-up in men [29.40 (95% CI 28.10 to 30.80) events per 1000 person-years].

Discrimination in people with type 2 diabetes was moderate to good [Harrell’s c-statistic: women, 0.741 (95% CI 0.722 to 0.760); men, 0.695 (95% CI 0.679 to 0.712)] and explained variance in the model was lower than for type 1 diabetes (i.e. 29.2% in women and 22.0% in men) (see Appendix 1, Table 29).

Ignoring competing risks (see Appendix 1, Figure 29, parts a and c), calibration in women with type 2 diabetes was good, but there was some overprediction in men at the highest predicted risk. Accounting for competing risks (see Appendix 1, Figure 29, parts b and d), there was progressively increasing overprediction in women with moderate to high predicted risk and in men in all but the lowest deciles of predicted.

Chronic kidney disease defined by Read code alone

There were 16,048 women with CKD defined by Read code alone potentially eligible for inclusion, of whom 4223 (26.3%) and 4897 (30.5%) were excluded because of prior CVD or prior statin prescribing for primary prevention, respectively. There were 15,784 men with CKD defined by Read code alone potentially eligible for inclusion, of whom 5645 (35.7%) and 3850 (24.4%) were excluded because of prior CVD or prior statin prescribing for primary prevention, respectively. Therefore, 6918 women (i.e. 43.1% of the potentially eligible women with CKD defined by Read code alone) and 5659 men (i.e. 35.9% of the potentially eligible men with CKD defined by Read code alone) were, therefore, included in analysis of CKD defined by Read code alone. The mean age of women was 63.0 years and mean age of men was 59.2 years.

During follow-up of the CKD defined by Read code alone cohort, there were 541 CVD events in 25,544 person-years’ follow-up in women [21.18 (95% CI 19.48 to 23.02) events per 1000 person-years] and 569 CVD events in 21,459 person-years’ follow-up in men [26.50 (95% CI 24.40 to 28.80) events per 1000 person-years].

Discrimination in people with CKD defined by Read code alone was good [Harrell’s c-statistic: women, 0.755 (95% CI 0.728 to 0.782); men, 0.734 (95% CI 0.708 to 0.760)] and explained variance in the model was 34.2% in women and 29.7% in men (see Table 30).

Ignoring competing risks (see Appendix 1, Figure 30, parts a and c), calibration in women with CKD defined by Read code alone was reasonable (allowing for small number of events and, therefore, a relatively noisy plot) and good for men (with some underprediction for both at higher predicted risk). Accounting for competing risks (see Appendix 1, Figure 30, parts b and d), there was overprediction in women at moderate and higher predicted risk and overprediction in men at higher predicted risk.

Chronic kidney disease defined by Read code and estimated glomerular filtration rate

Laboratory values were extracted for only people included in the CVD study cohort and so it is not possible to calculate the proportions excluded because of prior CVD or prior statin prescribing. There were 71,094 women and 33,699 men with CKD defined by Read code or eGFR included in analysis, with an older mean age than the cohort with CKD defined by Read code alone (mean age: women, 70.1 years; men, 69.1 years).

During follow-up of the CKD defined by Read code or eGFR cohort, there were 8877 CVD events in 348,982 person-years’ follow-up in women [25.44 (95% CI 24.92 to 25.96) events per 1000 person-years] and 5273 CVD events in 146,730 person-years’ follow-up in men [35.90 (95% CI 35.00 to 36.90) events per 1000 person-years].

Discrimination in people with CKD defined by Read code or eGFR was somewhat worse than for people with CKD defined by Read code alone. Discrimination was moderate to good [Harrell’s c-statistic: women, 0.705 (95% CI 0.699 to 0.712); men, 0.671 (95% CI 0.663 to 0.680)] and explained variance in the model was somewhat lower than for CKD defined by Read code alone (i.e. 24.9% in women and 17.4% in men) (see Table 30).

Ignoring competing risks (see Appendix 1, Figure 31, parts a and c), calibration in women with CKD defined by Read code or eGFR was excellent, but there was some underprediction in men. Accounting for competing risks (see Appendix 1, Figure 31, parts b and d), there was progressively increasing overprediction in women and men with moderate to high predicted risk.

QRISK3 external validation

At the whole-population level, QRISK3 has excellent discrimination (which is the ability of the model to distinguish people at higher or lower risk). However, as is expected when examining discrimination in subsets of the modelled population defined by strong predictors of the outcome,72 discrimination was poor to moderate when stratified by age and additionally worse when stratified by level of comorbidity (which was not a predictor in the model). Calibration is the extent to which predicted and observed event rates are similar, and it was excellent in the whole population when ignoring competing mortality risks; however, there was systematic underprediction after competing risks were accounted for. Calibration was considerably worse in older people and in people with higher levels of comorbidity, where QRISK3 systematically overpredicted risk, particularly after competing mortality risks were accounted for.

In people with diabetes, discrimination was excellent in type 1 diabetes and moderate to good in type 2 diabetes. Similar to the whole population, calibration was good, with some overprediction when ignoring competing risks, but there was more consistent overprediction once competing risks were accounted for. Similar findings were found for people with CKD, but it is important to recognise that the populations studied exclude people with prior statin prescribing, which excludes substantial numbers of people with either condition (based on statin prescribing in the primary prevention population, 27.0% of women and 49.9% of men with type 1 diabetes, 45.4% of women and 40.1% of men with type 2 diabetes, and 30.5% of women and 24.4% of men with CKD defined by Read code were excluded).

The published external validation of QRISK2 found excellent discrimination and calibration at the whole-population level when ignoring competing mortality risk (i.e. answering the question ‘what is the risk of CVD assuming this person does not die of anything else in the next 10 years?’). 15 This study found similar, but additionally found overprediction and poor calibration in people aged 75–84 years, and moderate calibration in people aged 65–74 years and in people with the highest levels of comorbidity (mCCI = 3).

Once competing mortality risk was accounted for (i.e. answering the question ‘what is the risk of CVD allowing for the risk of death from something else first?’), then there was greater overprediction at the whole-population level, and particularly in older people and in people with more comorbidity. These findings are consistent with other studies examining the impact of competing risks on estimated CVD risk in people without CVD34,73,74 and with established CVD. 75

QRISK2 has also been shown to systematically overpredict CVD risk in a contemporary population of people with type 2 diabetes, with increasingly poor discrimination with increasing age76 and underprediction in a contemporary population of people with type 1 diabetes. 60 This highlights that good performance at the whole-population level does not necessarily mean good performance in important subgroups,72 and also that models derived in populations excluding prior statin prescribing are likely to be increasingly unrepresentative as statin prescribing increases.

CRISK and CRISK-CCI derivation and internal validation

The two new competing risk models derived (i.e. CRISK and CRISK-CCI) has similar excellent discrimination for CVD events as QRISK3. CRISK was better calibrated than QRISK3 (after accounting for competing mortality in derivation, but without adding any new predictors to the model) and CRISK-CCI was better calibrated again (after adding the mCCI as an additional predictor).

Two studies74,77 in people aged ≥ 65 years have examined the impact of competing mortality risk on CVD prediction. Like this study, the two studies74,77 also found only moderate discrimination of whole-population CVD risk prediction tools in older adults, and that newly derived competing risks models were generally better calibrated than models derived using standard Cox regression. 77 In a UK study35 evaluating a new competing risk model against the QRISK2, differences between predicted and observed CVD risk were greatest among people with highest predicted risk, as was found in this study.

Limitations

Limitations of this study are largely those that are found in all studies using routine GP data, including the original QRISK3 derivation. 20 First, there is considerable missing data for key predictors. As with QRISK3 derivation, we used multiple imputation for missing data, but the assumption that data are missing at random is a strong one because risk factors are likely to be better recorded in people at higher CVD risk. 35 This weakness is balanced against the use of more representative population data than is found in individually recruited research cohorts where data are more complete. Second, we used a more recent index date (1 January 2014) than QRISK3 (1 January 1998), which likely means that we exclude more people with prior statin prescribing. Deriving clinical prediction tools on increasingly historical data is likely biased because CVD incidence is falling,73 but using more recent data with greater rates of exclusion because of prior statin initiation may also be biased.

Data sources, outcome definition, other variable definitions and missing data

The same data set used for QRISK3 external validation was used (see Chapter 2, Methods), with the same variation in methods applying (i.e. a cohort entry date of 1 January 2004 vs. 1 January 1998 for QRISK-Lifetime, our implementation did not allow the use of future cholesterol values and Townsend deprivation scores were fitted as the median of vigintiles of Townsend score).

Analytical methods

The lifetime (i.e. to age 95 years) and 10-year risk of experiencing a cardiovascular event was calculated for each patient using publicly available QRISK-Lifetime 2011 software (under GNU Lesser General Public Licence version 3) without recalibration. As lifetime risk is not observed in the validation data set, the performance of the risk score was assessed by examining discrimination and calibration of the model using the same methods as for QRISK3 external validation (see Chapter 2, Methods) over a 10-year time horizon, as was carried out in the original derivation and internal validation study. 27 For both men and women, calibration was evaluated in the whole population and in prespecified subgroups of age and mCCI. Calibration refers to how closely predicted risk and observed probabilities agree at group level. As QRISK-Lifetime accounts for competing mortality risk, we evaluated calibration using only the Aalen–Johansen estimator of observed risk in censored survival data (i.e. an extension of the Kaplan–Meier estimator, which allows for competing events, non-CVD death in this case). 84

Clinical guideline recommendations for primary preventative treatment of CVD classify patients in relation to thresholds of predicted risk. In England and Wales, NICE-recommend treatment if 10-year predicted CVD risk is ≥ 10%. Consistent with the validation of QRISK-Lifetime over a 10-year time horizon, we examined changes in which patients were recommended for treatment based on either a QRISK3 or QRISK-Lifetime 10-year predicted risk of ≥ 10%.

However, there is no recommended threshold of lifetime risk at which to offer treatment,10 although NICE has signalled that consideration of lifetime risk models is of interest in any future guideline update. 13 Therefore, we additionally calculated the proportion of men and women recommended for treatment by QRISK3 at the 10% threshold, and used a cut-off point of QRISK-Lifetime-predicted lifetime risk above which exactly the same proportion of participants lay.

For both comparisons (i.e. QRISK3 and QRISK-Lifetime 10-year prediction > 10%; and QRISK3 prediction > 10% and matched number of participants with the highest lifetime predicted risk), we examined the characteristics of patients recommended for treatment, the observed risk of CVD at 10 years, and the number needed to treat (NNT) to prevent one new CVD event assuming all people recommended for treatment actually took a statin and with a relative risk reduction of 25% for new CVD events. All models were fitted in R and Stata.

Reclassification using QRISK-Lifetime constrained to 10-year versus QRISK3 10-year prediction

In the first reclassification analysis examining people with 10-year predicted CVD risk > 10% using QRISK-Lifetime and QRISK3, QRISK-Lifetime classified fewer people as eligible to be offered a statin than QRISK3 (Table 4). QRISK-Lifetime classified 194,411 (15.4%) women as high risk, compared with 239,396 (19.0%) women classified as high risk by QRISK3. QRISK-Lifetime classified 276,369 (22.6%) men as high risk, compared with 341,962 (28.0%) men classified as high risk by QRISK3. For women, 15.1% were classified as high risk by both tools, 3.9% were classified as high risk by QRISK3 only and 0.3% were classified as high risk by QRISK-Lifetime only, over 10 years. For men, 21.9% were classified as high risk by both tools, 6.1% were classified as high risk by QRISK3 only and 0.7% were classified as high risk by QRISK-Lifetime only, over 10 years.

Based on 10-year risk prediction, the characteristics of people classified as high risk by each tool were similar (Table 5). Fewer people were recommended for treatment by QRISK-Lifetime and there were fewer observed events in people recommended for treatment by QRISK-Lifetime (women, 25,461 vs. 28,373; men, 33,450 vs. 37,026), but the percentage of people experiencing an event was higher (women, 13.2% vs. 11.9%; men, 12.1% vs. 10.8%). Among people recommended for treatment with a statin, the estimated NNT from statin prescription to prevent one event was 30 and 34 in women, and 33 and 37 in men, for QRISK-Lifetime and QRISK3, respectively (see Table 5).

Reclassification using QRISK-Lifetime lifetime risk versus QRISK3 10-year prediction

By design, the comparison with QRISK-Lifetime predicting to age 95 years (i.e. lifetime risk) was constrained to include 19.0% of women and 28.0% of men at the highest predicted lifetime risk (i.e. the same proportion of people identified as high risk based on QRISK3 10-year prediction) (see Table 4). Only 5.3% of women were identified as high risk by both QRISK3 and QRISK-Lifetime predicting to age 95 years, with a different 13.7% of women identified as high risk by one or other of the prediction tools. For men, 8.9% were identified as high risk by both prediction tools and a different 19.1% of men by one or other of the tools.

Compared with people identified as high risk by QRISK3 10-year prediction, people with the highest predicted lifetime risk were much younger, had a lower mean SBP and had a somewhat higher mean TC : HDL and BMI. In addition, a lower proportion of people with the highest predicted lifetime risk had treated hypertension and a much higher proportion had a family history of premature CVD and were from a minority ethnic background (see Table 5). Compared with people recommended for treatment based on 10-year predicted risk, there were fewer CVD events observed in people at the highest predicted lifetime risk [women, 9652 (4.0%) vs. 28,373 (11.9%); men, 14,725 (4.3%) vs. 37,026 (10.8%)]. For QRISK-Lifetime predicting to age 95 years, the estimated NNT to prevent one CVD event from statin treatment was 99 and 100 in women and men, respectively, compared with 34 and 37 in women and men, respectively, among those with > 10% 10-year QRISK3-predicted risk (see Table 5).

QRISK-Lifetime external validation

It is essentially impossible to validate lifetime risk models because lifetime follow-up of observed events either is not available or would be so historical that it would be a poor guide to performance of the model in a contemporary population. 60,85 Similar to the original internal validation,27 we evaluated QRISK-Lifetime over a 10-year prediction horizon. Within this constraint, QRISK-Lifetime had excellent discrimination in the whole population; however, as with QRISK3, discrimination was poor to moderate within the age strata, and moderate within the comorbidity strata. Calibration plots showed some underprediction in the whole population, with large underprediction in older people and people with multimorbidity at higher levels of predicted risk.

Comparing people recommended for treatment at a 10-year 10% risk threshold, QRISK-Lifetime (predicting over 10 years) recommended fewer people for statin treatment (i.e. 15.4% of women and 22.6% of men) than QRISK3 (i.e. 19.0% of women and 28.0% of men), although the people recommended experienced slightly more CVD events and the estimated NNT to prevent one CVD event was slightly lower for QRISK-Lifetime. Characteristics of people recommended for treatment over a 10-year prediction horizon were broadly similar.

Comparing people recommended for treatment by QRISK3-predicted 10-year risk ≥ 10% compared with the same proportion at highest estimated lifetime risk by QRISK-Lifetime, there was only a small overlap between the populations at highest predicted risk by the different tools. By design, each tool ‘recommended’ 19.0% of women and 28.0% of men for treatment. Only 5.3% of women and 8.9% of men were recommended for treatment by both tools, and the people recommended for treatment were considerably different, as other studies have found. 27,78 People with highest predicted lifetime risk were considerably younger, were more likely to have a family history of premature CVD and were more likely to be from a minority ethnic background. Over 10 years, people with highest predicted lifetime risk experienced many fewer CVD events (as expected given age differences), with people at highest predicted lifetime risk having an estimated NNT (with a statin to prevent one CVD event) approximately three times larger than for people recommended for treatment by QRISK3 (in women, QRISK-Lifetime highest risk NNT = 99 vs. QRISK3 risk > 10% NNT = 34; in men, QRISK-Lifetime highest risk NNT = 100 vs. QRISK3 risk > 10% NNT = 37). There is, therefore, a considerable leap of faith involved in treating based on lifetime risk, as the medium-term (10-year) benefit is considerably lower.

Limitations

This study has the same limitations as those already described in Chapter 2, Limitations, for any routine data study, but two particular limitations specifically apply. First, as with previous studies,35,67 we used multiple imputation, but the assumption that data are missing at random may be incorrect, and this is probably more likely to apply in younger people (i.e. the key target population for lifetime CVD risk prediction), as CVD risk assessment (particularly measurement of cholesterol) is likely to be carried out in people who are already suspected to be high risk. 35 Second, evaluating lifetime risk in a study with relatively short follow-up is intrinsically problematic. In this study, median follow-up of study participants was 5.7 years (interquartile range 2.2–10.2) years in women and 5.2 (interquartile range 2.0–9.3) years in men. The QRISK-Lifetime derivation paper does not state follow-up time,27 but follow-up time in this study is similar to follow-up in QRISK2 derivation, which uses a similar data set to QRISK-Lifetime. 86 Lifetime risk is, therefore, being estimated and evaluated from relatively short periods of observation. However, lifetime risk is estimated by assuming that future risk beyond the period of observation will be the same as that observed for older people during the period of observation. As QRISK-Lifetime systematically underpredicts CVD risk over the period of observation, then lifetime estimates must also underpredict. More generally, however, it is a very strong assumption that age-specific CVD incidence will be stable over the next few decades, given falling CVD incidence over the last few decades and large increases in obesity, sedentary behaviour and diabetes in the last two decades.

Data sources

We used the same data from CPRD GOLD as for the CVD risk modelling. 66 These data are similar to data in the QFracture-2012 derivation data set in terms of their inclusion of linked primary care and mortality data, but the CPRD GOLD data used also has linked hospital admission data, which we used for fracture ascertainment.

Identical to QRISK3 derivation and internal validation, patients were eligible for inclusion if they:

• were permanently registered with a general practice, contributing up-to-standard data for at least 1 year and with consent to link GP data to hospital discharge (HES admitted patient care) and mortality (ONS mortality registration) data

• were aged ≥ 30 years and < 100 years.

Cohort entry was defined as the latest date of an individual’s date of registration plus 1 year, the individual’s 30th birthday or 1 January 2004. Cohort exit was defined as the earliest of:

• the first MOF or first hip fracture (two distinct outcomes are modelled)

• deregistration from their participating general practice

• date of last data collection from their participating general practice

• end-of-study follow-up on 31 March 2016.

All outcomes and predictors are recorded blind to the study hypothesis because they are recorded as part of routine clinical care.

Sample size

The sample size is fixed by the size of the CPRD GOLD data set. No formal power calculation was, therefore, carried out, as the calculation could not alter study design and the available sample size was considered more than sufficient for the purpose. 54

Outcome definition

QFracture-2012 has separate models to predict two outcomes (i.e. MOF and hip fracture)26 and both models were validated. In this study, MOF was defined as hip, vertebral, wrist, proximal humeral or osteoporotic fractures leading to hospital admission ascertained from codes in the GP electronic health record (using Read codes), HES discharge diagnoses (ICD-10 codes) or ONS death registration (ICD-10 codes). QFracture-2012 does not publish lists of Read codes used to define these outcomes and, therefore, we used published codesets where available98 and otherwise derived our own (see Appendix 3, Tables 32 and 33).

Other variable definitions

We implemented the published QFracture-2012 risk model (under GNU Lesser General Public Licence version 3) and calculated QFracture-2012-predicted 10-year risk of a MOF and the risk of a hip fracture for all patients in our cohort. As with fracture outcomes, we derived codesets for each predictor (see Appendix 3, Tables 34–36), which were based on published codesets where available. 98

There were two differences in our implementation from QFracture-2012 derivation. First, in our data set, Townsend deprivation score was available as vigintiles (twentieths) only, rather than as raw values (to minimise disclosure), and so for this variable we fitted the median Townsend score of each vigintile (as with the CVD models described previously). Second, QFracture-2012 allowed values such as BMI recorded after the date of study entry but before any fracture outcome to be used in prediction, whereas in this analysis we restricted predictor values to those recorded before study entry to avoid using future information in prediction.

In addition to QFracture-2012, we calculated the CCI for each patient at baseline based on primary care Read codes. As with the CVD models, CCI was used only to stratify the analysis of discrimination and calibration by level of comorbidity (i.e. CCI score grouped into 0, 1, 2 and ≥ 3).

Missing data

The extent and management of missing data are detailed in Appendix 4, Table 37. As with QFracture-2012 derivation, patients with missing Townsend deprivation score were excluded from the cohort and patients with missing ethnicity were assumed to be white. For missing BMI, smoking status and alcohol status, multivariate imputation via chained equations68 was used to generate five imputed data sets. These data sets were then combined using Rubin’s rules69 to give summary point estimates, with confidence limits reflecting the uncertainty associated with imputing missing values. Reflecting that morbidity and prescribing recording in CPRD is generally good, and consistent with QFracture-2012 derivation, morbidities and prescribing variables used in the prediction score were assumed to be absent if not recorded. 66,99

Analytical methods: external validation

The performance of the QFracture-2012 risk score was assessed by examining discrimination and calibration in the same way as for the QRISK3 external validation (see Chapter 2, Methods). Discrimination was primarily evaluated using the truncated version of Harrell’s c-statistic to include only pairs where the earliest survival time is no later than 10 years after entry. A c-statistic of 0.5 indicates discrimination no better than chance, whereas a c-statistic of 1 indicates perfect discrimination. Evaluating how good discrimination is for values between 0.5 and 1 is arbitrary and involves judgement. We considered c-statistic values of 0.5–0.599 poor, values of 0.6–0.699 moderate, values of 0.7–0.799 good and values ≥ 0.8 excellent. In addition, we calculated the D-statistic of Royston and Sauerbrei (where higher values indicate better discrimination, and a difference of ≥ 0.1 is suggested as indicating meaningfully different discrimination)70 and a related R²-statistic appropriate for estimating explained variation for censored survival data. 71 Calibration was assessed for 10 equal-sized groups (deciles) of participants ranked by predicted risk, by plotting observed proportions versus predicted probabilities. Plots were generated separately by sex for all patients and for subgroups of age and CCI based on summary statistics pooled across the imputed data sets. When examining calibration, we estimated observed risk for censored data in two ways: (1) using the standard Kaplan–Meier estimator (which is consistent with the assumptions made in QFracture-2012 derivation in that it ignores competing risks) and (2) using the Aalen–Johansen estimator (an extension to allow for competing events, non-fracture death in this case). 84 All models were fitted in R and Stata.

Analytical methods: competing risk model derivation and internal validation

Competing risk model derivation and internal validation was carried out in the same data set as QFracture-2012 external validation. For this purpose, participants were randomly allocated to distinct derivation and test data sets in a 2 : 1 ratio, with allocation balanced in terms of age and final event status. The derivation data set was used to derive CFracture, a Fine–Gray model to predict the 10-year risk of experiencing a CVD event, accounting for the competing risk of non-CVD death. Separate models were estimated for men and women. Reflecting the overall aim of the project where we wished to explicitly compare prediction in models accounting for competing risk compared with ignoring competing risk, we included all the same main effects (predictors) and age interactions as QFracture-2012 modified, as follows. we accounted for non-fracture death as a second (competing) outcome using the Fine–Gray model and we additionally included the CCI score in the model (categorised as 0, 1, 2 or ≥ 3), as CCI is a well-validated predictor of total mortality. 12 Fine–Gray models allow the cumulative incidence function or probability of a CVD event occurring over time to be directly predicted, but the subdistribution HRs in the Fine–Gray models do not have a straightforward interpretation, as they describe the direction but not the magnitude of the effect of predictors on the cumulative incidence function. The use of fractional polynomials and the inclusion of complex interactions with age further complicate direct interpretation of model coefficients. Model coefficients are, therefore, not straightforwardly interpretable, but the derived model is provided in Appendix 4, Tables 38–41, to allow replication.

The performance of CFracture and QFracture-2012 was evaluated in the independent validation data set by examining discrimination and calibration, as described above. R 4.0.0 was used for all analyses.

Performance of CFracture for predicting major osteoporotic fracture

In the internal validation cohort, discrimination of CFracture for MOF in the whole population was excellent in women (Harrell’s c-statistic 0.813, 95% CI 0.810 to 0.816) and good in men (Harrell’s c-statistic 0.738, 95% CI 0.732 to 0.743), which is similar to QFracture-2012 for MOF in the same cohort (Table 8). Stratified by age, discrimination of both CFracture and QFracture-2012 was lower in all age groups and declined with age in women, although was similar in all age groups in men. Stratified by comorbidity, discrimination in the strata was somewhat worse than in the whole population, but there was no clear pattern of change with increasing comorbidity.

Calibration of CFracture for MOF was better than for QFracture-2012 in the whole population. In women, there was some underprediction at higher levels of predicted risk and a in men there was a similar calibration, but with some overprediction in the middle range of predicted risk (Figure 11). Stratified by age, calibration was good in women aged 30–64 years, with some underprediction in women aged 65–74 years and 75–84 years (Figure 12). In men, there was some overprediction a higher levels of predicted risk in those aged 30–64 years and 75–84 years, and some overprediction in those aged 65–74 years. In both men and women aged 85–99 years, calibration was poor, although considerably better than QFracture-2012 in this age group, as well as all others. Stratified by CCI, calibration in women was good, although with some underprediction in women with a CCI score of 0 and some overprediction in women with a CCI score ≥ 3 in the highest decile of predicted risk (Figure 13). Calibration in men stratified by CCI was more variable in those with a CCI score of 0 or 1 but was good in those with a CCI score of 2 and in those with a CCI score ≥ 3, for whom it was good apart from some overprediction in the decile of highest predicted risk.

FIGURE 11.

Whole-population calibration of CFracture and QFracture-2012 for MOF in the internal validation data set in the whole population. (a) Women; and (b) men. Observed risk is based on the Aalen–Johansen estimator, which accounts for competing mortality risk. Source: Livingstone et al. 51

FIGURE 12.

Calibration of CFracture and QFracture-2012 for MOF in the internal validation data set stratified by age group. (a) Women aged 30–64 years; (b) men aged 30–64 years; (c) women aged 65–74 years; (d) men aged 65–74 years; (e) women aged 75–84 years; (f) men aged 75–84 years; (g) women aged 85–99 years; and (h) men aged 85–99 years. Observed risk is based on the Aalen–Johansen estimator, which accounts for competing mortality risk. Ideal calibration lies on the reference line. Below the reference line is overprediction and above the reference line is underprediction. Source: Livingstone et al. 51

FIGURE 13.

Calibration of CFracture and QFracture-2012 for MOF in the internal validation data set stratified by mCCI. (a) Women, mCCI = 0; (b) men, mCCI = 0; (c) women mCCI = 1; (d) men, mCCI = 1; (e) women, mCCI = 2; (f) men, mCCI = 2; (g) women, mCCI ≥ 3; and (h) men, mCCI ≥ 3. Observed risk is based on the Aalen–Johansen estimator, which accounts for competing mortality risk. Ideal calibration lies on the reference line. Below the reference line is overprediction and above the reference line is underprediction. Source: Livingstone et al. 51

Performance of CFracture for predicting hip fracture

In the internal validation cohort, discrimination of CFracture for hip fracture in the whole population was excellent in both women (Harrell’s c-statistic 0.914, 95% CI 0.908 to 0.883) and men (Harrell’s c-statistic 0.886, 95% CI 0.877 to 0.895), similar to QFracture-2012 for hip fracture in the same cohort (Table 9). Stratified by age, discrimination of both CFracture and QFracture-2012 was lower in all age groups and declined with age in both men and women. Stratified by comorbidity, discrimination was good to excellent in all strata, although declined with increasing comorbidity.

Calibration of CFracture for hip fracture was better than for QFracture-2012 in the whole population. Calibration was good in women across all levels of predicted risk. Calibration was good in men, except in the highest decile of predicted risk where there was underprediction (Figure 14). Stratified by age, calibration was good in women aged 30–64 years and 75–84 years, with some overprediction in women aged 65–74 years (Figure 15). In men, calibration was reasonable in those aged 30–64 years, 65–74 years and 75–84 years, with some underprediction at the highest level of predicted risk. In both men and women aged 85–99 years, calibration was poor, with overprediction at most levels of predicted risk, although was considerably better than QFracture-2012 in this age group, as well as all others.

FIGURE 14.

Whole-population calibration of CFracture and QFracture-2012 for hip fracture in the internal validation data set in the whole population. (a) Women; and (b) men. Observed risk is based on the Aalen–Johansen estimator, which accounts for competing mortality risk. Source: Livingstone et al. 51

FIGURE 15.

Calibration of CFracture and QFracture-2012 for hip fracture in the internal validation data set stratified by age group. (a) Women aged 30–64 years; (b) men aged 30–64 years; (c) women aged 65–74 years; (d) men aged 65–74 years; (e) women aged 75–84 years; (f) men aged 75–84 years; (g) women aged 85–99 years; and (h) men aged 85–99 years. Observed risk is based on the Aalen–Johansen estimator, which accounts for competing mortality risk. Ideal calibration lies on the reference line. Below the reference line is overprediction and above the reference line is underprediction. Source: Livingstone et al. 51

Stratified by CCI, calibration in women was good, except in the highest decile of predicted risk where there was some overprediction in CCI scores 1, 2 and ≥ 3 (Figure 16). Calibration in men stratified by CCI was good, apart from underprediction in the highest decile of predicted risk for CCI scores 0 and 1, underprediction in the middle of the predicted risk range for CCI score 2 and some overprediction in the highest decile of predicted risk for CCI score ≥ 3.

FIGURE 16.

Calibration of CFracture and QFracture-2012 for hip fracture in the internal validation data set stratified by CCI. (a) Women, CCI = 0; (b) men, CCI = 0; (c) women, CCI = 1; (d) men, CCI = 1; (e) women, CCI = 2; (f) men, CCI = 2; (g) women, CCI ≥ 3; and (h) men, CCI ≥ 3. Observed risk is based on the Aalen–Johansen estimator, which accounts for competing mortality risk. Ideal calibration lies on the reference line. Below the reference line is overprediction and above the reference line is underprediction. Source: Livingstone et al. 51

QFracture-2012 external validation

QFracture-2012 was found to have very good to excellent discrimination in the total population, with higher discrimination observed in the hip fracture model than in the MOF model. However, discrimination was typically only poor to moderate in important subgroups, including older patients (as expected given that age is such a strong predictor of risk72,100) and moderate to good those with higher levels of multimorbidity. However, calibration was very poor, irrespective of how evaluated. When evaluated in its own terms (ignoring competing risk), QFracture-2012 showed consistent underprediction for both MOF and hip fracture. The most likely explanation for this underprediction is that fracture ascertainment in this study is more complete, as it includes fractures recorded during hospital admission in addition to fractures recorded in GP electronic health records and mortality registration data. For hip fracture, the inclusion of hospital data for fracture ascertainment likely explained much of the observed difference in hip fracture incidence, but less than half of observed differences in MOF incidence. Other reasons for higher observed incidence are the use of different index dates for first study entry (1 January 2004 in this study vs. 1 January 1998 in QFracture-2012), as recording of fractures in GP data is likely to have improved over time, as well as the use of different codesets to identify fracture (although QFracture-2012 does not have a published list of Read codes for predictors and so direct comparison is not possible).

When evaluated against observed fractures estimated accounting for competing risk, in general, underprediction reduced (i.e. failing to account for competing risk causes overprediction); however, there was very large overprediction at higher levels of predicted risk in older people and in people with more complex multimorbidity. Notably, in people aged 85–99 years and with CCI score ≥ 3, calibration was extremely poor, with observed risk flat or even declining across deciles of increasing predicted risk. Therefore, in summary, QFracture-2012 has two causes of poor calibration that operate in different directions. First, QFracture-2012 considerably underpredicts in all patients because derivation is based on incomplete ascertainment of fracture. Second, QFracture-2012 considerably overpredicts in people with high competing risk of death (primarily older people and in people with more complex multimorbidity).

CFracture derivation and internal validation

In the internal validation cohort, discrimination of CFracture for MOF and hip fracture in the whole population was good to excellent for MOF and excellent for hip fracture, and was similar to QFracture-2012 performance in the same data set. Stratified by age, CFracture discrimination was less good in all age groups, with worse discrimination with increasing age, except in men with MOF, for whom there was no obvious relationship with age. Stratified by CCI, discrimination for MOF was somewhat worse in all strata than in the whole population, but there was no clear pattern with age. Stratified by CCI, discrimination for hip fracture was similar to the whole population for a CCI score of 0, but discrimination declined somewhat with increasing comorbidity (although was good to excellent in all strata).

CFracture was better calibrated than QFracture-2012 in the whole population and in every strata, although was better calibrated in women than men and was better calibrated for hip fracture than for MOF. Strikingly, however, CFracture was poorly calibrated in women and men aged 85–99 years (although much less so than QFracture-2012).

These findings are somewhat different from previous external validations. 18,19 The first version of QFracture-201217 was independently, externally validated in The Health Improvement Network data set, which is a similar (and partially overlapping) set of practices to CPRD. This study found QFracture-2012 to have excellent discrimination and calibration in the whole population. 18 The updated QFracture-2012 (as evaluated in this study)26 was externally validated in CPRD by the QFracture-2012 derivation team, and this study, again, found excellent discrimination and calibration in the whole population. 19 This study differs from both previous external validations in two ways. First, we also identified fractures recorded during hospital admission, whereas the two previous studies identified fractures recorded in GP records and mortality registration data only. Better ascertainment of fractures would be expected to lead to underprediction of risk, as observed in this study. An Israeli external validation using both community and hospital data for ascertainment also observed considerable underprediction by QFracture-2012, consistent with our findings. 92 Second, this study examined calibration against observed outcomes estimated in the same way as previous external validations (using the Kaplan–Meier estimator, which ignores competing risk) and additionally accounting for competing risk (using the Aalen–Johansen estimator). As not accounting for competing risk leads to overprediction of fractures, this partially, but not completely, compensated for the observed underprediction in the whole population. However, as expected,96,97,100 accounting for competing risks led to large changes in observed risk in older people and in people with more multimorbidity, consistent with systematic overprediction by QFracture-2012 in groups with high competing mortality (despite systematic underprediction due to incomplete fracture ascertainment in QFracture-2012 derivation).

Limitations

This study has the same limitations as any routine data study (see Chapter 2, Limitations), but several particular limitations specifically apply. First, the observed incidence of fracture is higher in this study than in QFracture-2012 derivation and the two previous external validations. 18,19,26 There are three possible reasons for this. First, this study used hospital data to identify fractures (which accounts for most of the difference for hip fracture and some of the difference for MOF). Second, this study has an earliest study entry date of 1 January 2004 compared with the original QFracture-2012 external validation date of 27 June 1994 and the updated version external validation date of 1 January 1998, and recording of fracture outcomes in GP data are likely to have improved over time. Third, the codes used to identify fractures in GP data may be different; however, because none of the previous studies published their codesets we cannot explore this further. Within our own data, one issue is that humeral fractures were most commonly recorded in the GP data without specifying whether proximal or more distal. Therefore, we included non-site-specific humeral fractures as ‘proximal humerus’ fractures, which may lead to some misclassification (i.e. some false positives). However, registry data show that ≈ 80% of humeral fractures are proximal (and > 80% of non-proximal fractures are also low energy)101 and only including humeral fractures specified as proximal would lead to larger misclassification (i.e. a larger number of false negatives).

Second, although we explicitly accounted for censoring due to death in this study, the analysis still assumes that people who deregister from a CPRD practice have the same fracture risk as people who do not. This assumption is likely to be strong in older people for whom deregistration due to moving into extra-care housing or a care home is likely to be associated with higher falls and fracture risk.

Selection of exemplar medicines

We quantified the DTD of statins for the primary prevention of heart disease and bisphosphonates for the primary prevention of fractures using two medicine-taking case studies. Statins were chosen as an example of a class of medicines that are perceived by professionals to be benign, but which some people perceive as harmful. Bisphosphonates were selected because they are medicines we thought had an obvious influence on daily life (i.e. people who take a bisphosphonate are required to drink a large glass of water before taking the medication, remain upright for 30 minutes after taking it and avoid food and drink for 2 hours thereafter).

Selection of elicitation methods

We sought to get an accurate description of the medicine-taking health states within our survey so that our respondents could appropriately value them. In addition, we wanted to check that the methods we adopted were appropriate and robust enough to estimate small disutilities associated with the ongoing use of a medicine. Therefore, we reviewed methods from previous studies47,48,103 that had attempted to estimate DTD values and integrated those methods with the experience and views of our clinical research team, as well as our two patient representatives who were members of the research team. This review suggested that TTO should be the predominant method, but that there is emerging interest in using BWS to elicit DTD.

The context

The elicitation surveys were designed to take account of the specific context in which a respondent is taking a medicine for primary prevention. Two distinct surveys for each method (i.e. TTO and BWS) were designed to focus on each of the two exemplars of statins and bisphosphonates.

For each exemplar, respondents were asked to consider taking medicine A, which was a ‘one-off pill’ assumed to have no ongoing inconvenience, or medicine B, which was a daily pill for 10 years. In the TTO exercise, we also decided to include four scenarios to understand whether or not DTD values differed based on how the benefits and harms of the medication were framed. Therefore, we asked our respondents to consider medicines A and B in the context of the pills having (1) no side effects, (2) minor side effects (MSEs), (3) severe side effects (SSEs) and (4) reduced effectiveness. Finally, we also wanted to explore if there was any systematic difference between how different groups might value DTDs. In particular, we wanted to understand whether or not patients with experience of taking pills valued DTD differently from people with little or no experience. In addition, we thought it was important to explore other factors, such as age and sex.

Time trade-off exercise

The first method we selected to value the disutility of medicine-taking health states was the TTO method. TTOs are a widely used approach for eliciting utility values and were used to generate the valuations for the EuroQol-5 Dimensions, three-level version (EQ-5D-3L) health states, which are part of the current NICE reference case. In general, the TTO method involves asking respondents to consider the relative amounts of time (e.g. number of life-years) they would be willing to sacrifice to avoid a poorer health state. 110 In this study, the TTO method followed the approach taken by Hutchins et al. 48,103 and asked respondents the maximum amount of time they are willing to give up at the end of their life to avoid having to take a medicine. Per respondent, for each of the four questions, the estimated utility was calculated as the ratio x/t, where x is the final selected time period for the medicine A option (i.e. one pill taken once) and t is the full life-years assumed for the medicine B option (i.e. a pill taken every day for 10 years). 111

Best–worst scaling experiment

Best–worst scaling experiments are an extension of discrete choice experiments. 112 There are three types of BWS: case 1 (object case), case two (profile case) and case three (multiprofile case). 113 This study used a case 2 (profile case) BWS experiment. A profile case BWS experiments ask respondents to select their most preferred and least preferred items (defined by attributes and levels) in a question.

An argued advantage of profile case BWS over standard discrete choice experiments is that the choices made reveal more information about the relative strength of people’s preferences for each attribute in the design, using fewer questions, which could, in turn, reduce the response error. Importantly, using a BWS allows a rank ordering of the attributes in the experiment together with utility weights to be estimated. In this study, the BWS experiment was framed around the choice question ‘we want you to indicate which of the listed features of medicine A you think are the most and least likely to make you want to avoid taking the tablet’. The BWS experiment contained three attributes: (1) inconvenience (levels: no inconvenience and inconvenience), (2) probability of a MSE (levels: 1%, 5%, 9% and 13%) and (3) probability of a SSE (levels: 0.1%, 0.3%, 0.5% and 0.7%). Rapid reviews of the relevant published literature and input from the research team, including patient involvement, were used to generate a list of potential attributes and assigned levels. The BWS exercise was created using a full factorial design of 32 choice sets, which were split into four randomised blocks. Each respondent was assigned to one of the four blocks, comprising eight choice sets, with each choice set displaying three features of medicine A (i.e. inconvenience, MSE and SSE). The levels for these three features varied across choice sets.

Design of training materials

Consistent with emerging good practice in the design of a stated-preference study, training materials, that introduced the background for each case study (i.e. statin or bisphosphonate) and the attributes and levels used in the BWS exercise were created using a storyboard approach. 114 The same training materials were used for the TTO and BWS experiment. In addition, we used visual arrays to communicate absolute risk consistent with best practice. 115

Survey format and content

Online surveys were designed for each method (i.e. TTO and BWS), using the same approach and format for each exemplar medicine (i.e. statins and bisphosphonates). The surveys were formatted and administered online using Sawtooth software (Sawtooth Software, Inc., Provo, UT, USA). Respondents were sent a secure link to complete one of the surveys (no reminders were used). There were three parts to each survey. The first part of the survey consisted of the EQ-5D-3L questions, training materials and questions about the respondent’s attitude towards taking a medicine for the first time. The second part of the survey consisted of the main elicitation exercise (i.e. TTO or BWS). The third part of the survey included questions about whether respondents were currently taking or had ever taken one of the branded bisphosphonates (or statins), respondents’ perceived benefits and harms of bisphosphonates (or statins), whether or not respondents had experienced any side effects from bisphosphonates (or statins) and whether or not respondents found it inconvenient to take the medicine. Participants were also asked about their opinion on taking the medicine (i.e. whether they mind/dislike taking it), whether or not they took any other medicines, the number of medicines taken on a regular basis, the number of times medicines are taken daily and sociodemographic questions, such as age, gender, qualifications, employment status, ethnicity and religion. There were also non-compulsory probability questions to understand the respondents’ understanding of risk. Feedback questions towards the end of the survey included asking respondents about their confidence in making similar choices in real life, their perceived difficulty in making choices between alternatives and in understanding the survey, as well as general feedback comments to improve the clarity of the survey.

Piloting of experiments

Two pilot phases were conducted for the design of each experiment. For the TTO and BWS surveys, early piloting involved ‘think-aloud’ interviews with a sample of 19 patients recruited from a general practice in Greater Manchester. The intention of the early pilot was to understand whether or not the draft surveys, training materials and valuation exercises were sufficiently clear for respondents. We also wanted to understand how our respondents were interacting with the material presented. Following on from this, a few minor changes were made to the training materials and valuation exercise. The survey was then tested again in quantitative pilot studies involving members of the public to assess whether or not the data could be analysed from the survey design. No changes were made following the quantitative pilot study. The final elicitation studies were then launched and involved a sample of members of the public and a sample of people with experience of taking a statin or bisphosphonate.

Data samples

People with experience of taking a statin or bisphosphonate were recruited from general practices via the NHS Research Scotland Primary Care Network and the Scottish Health Research Register (SHARE which is a register of people living in Scotland, allowing recruitment after a search of their medical records). To be included, patients needed to have been prescribed a statin or bisphosphonate in the previous year, needed to be aged ≥ 30 years and should not have been diagnosed for dementia or be taking a drug for dementia. Members of the public for the valuation study were recruited online using the panel company Dynata (Shelton, CT, USA). Within this sample, we also identified members of the public with experience of taking a statin or bisphosphonate. Respondents to the online survey needed to be aged over ≥ 30 years, but otherwise the sample should be a demographically balanced representation of the general public willing to take an online survey.

Time trade-off

Analysis characteristics for respondents to the TTO for statins (n = 514) and respondents to the TTO for bisphosphonates (n = 365) are reported in Table 10. Statin patients tended to report marginally higher mean TTO values than public respondents [difference 0.007 (SE 0.003); i.e. the amount of life expectancy respondents were willing to sacrifice to avoid taking statins was 0.7 percentage points greater in the public respondents than in people with experience of statins], although this finding was not statistically significant. Bisphosphonate patients reported much higher TTO values – indicating less DTD – than public respondents [difference 0.024 (SE 0.006)], with this difference being statistically significant. For both statins and bisphosphonates, changing the question context did not alter the mean TTO scores by more than 0.01. Irrespective of the type of question or respondent, there was a clear difference between statins and bisphosphonates survey results (Figure 17). Mean TTO values for the entire statin sample (0.967) were higher than TTO scores for the bisphosphonate sample (0.933), meaning that the average respondent was willing to trade twice as much life expectancy to avoid bisphosphonates as they would to avoid statins [0.033 vs. 0.067; absolute difference 0.034 (SE 0.004)] (Table 11). Respondents who had experience of taking medications more than three times a day provided a much lower DTD value than respondents who had no experience of daily medicine use. None of the other explanatory variables had a statistically significant association with DTD size, including the number of pills taken per day.

FIGURE 17.

Kernel density plots showing distribution of TTO responses stratified by medicine, question context and respondent type. (a) Statin questions 1–4 utility values in patients; (b) statin questions 1–4 utility values in the public; (c) bisphosphonates questions 1–4 utility values in patients; and (d) bisphosphonates questions 1–4 utility values in the public.

Best–worst scaling experiment

Analysis characteristics for respondents completing the BWS for statins (n = 319) and respondents completing the BWS for bisphosphonates (n = 312) are reported in Table 12. Appendix 5, Table 46, shows the count data [normalised based on the number of levels for each attribute: in/convenience (two levels), MSEs (four levels), SSEs (four levels)] for the number of times a respondent chose the attribute level as ‘best’ and ‘worst’, and the difference between ‘best and worst’. Appendix 5, Figures 34 and 35, show the distribution of the count data for statins and bisphosphonates, respectively.

The data show that all the respondents had the strongest preference for taking a medicine with no inconvenience, which is the result expected a priori. This result was consistent between respondents with and without previous experience of taking a statin or bisphosphonate. The results from the pooled sample indicated that respondents had a strong dislike for experiencing the highest level of risk of a SSE. Overall, the respondents with experience of taking a statin or bisphosphonate indicated a greater dislike for a SSE and stronger preference for no inconvenience.

There was significant preference heterogeneity in the results. A fully correlated mixed logit model was, therefore, used to estimate preference weights for each attribute level relative to the reference level of 0.7% risk of a SSE (Table 13). The signs for all estimated coefficients in the fully correlated mixed logit model were consistent with a priori expectations based on the normalised best–worst scores. The estimated coefficients from the fully correlated mixed logit model were ‘re-scaled’ by setting no inconvenience at a value of 1 and 0.7% risk of SSEs at zero to calculate an indicative utility score for no inconvenience. Using this approach, the DTD for a statin was 0.2 and 0.46 for a statin and bisphosphonate, respectively (i.e. compared with a life free of the medications, life lived with medications should be seen as 80% or 54% as desirable, respectively).

Limitations

There are some limitations to our study that need to be considered. Although we endeavoured to communicate complex concepts clearly within the study, this had a set of clear trade-offs for participants. The length of the survey, the cognitive burden and the time required to complete the survey were noted as challenges. Some of our respondents reported that they had difficulty understanding the TTO methods, whereas others provided inconsistent values across the survey questions. In addition, participants who took part in our survey ultimately self-selected and this could be related to their ability to complete the survey, as well as their self-reported health. Applicability for a different patient or general population should be made based on a careful judgement of the self-reported characteristics summarised for this study cohort.

It is not immediately clear which values elicited for DTD should be used. In terms of face validity, the DTD elicited from the TTO seemed ‘better’. For the decision-analytic modelling reported in Chapters 6 and 7, we used the DTD values from the TTO study. This decision was made because the TTO method is consistent with how utility values are typically elicited to attach to health states in CEA. 117

Conclusion

This study indicated that DTD does exist. It is not clear whether or not the DTD elicited from the TTO or BWS should be used; however, using face validity would suggest that the TTO method produced more realistic DTD values. We suggest that DTD should be included in model-based CEA of long-term preventative medicines, and Chapters 6 and 7 explore the impact of using the elicited DTD values from the TTO method.

Selection of the statins case study

The decision-analytic model we adopted and developed for this decision problem takes a cohort-level state-transition (Markov) approach to predicting events of interest and estimating costs and effects. 10

Decision-analytic models of this type account for competing risks by design. For example, consider a Markov model with the three states of (1) well, (2) myocardial infarction (MI) and (3) dead (for our purposes, we will define the latter two as absorbing states from which further transitions are not of interest) and let per-cycle transition probabilities p_well→MI and p_well→dead both be 0.1. Assuming the whole cohort starts in the well state, after one cycle of the model 10% of people will be in the dead state and 10% will be in the MI state, leaving 80% in the well state. Note that, via the transition to death, 10% of the people who were initially at risk of a MI have been taken out of the at-risk state before cycle 2 without experiencing a MI. After the second cycle, we will have 18% of the cohort in the MI state (1 × 0.1 + 0.8 × 0.1); however, if we had a two-state well MI model with no competing risk of death, then it would have been 19% (1 × 0.1 + 0.9 × 0.1). If we run the model on for 10 cycles, and we find that 44.6% of people have had a MI, then without the competing transition to death the number would be 65.1%.

This simple example demonstrates that a competing risk structure is hardwired into the kind of state-transition models that are very commonly used to assess the cost-effectiveness of healthcare interventions. The example also shows that, to parameterise any given transition within such a structure, it is appropriate to rely on evidence that censors for competing events. Were we to base our inputs on time-to-event models that, themselves, adjust for competing risks (e.g. models explored in Chapters 2 and 3), our model would underestimate the probability of events of interest, as we would be effectively double counting the competing risks.

Clearly, however, it is important that the time-to-event evidence on which a decision model relies should be valid for the population to which the decision problem pertains. It should be clear from the above that if we have an inappropriate estimate of any one of our transitions, then it is not only with respect to that outcome that our model will be biased, as it will also compromise our estimate of any events with which it competes.

In an epidemiological risk prediction model, failing to account for competing risks will result in an overestimate of event rates for the whole cohort. In a state-transition decision model, the danger is that we will fail to account accurately for competing risks, which will lead to overprediction of events in people with higher-than-average competing risk and underprediction of events in people with lower-than-average competing risks.

General model updates

Except where otherwise described, our version of the model retains the inputs used for the CG181 analysis (see guideline documentation for full details10).

New model features specific to this project

Accounting for competing risk of non-cardiovascular death

The CG181 model10 accounts for background non-cardiovascular mortality using general population lifetables. A key motivation of this project (see Background) was to improve on this by ensuring that, for any modelled cohort, the competing risk of non-cardiovascular death reflects the life expectancy of people with the specified level of cardiovascular risk. Official lifetables stratify by age and sex; however, the lifetables do not tell us how all-cause mortality might be associated with cause-specific factors, such as cardiovascular risk. Therefore, we needed a method that will enable us, in effect, to estimate cardiovascular risk-specific non-cardiovascular lifetables.

To do this, we fitted a relative survival model to the cardiovascular CPRD data set used in the cardiovascular element of this project (see Chapter 2), using the multiplicative method of Andersen et al. 142 This approach is a variant of a Cox proportional hazards model that uses the life expectancy of a reference group (in this case, non-cardiovascular survival in the general population) as a time-varying covariate. For a formal definition of the model, see Appendix 6.

Our rationale for using a relative survival model is threefold. First, we needed some way of not only characterising observed survival in the CPRD data set (which has a median follow-up of 5 years), but also extrapolating to the lifetime horizon of the decision model. Decision-analytic economic models commonly use parametric models to accomplish this type of task; however, such models invariably characterise a discrete cohort of which we can plausibly assume that the members share a survival distribution (e.g. people with newly diagnosed cancer of a given type). In this case, we want to simulate heterogeneous groups of people (e.g. 40-year-old women with a 10-year cardiovascular risk of 2%; 80-year-old men with a 10-year cardiovascular risk of 40%) and it is implausible to assume that we can describe the expected survival of all such groups with a common parametric foundation. With a relative survival model, we do not need to define a functional form, we simply have to estimate a model that we can apply to a baseline drawn from empirical lifetable data, as this will naturally produce survival curves that vary in shape and scale, even though the relative difference between the general population and the modelled cohort is a constant function of the covariates. Second, we wanted an approach we could use to estimate the life expectancy of cohorts for whom the decision to initiate statins takes place in the present day. Characterising the relative relationship between observed survival and expected population survival from contemporaneous lifetables provides us with a model we can then apply to present-day lifetables. In doing this, we make the assumption that, although absolute life expectancy may have changed over time, the extent to which factors of interest affect it generalise across time with minimal bias. Third, relative survival models provide a parsimonious way of summarising complex data. If, as discussed above, we relied on parametric models to extrapolate to lifetime expectation, then we would inevitably need to fit multiple models to discrete subsets of the data; however, if we had chosen too few, then we would have inappropriately lumped people with heterogeneous expectations together, and if we had chosen too many, then we would have generated high-variance estimates with a tendency to overfit to artefacts in uncertain data. In contrast, the relative survival approach enables us to make use of all available data while adjusting for the things that might define groups with different expectations. Another reasonable approach would have been to fit a multistate model to the CPRD data; however, this would not solve the problem of extrapolation beyond the observed evidence and may have given us less flexibility in simulating cohorts with different baseline characteristics and risk profiles.

The output of the relative survival model is the multiplicative difference in time to event (in our case, non-cardiovascular death) between a person in a specific population and someone of the same age and sex in the general population. The coefficients show how this difference varies with each covariate.

For our models, the critical covariate is ΔQ, which we define as the difference between an individual’s predicted 10-year QRISK3 score and the average score for a person of the same age and sex. Because QRISK3 predicted risk is a probability it is bounded between 0 and 1 and, therefore, it becomes mathematically convenient to transform the estimates to a log-odds scale and so we work with Δlogit(Q). Because population lifetables do not provide information on QRISK3 scores we estimated average 10-year risks by age and sex using a regression on the CPRD data set (see Appendix 6, Table 54 and Figure 39). Owing to computational constraints, we were unable to fit a single relative survival model to the whole CPRD extract and, therefore, we developed separate models for men and women.

To begin with, we fitted relative survival models with a single coefficient [i.e. Δlogit(Q)].

As lifetables are already stratified by age, it is theoretically not necessary to include a term for it in a relative survival model. However, it is plausible that the effect of other covariates – in this case, Δlogit(Q) – varies with age. Therefore, we also fitted a more complicated model including age and its interaction with Δlogit(Q). At the same time, we explored introducing polynomial terms to account for non-linearity of effect. We found that using quadratic terms for both Δlogit(Q) and age improved the fit of the model (reducing AIC by several hundred points); however, the inclusion of higher-order polynomials provided diminishing returns. Therefore, the linear component of our more complex model takes the form:

⁽¹⁾ Δ l o g i t ( Q ) + Δ l o g i t ( Q ) 2 + a g e + a g e 2 + Δ l o g i t ( Q ) × a g e + Δ l o g i t ( Q ) 2 × a g e + Δ l o g i t ( Q ) × a g e 2 + Δ l o g i t ( Q ) 2 × a g e 2 .

Table 15 shows the results of simple and polynomial models fitted for men and women separately.

TABLE 15 - Relative survival models: non-cardiovascular mortality as a function of predicted cardiovascular risk and age in men and women

Parameter	Simple model, HR (95% CI)	Polynomial model, HR (95% CI)
Men
Δlogit(Q)	1.9314 (1.8747 to 1.9898)	2.7540 (2.5997 to 2.9174)
Δlogit(Q)2	–	0.8787 (0.8336 to 0.9262)
Age	–	1.0111 (1.0090 to 1.0133)
Age2	–	0.9998 (0.9998 to 0.9999)
Δlogit(Q) × age	–	1.0229 (1.0183 to 1.0275)
Δlogit(Q)2 × age	–	0.9858 (0.9833 to 0.9883)
Δlogit(Q) × age2	–	0.9989 (0.9987 to 0.9990)
Δlogit(Q)2 × age2	–	1.0001 (1.0000 to 1.0003)
AIC	518,564	517,989
Women
Δlogit(Q)	1.7061 (1.6789 to 1.7338)	2.2100 (2.0959 to 2.3304)
Δlogit(Q)2	–	0.9809 (0.9473 to 1.0158)
Age	–	0.9978 (0.9954 to 1.0001)
Age2	–	1.0002 (1.0001 to 1.0002)
Δlogit(Q) × age	–	1.0235 (1.0190 to 1.0281)
Δlogit(Q)2 × age	–	0.9918 (0.9899 to 0.9938)
Δlogit(Q) × age2	–	0.9994 (0.9992 to 0.9995)
Δlogit(Q)2 × age2	–	0.9999 (0.9998 to 1.0000)
AIC	536,994	536,427

To give a worked example, consider a 60-year-old man with a 10-year QRISK3 predicted CVD risk of 20%. First, using Appendix 6, Table 54, we calculate the average logit(QRISK3) for a person with those characteristics:

⁽²⁾ 12.90 + 0.353 × 60 − 5.60 × 10 − 3 × 60 2 + 6.31 × 10 − 5 × 60 3 − 4.18 × 10 − 7 × 60 4 + 1.42 × 0 − 9 × 60 5 − 11.02 + 0.986 × 60 + 0.0320 × 60 2 + 5.20 × 10 − 4 × 60 3 − 4.28 × 10 − 6 × 60 4 + 1.41 × 10 − 8 × 60 5 = − 1.815.

On a probability scale, this is a 10-year risk of 14% (i.e. our hypothetical man is subject to higher-than-average risk). Next, we calculate Δlogit(Q), that is the difference between observed and expected risk, and in this case it is logit(0.2) minus −1.815 = 0.429. Plugging this value into our polynomial model gives:

⁽³⁾ 1.013 × 0.429 − 0.129 × 0.429 2 + 0.01106 × 14.57 − 0.000176 × 14.57 2 + 0.02262 × 0.429 × 14.57 − 0.014348 × 0.429 2 × 14.57 − 0.001134 × 0.429 × 14.57 2 + 0.00014 × 0.429 2 × 14.57 2 = 0.540.

(Note that for summation we use the natural logarithm of the coefficients given in Table 15 and age is centred around the population mean of 45.43. Therefore, in this case, 60 – 45.43 = 14.57.)

When exponentiated, this value provides a HR of 1.72. Therefore, we estimate that a 60-year-old man with a 10-year QRISK3-predicted CVD risk of 20% is, alongside his raised cardiovascular risk, subject to a hazard of non-cardiovascular death that is more than two-thirds higher than that of an average person of that age and sex. This is the value that we apply to general population lifetables (adjusted to remove cardiovascular deaths) when we simulate the non-cardiovascular life expectancy of a cohort with those characteristics.

Figure 19 shows estimates of non-cardiovascular mortality generated in the way described above, compared with empirical data across different categories of cardiovascular risk and age. The empirical Kaplan–Meier curves represent observed non-cardiovascular mortality (censored for cardiovascular death) in the CPRD extract for people of the stated age, with baseline QRISK3 predictions in the specified brackets. For our modelled estimates, we start from ONS 3-year lifetables for England and Wales (for this illustration we use 2009–11 tables, as this is in the middle of the period covered by the CPRD data). We adjust these data to remove cardiovascular deaths, estimated using proportions recorded under ICD-10 codes I00–I99 in ONS’s ‘Deaths registered in England and Wales’ series. For each combination of risk and age bracket, we calculate HRs, as described above (for comparability, we fit at the mean QRISK3 prediction and age observed within that category in the CPRD data), and apply the HRs to the ONS curves, and this produces the adjusted curves shown in Figure 19.

FIGURE 19.

Fitted relative survival models for non-cardiovascular mortality compared with empirical data. (a) Men; and (b) women.

This exercise shows that the approach produces an excellent fit to the observed data. Even though the simple models (i.e. short dash) depend on a single variable, the simple models capture most of the observed variability in non-cardiovascular death. The polynomial models with age interactions (i.e. long dash) improve this fit somewhat further, at least falling within the 95% confidence limits of the observed survival functions in almost all cases. (One exception is 45- to 54-year-old women with a 10-year QRISK3 CVD prediction of 25% ± 2.5%, as these women appear to have better survival than women of the same age with lower cardiovascular risk; however, we consider this an implausible artefact of the relatively small sample of women in the category.)

Note that within each age bracket the unadjusted general population estimate (see Figure 19) remains the same across cardiovascular risk categories. Estimating competing-cause mortality on the basis of age and sex alone, without accounting for the way in which it correlates with the risk of interest, is the approach the CG181 model takes,10 but this analysis shows that such an approach can misestimate competing mortality risk to a potentially substantial degree.

For example, for women aged 65–74 years the unadjusted general population estimate suggests that women of this age have around an 84% chance of surviving 10 years without experiencing non-cardiovascular death. However, if a woman in this category had a predicted cardiovascular risk between 7.5% and 12.4%, then she would be uncommonly healthy and so her true 10-year non-cardiovascular survival would be around 91%. Conversely, if a woman had a high QRISK3-predicted risk of the order 27.5%–32.4%, then she would have below-average health and her 10-year non-cardiovascular survival probability would be around only 70%. In the former case, a model relying on unadjusted competing-cause mortality estimates would underestimate likely survival by a non-trivial amount. In the latter case, the bias is reversed. Note that in both cases the unadjusted population expectation is not within the confidence limits of the observed Kaplan–Meier estimate, whereas our adjusted estimates pass closely through them.

Model verification

Before making any updates or amendments, we verified that our version of the model exactly reproduced reported results from CG181. 10 We also performed technical verification to confirm that the model results changed in the expected direction when extreme values were used for specific model inputs.

Base-case analysis

Clinical Guideline 18110 presented cost–utility results for no preventative treatment and for statins at three levels of intensity (i.e. low, medium and high). All of the statin arms assumed class-level effects and used the acquisition costs of a single representative drug (i.e. simvastatin 10 mg/day, simvastatin 20 mg/day and atorvastatin 20 mg/day for low-, medium- and high-intensity therapy, respectively). An additional arm used identical effects to the high-intensity arm, but adopted the costs of atorvastatin 80 mg/day. Therefore, there were five arms in total.

Our model retains these five arms, and full incremental results are presented in Appendices 7 and 8. However, to simplify interpretation, the main analyses presented below focus on a pairwise comparison between high-intensity therapy with atorvastatin 20 mg/day and no treatment, as this is the approach the guideline ultimately recommended.

Sensitivity analyses

We use three types of sensitivity analyses (i.e. probabilistic sensitivity analysis based on 1000 iterations per output set, deterministic one-way sensitivity analysis and scenario analyses) to look at the impact of different assumptions about DTD.

General model updates

Before accounting for competing risks or DTD, in Table 16 we break down the independent and cumulative impact of each of the steps detailed in General model updates (above). The updates with the largest impact serve to improve the estimated cost-effectiveness of statins. Correcting the miscalculation by which statins were erroneously associated with increased incidence of some cardiovascular effects has the largest effect and, unsurprisingly, this increases the incremental QALYs conferred by statins and reduces their net costs. Updating costs relating to cardiovascular events also has a large effect. Because we now use higher estimates of cost across most categories (substantially so in some cases), the value of statins in preventing events increases. For both men and women, the net effect of the updates is to increase incremental QALYs and decrease incremental costs, although this does not lead to large changes in estimated cost-effectiveness.

Appendix 7 provides full results for all five simulated strategies at varying levels of cardiovascular risks, as well as a range of deterministic and probabilistic sensitivity analyses.

Accounting for competing risk

Before examining the impact that introducing adjustment for competing risks has on the cost-effectiveness of statins, it is useful to understand how it affects the modelled natural history of people at risk of CVD without accounting for any treatment effect. Appendix 8, Figure 47, provides model state occupancy graphs in two untreated cohorts with and without adjustment for competing risk of non-cardiovascular mortality. From these examples (see Appendix 8, Figure 47), we can see that the impact of adjusting for competing risk depends on a person’s cardiovascular risk compared with an average person of their sex and age:

• A 60-year-old woman with a 10-year QRISK3-predicted risk of 10% has above-average risk. In Figure 19, we can see that people in this category have slightly shorter non-cardiovascular life expectancy than people of the same age and sex in the general population. Consequently, when we include this adjustment in the cost-effectiveness model, the cohort experiences somewhat diminished life expectancy. The unadjusted lifetables tell us that a woman of this age has around a 35% chance of surviving until her 90th birthday. However, once we adjust for the raised risk of non-cardiovascular death that is associated with a 10-year cardiovascular risk of 10%, that figure drops to around 25%. A direct corollary of this is that the average woman with these characteristics will spend somewhat less time at risk of cardiovascular events, and we can see that the areas reflecting fatal and non-fatal cardiovascular history in the state occupancy graph (see Appendix 8, Figure 47) shrink between the unadjusted and adjusted versions.

• Conversely, a 60-year-old man with a 10-year QRISK3-predicted risk of 10% has below-average risk (note that in Figure 19 empirical and modelled data show a slightly longer non-cardiovascular life expectancy than would be expected for an average person of that age and sex). Now we have the opposite result, that is the man spends longer alive in the model than the unadjusted lifetable would suggest. The man’s chance of reaching 90 years of age rises from around 28% to around 38%. In addition, the man spends longer at risk of cardiovascular events and so the relevant areas grow instead of shrinking.

When we layer treatments effects, costs and QALYs on these dynamics, we derive cost-effectiveness results with similar characteristics (Table 17). When estimating the cost-effectiveness of statins, the impact of adjusting for competing risk depends on the cohort’s cardiovascular risk compared with an average person of their sex and age:

• For the 60-year-old women with a 10-year QRISK3-predicted risk of 10%, expected QALYs in both intervention and control arms decrease by around 0.65, reflecting the cohort’s below-average health. This means that a proportion of the cohort will die before they accrue all of the benefits of statins previously predicted for them. Consequently, health gains associated with high-intensity therapy (i.e. atorvastatin 20 mg/day), compared with no treatment, reduce by around 0.04 QALYs compared with the unadjusted model.

• Our 60-year-old men with a 10-year QRISK3-predicted risk of 10% were more healthy than their average contemporaries and so our cohort accrues around 0.6 QALYs per person more in both arms when adjusting for competing risk. This additional life expectancy potentially increases each man’s capacity to benefit from statins. Therefore, we see the incremental benefits of high-intensity therapy (i.e. atorvastatin 20 mg/day), compared with no treatment, rise by a little under 0.04 QALYs per person compared with the unadjusted model.

However, these adjustments affect treatment and non-treatment arms to a fairly similar degree. As a result, the cost-effectiveness of statins is not materially changed. The incremental cost-effectiveness ratio (ICER) for 60-year-old men with a 10-year QRISK3-predicted risk of 10% goes down by around £500 per QALY. For women with the same age and risk, the ICER goes up by a similar amount. In all cases, ICERs remain far below conventional thresholds representing an effective use of NHS resources.

Appendix 8, Figure 54, illustrates the same comparisons when analysed probabilistically. The same features are evident as in the deterministic results. For 60-year-olds with a 10% 10-year cardiovascular risk, adjusting for competing risk of non-cardiovascular death increases QALYs for men and decreases QALYs for women, although resulting cost-effectiveness estimates are not materially altered.

In Appendix 8, Figure 48, we depict the cost-effectiveness of statins for people of different ages and cardiovascular risks, and how adjusting for competing risk of non-cardiovascular death affects these results. Even before adopting this adjustment, the model suggests that statins represent a good use of resources for almost everyone. It is only for people aged ≥ 60 years with the lowest cardiovascular risk that statins represent poor value for money. However, adjusting for competing risk of non-cardiovascular death removes even this small subgroup. In practice, the distinction is moot if QRISK3 is used to predict cardiovascular risk, as it is essentially impossible for people in those age brackets to have 10-year risks low enough to enter the cost-ineffective zone. If such people did exist, then they would have extraordinary life expectancy, which is why the adjusted model concludes that it would still be good value to offer them statins, as there is every chance that even the oldest people would live to realise their benefit.

Appendix 8 provides full results from the updated model, adjusted for competing risk for all five simulated strategies at varying levels of cardiovascular risks, as well as a range of deterministic and probabilistic sensitivity analyses.

Incorporating direct treatment disutility

Introducing DTD to the analysis has a substantial impact on the results. Table 18 shows cost-effectiveness results for 60-year-old men and women under our four DTD scenarios. For the permanent scenario (with a utility multiplier of 0.966), DTD throughout the period people are taking statins for primary prevention is enough to render the statins net harmful. Although the intervention provides around 0.3 QALYs, the small day-to-day disutility of taking the tablets amounts to more than 0.4 discounted QALYs over a lifetime, and so the tablets end up doing more harm than good. Assuming that DTD lasts for no more than 10 years (with either constant or diminishing detriment over that period) attenuates but does not eradicate the benefits of statins in these cohorts, and so the statins remain reasonable value for money, assuming that QALYs are valued at conventional levels.

We see the effect of differing DTD assumptions on the cost-effectiveness of statins across a range of ages and baseline cardiovascular event risks in Appendix 8, Figure 50.

Even though a lifetime’s unremitting DTD represents a substantial decrement to expected quality of life, statins remain an effective use of resources for younger people with higher levels of risk. This is not only because people with high levels of risk stand to gain more from taking preventative statins (as they are more likely to have an event to prevent), but also because they spend less long in the primary prevention state (as they move to the event states more swiftly and they are also, now we adjust for competing risk of non-cardiovascular mortality, more likely to die of other causes). In contrast, people with lower levels of risk may need to take statins for a long time before realising benefits, and this becomes a poor trade-off when DTD is permanent. Older people, too, experience net harm from statins if they cannot get used to taking them, and this is because the lower life expectancy of older people gives them limited time in which to experience potential benefits of the intervention, while its harms are unavoidable.

For time-limited (10-year) DTD, the age–risk threshold at which statins become cost-effective is very close to the mean QRISK3-predicted risk observed in the population. Under this scenario, anyone with above-average risk stands to gain from treatment, whereas the model suggests that anyone with below-average risk will experience net harm.

If we can assume that DTD declines over time, then statins remain good value for everyone except people with a theoretically possible – although practically unlikely – profile reflecting extraordinarily good health for their age.

We show probabilistic versions of these calculations for a representative range of age–risk profiles in Appendix 8, Figure 58. In Appendix 8, Figure 58, the ‘clouds’ converge as age and risk rises, showing that DTD has less of an absolute effect on model outputs (owing to shorter time on preventative treatment with lesser life expectancy and higher event rates). However, net value for money remains relatively unaffected. In the cost-effectiveness acceptability curves, statins are always associated with less than 20% chance of cost-effectiveness when we assume permanent DTD (assuming we value QALYs at ≥ £20,000). Conversely, without DTD, statins are certain to be cost-effective as long as we value QALYs at ≥ £6000. The temporary DTD scenarios also result in a high probability of cost-effectiveness, with the exception of 10-year DTD in 80-year-olds with a 30% risk of cardiovascular event. Under this scenario, there is around a 50 : 50 chance that statins are worth paying for, and this is fairly invariant to the value we place on QALYs.

The interaction between life expectancy, level of risk and duration of DTD discussed above begins to suggest another way of conceptualising benefit–harm trade-offs. We have previously argued that the ‘pay-off time’, that is the minimum time until people can expect net benefit from a course of action with up-front harms and long-term benefits, may provide a useful heuristic for thinking about these issues in shared decision-making. 41 Appendix 8, Figure 59, shows cumulative incremental QALYs over time for four example profiles across our four DTD scenarios (we do not discount outcomes in this instance, as we take the view that it would be impossible for a patient to arrive at their own conception of the interaction between benefit and time while mentally adjusting for the fact that an interaction of that type is hardwired into the calculation).

Under the permanent DTD scenario, for a 50-year-old with a 10-year CVD risk of 5%, a decision to take statins would confer net QALY gains only after 48 years of treatment. For 60-year-olds at 10% risk, 70-year-olds at 20% risk and 80-year-olds at 30% risk, the expected pay-off times are 36 years, 25 years and 22 years, respectively. When we assume temporary DTD of one form or another, the expected pay-off times are typically in the range of 10–20 years.

Interaction of competing risk and direct treatment disutility

Bringing the analyses above together, we can look at the combined effect of competing risk of non-cardiovascular death and DTD on the effectiveness and cost-effectiveness of statins for the primary prevention of CVD.

To do this, we present a 4 × 4 cross-categorisation of scenarios, within each of which we derive results for people at different ages and different levels of baseline risk. For competing risk calculations, we present two additional scenarios in which we halve and double the cohort’s hazard of non-cardiovascular death. These scenarios are representative of a decision problem in which people are atypically healthy or atypically well to a degree over and above what we would expect by adjusting for their cardiovascular risk alone.

Figure 20a visualises the expected incremental QALY gains associated with high-intensity statins compared with none in each of these scenarios. Without DTD, statins are associated with very little disutility in the model and, therefore, statins generate some degree of QALY benefit for people at any age and any level of cardiovascular and competing risk (see first column Figure 20a). However, when we introduce DTD, the threshold at which treatment confers net QALY gains begins to depend on age. For someone who is persistently averse to pill-taking (implying lifelong DTD; see last column of Figure 20a), many combinations of age and risk result in statins doing more harm than good (see black area in Figure 20a). For example, people with such preferences in their mid-70s and older would need a 10-year QRISK3-predicted risk of over 30% before statins would be net beneficial. If people had additional long-term conditions leading to a doubled risk of non-cardiovascular death, then that threshold would rise to 40%. Even if we assume taking pills is something people get used to (i.e. time-limited DTD), it is easy to find combinations of age and risk where statins’ DTD outweighs their cardiovascular benefits.

FIGURE 20.

Clinical effectiveness and cost-effectiveness of high-intensity statins, as a function of age and cardiovascular risk, with different levels of DTD and competing risk of non-cardiovascular death. (a) Effectiveness (incremental QALYs); and (b) cost-effectiveness (when QALYs are valued at £20,000–30,000).

Figure 20b provides a similar analysis but introduces expected costs to the equation to establish the circumstances in which statins would be judged to provide an effective use of NHS resources (when QALYs are valued at NICE-recommended levels of £20,000–30,000). The combinations of characteristics for which statins do not provide good value for money are similar to characteristics for which they generate QALY loss (except under unusual circumstances, a technology will not be cost-effective if it is not associated with QALY gains). However, there are a small number of combinations where, although the model predicts positive QALY gains for statins, the positive QALY gains are so small that the cost of the drugs outweighs the gains (e.g. under the permanent DTD scenario, 50-year-olds with a 10% risk of cardiovascular event and an adjusted but unaltered hazard of non-cardiovascular death, as such people expect a tiny lifetime benefit of 0.006 QALYs, but at a cost that equates to an ICER of around £44,000 per QALY).

As a scenario analysis, we also explored applying DTD as an absolute decrement rather than a relative multiplier, and this results in a slightly greater impact for DTD. For example, for 80-year-olds under time-limited and permanent DTD scenarios, statins do more harm than good at all levels of risk we analyse (see Appendix 8, Figure 60).

Main findings

When it comes to state-transition decision models, the question is not whether or not to account for competing risks, but how to do so. Our analysis shows that, where competing risks are handled in a naive manner (e.g. assuming population lifetables for all-cause death), analyses overestimate the clinical effectiveness and cost-effectiveness of preventative interventions in people with above-average risk and underestimate the clinical effectiveness and cost-effectiveness in people with below-average risk. In the case of statins for the primary prevention of CVD, the intervention appears to be almost universally cost-effective in a way that renders these biases of limited impact. However, we have shown that the potential for meaningful differences exists. Indeed, one only has to imagine high-intensity statins cost £18 per pack instead of £1.10, and appropriate adjustment for competing risk of non-cardiovascular death would make the difference between the intervention meeting and not meeting conventional thresholds of cost-effectiveness for some people, for example 75-year-old men with a 10-year cardiovascular risk of 15% (unadjusted ICER compared with no treatment, £25,700/QALY; adjusted ICER, £18,100/QALY).

Incorporating DTD has a more immediately obvious effect on the cost-effectiveness outputs of the model. However, although we see no reason for decision models not to include the most robust estimates of competing risks possible in their base cases, we would hesitate to recommend that modellers should include population-average DTD in their base-case cost–utility results, as this would result in statins being recommended for only people who appear in the dark blue areas in the bottom-right of Figure 20b for whichever DTD scenario decision-makers prefer. In Chapter 5, we argued that, rather than assuming blanket disutility, we should provide decision-makers with evidence of circumstances under which, and people for whom, DTD might tip the balance of benefits and harms (if any). In turn, decision-makers should use this information to provide guidance to prescribers that can inform shared decision-making. This case study illustrates the point well.

To return to our example of a 75-year-old man with a 10-year cardiovascular risk of 15%, without accounting for DTD, we would expect the man to achieve net QALY gains from statins (i.e. an expected value of 0.19 QALYs, which is equivalent to over 2 months in perfect health). However, once we introduce DTD, that expectation reduces. Under the diminishing DTD scenario, the man can expect fewer than 0.06 QALYs, whereas the time-limited and permanent DTD scenarios predict net QALY loss.

Nevertheless, even if we think the permanent DTD scenario is the best way to capture population-level average disutility, we would not argue in favour of adjusting everyone’s expectation to account for it (which would have the effect of rendering statins poor value for money in people with these characteristics). This is because – as clearly shown in Chapter 5 – some people anticipate negligible disutility from taking statins and, indeed, some people may even anticipate benefit over and above the directly clinical effects (i.e. ‘peace of mind’). To deny such people access to an effective treatment would clearly not maximise societal welfare, even if the average person would not feel the same. Therefore, we suggest that the appropriate way to make use of this information is as an objective basis on which decision-makers can identify what NICE refers to as ‘preference-sensitive decision-points’. 125 In this paradigm, guidelines should help prescribers to summarise benefits and harms of possible options for people for whom they are indicated, and this should include the process characteristics underlying DTD.

Precisely how such information makes it into shared decision-making is an unanswered question. It is theoretically possible to envisage ex ante preference elicitation, along similar lines as our experiment, which could be integrated into a personalised QALY-maximising analysis to provide recommendations for individual patients, and, unquestionably, this would be practically onerous, but, arguably, it would be unnecessarily prescriptive. However, we think it is useful to augment prescribers’ understanding of the strength of preference for avoiding a technology that would be necessary to outweigh its expected benefits, and our findings give an example of how this might be possible.

Since CG181, there have been two CEAs143,144 of statins for the primary prevention of CVD that predicted cardiovascular events and non-cardiovascular death based on a common set of risk factors. To do this, the analyses used either separate models144 or cause-specific survival models. 143 Like our relative survival model, these approaches should accurately predict time spent living with CVD and, therefore, QALYs because they capture the crucial correlation between cardiovascular risk and competing risk of non-cardiovascular death. However, neither analysis tested whether or not its results were sensitive to this competing risk adjustment, that is neither analysis modelled scenarios where risk of non-cardiovascular death was assumed to be independent of cardiovascular risk factors for comparison. In addition, one of these analyses found that cost-effectiveness was sensitive to DTD in the form of an additive decrement ascribed to ‘pill burden’, which was applied to each QALY. 144 Another primary prevention model124 has also reported that the cost-effectiveness of statins is sensitive to an arbitrary disutility for DTD.

There is recent, compelling evidence that many self-reported muscle symptoms in people taking statins can be explained as ‘nocebo’ effects. 43 Although this finding may represent a fascinating psychological insight, it is of limited practical value to prescribers, as few people experiencing pain when taking a statin will easily accept that their symptoms are ‘all in the mind’. However, we argue that DTD provides a helpful way to conceptualise and quantify the authentic and predictable harm that such people experience. Authors usually cite process and inconvenience factors to explain why people ascribe disutility to taking a preventative medicine they have been assured is benign. 40 However, it is also likely that part of the internal calculus reflects (1) a degree of mistrust that the substance truly has no adverse effects and (2) a reluctance to become ‘a patient’ who – among other inconveniences – will henceforth be on-guard for unpleasant symptoms. 145 Although few people would describe it in such terms, no one wants to volunteer for a nocebo effect. Therefore, by quantifying the things that people ex ante want to avoid from preventative treatment, we believe that we capture some or all of the subtler harms that people taking them report ex post. As explored in this chapter, this enables us to balance these harms against the expected benefits and costs of treatment to estimate whether or not we expect net benefit overall.

Limitations

First, we have not been able, as originally intended, to take account of the extent to which specific co-existing long-term conditions might complicate the issues explored in this chapter. In particular, we wanted to look at type 1 diabetes and CKD, as NICE identified these conditions as areas of interest in the surveillance proposal in which it announced its intention to update CG181. 13 However, we found that including additional binary covariates for one or both of these factors did not improve the fit of our relative survival models. Therefore, this implies that the impact of type 1 diabetes and CKD on cardiovascular risk is well captured by the QRISK3 algorithm, and the factors do not differentially confound the relationship between cardiovascular risk and non-cardiovascular mortality. We also found that type 1 diabetes and CKD had no meaningful effect on the type of first cardiovascular event people experience. Having found no clear point of differentiation, it would not have been meaningful to stratify our analyses. Of course, it would also be possible to tailor other model inputs (e.g. baseline utility, event disutility, cost of events) to represent people with particular long-term conditions, and this might affect outputs in non-trivial ways.

Second, our relative survival models use predicted cardiovascular risk (i.e. QRISK3 10-year predicted risk) as an overarching indicator of non-cardiovascular life expectancy. We believe that we achieved a very good fit to observed data in this way (see Figure 19). Nevertheless, it would be possible to construct a more sophisticated model that, instead of using a summary risk prediction measure, estimates relative survival as a function of each of the individual covariates on which the prediction itself relies. For example, an 80-year-old non-smoking man with type 2 diabetes and SBP of 140 mmHg has an identical 10-year QRISK3-predicted CVD risk to an 80-year-old moderate-smoking man with no diabetes and SBP of 166 mmHg. If the factors that distinguish these two individuals affect their non-cardiovascular life expectancy differentially, then our relative survival model will fail to capture this. By accounting for these factors, it might be possible to achieve an even more accurate prediction of life expectancy. However, as noted above, we found that including terms for type 1 diabetes and/or CKD did not improve the fit of our relative survival models. Nevertheless, we cannot rule out the possibility that other variables would explain some residual heterogeneity between observed and modelled survival expectation.

Third, because secondary prevention of CVD is beyond the scope of our project, we did not review any evidence about the natural and treated history of people who have had a cardiovascular event. For all secondary transitions, the model relies on the same inputs that the developers identified for CG181. 10 Although this part of the pathway is outside our decision problem, it has an important impact on our estimates of lifetime costs and effects, as it defines some of the value that an intervention provides by preventing first events. We suggest that any future updates to the model should review evidence in this area, even if secondary prevention is not the analysts’ focus.

Fourth, we do not make any adjustments to the model’s effectiveness inputs to account for imperfect patient adherence to statins, and this means we effectively assume adherence is identical to that observed in the intention-to-treat trials that underpin the effect estimates. If adherence in practice is worse than in trials, then the model is likely to overestimate, to some degree, the clinical effectiveness and cost-effectiveness of statins.

Finally, the model compares all-or-nothing strategies. In reality, the counterfactual to offering statins to a given cohort is not never offering them statins, it is not offering statins for the time being, with the option of revisiting the decision later. To simulate this decision problem, we would require detailed longitudinal data with which to project the development of risk factors over time. It would also probably necessitate moving away from a cohort model to a stochastically evaluated individual patient simulation. However, we would argue that this analysis highlights some of the major advantages of analytically evaluated cohort models. The kind of many-way scenario analyses we have been able to generate would be computationally intractable and muddied by Monte Carlo error in a patient-level simulation, and we feel sure that the kind of analytical flexibility from which our analysis benefits outweighs the flexibility in representing the pathway that would come with a patient-level simulation.

Conclusion

This analysis has shown that updating the CG181 model with newer model inputs did not materially change the conclusion that statins are an effective use of healthcare resources. In the case of statins for the primary prevention of CVD, the intervention appeared to be almost universally cost-effective in a way that renders the impact of including competing risk inconsequential. However, we argue that failure to account for competing risk in CEAs of primary prevention strategies will produce biased results and the impact of this omission is likely to become more apparent for interventions with relatively higher costs compared with their health benefits. Incorporating DTD had a more immediately obvious effect on the cost-effectiveness of statins. We advise that the impact of including DTD in sensitivity analyses on the outputs of a decision-analytic model should be an integral component of CEAs of primary prevention strategies.

Selection of the case study

The decision to offer oral bisphosphonates for the primary prevention of osteoporotic fragility fracture is a second potentially illuminating case study for examining the impact of competing risks and DTD on CEAs. The risk factors for fragility fracture include many behavioural (e.g. smoking) and pathological (e.g. type 1 diabetes) factors that are also associated with shorter life expectancy. Failing to account for these factors could lead to inflated estimates of the value of medicines that reduce the risk of fracture. Furthermore, although there is little doubt that oral bisphosphonates are effective in reducing incidence of fragility fracture, bisphosphonates are inconvenient to take, and this raises the possibility that DTD may attenuate or even outweigh the gains that previous CEAs have estimated.

NICE initially recommended two oral bisphosphonates (i.e. alendronic acid and risedronate sodium) in in its Technology Appraisal 160 in 2008. 146 A subsequent appraisal (i.e. TA46412) confirmed this decision and added a third oral medicine in the class (i.e. ibandronic acid) to the recommended options. The evidence available to the decision-making committee for TA46412 included a CEA, synthesising evidence on the benefits, harms and costs of preventative treatments, including oral bisphosphonates. The analysis concluded that oral bisphosphonates are likely to represent an effective use of NHS resources for anyone whose 10-year risk of MOF exceeds 1.5%. 147

The same authors updated their model as part of a further health technology assessment (HTA) for a subsequent NICE appraisal, and this focused on antiosteoporosis medicines other than bisphosphonates, but included bisphosphonates as comparators, although NICE subsequently chose to suspend the appraisal. 148 The updated model found that a somewhat higher level of risk, than in TA464,12 would be necessary for oral bisphosphonates to provide reasonable value for money, although the precise threshold is not quantified. 149

The model underpinning these analyses was a discrete-event simulation. Unlike the cohort model explored in Chapter 6, discrete-event simulations work by generating a virtual population of people and simulating their lives one by one. Discrete-event simulations handle competing events by randomly generating next-event times for all outcomes of interest, processing the event that will occur first and then, depending on event type and model structure, moving on to the next event or recalculating some or all next-event times to reflect the simulated person’s updated history. The simulation may also terminate (e.g. if the event is death).

Given adequate data and implementation, discrete-event simulation is an appropriate way of accounting for competing events and, because the model handles them simultaneously, it is appropriate to parameterise each event distribution using time-to-event methods that censor for competing events. This is because, in a patient-level simulation as in clinical reality, events compete naturally, and the occurrence of some events will preclude the occurrence of others. In the case in hand, some simulated people will experience hip fracture before death; however, some people will experience death first, thereby precluding a fracture event (although the model ‘knows’ when that person was destined to experience a fracture had death not intervened).

Therefore, in contrast to the situation where we want to estimate the probability that a single event in isolation will occur (as in the epidemiological analyses in Chapters 2 and 3), modellers developing simulations with multiple events need inputs for each that censor for intervening events. In the case in hand, we actively want our estimates of time to hip fracture to reflect a world in which death is not possible (e.g. relying on a Kaplan–Meier estimator censoring for death), as the model will simultaneously and independently apply an estimate of time to death. So long as both inputs represent valid estimates for the modelled population, the incidence of hip fracture observed in the model will replicate the Aalen–Johansen estimator in a real population.

However, as with a cohort model, if the estimate for any one of our competing events is biased, then we will have biased results for all outcomes. Again, a commonly overlooked question is whether or not unadjusted general population lifetables represent an appropriate estimate of life expectancy for the people the model simulates. In a patient-level model, we have information about some characteristics for each person simulated. It is common to use the age and sex of the simuland to inform the distribution from which the model draws the time to all-cause death. However, other characteristics may also be associated with greater or lesser life expectancy. In the case in hand, we are particularly interested in predicted fracture risk. If people with high fracture risk also have an increased risk of other-cause death and we fail to account for this, then we will overestimate the number of fractures that occur and this, in turn, will overstate the clinical effectiveness and cost-effectiveness of any intervention that has the potential to reduce fracture incidence.

Previous models have estimated that the benefits of bisphosphonates for the average person are, in absolute terms, small. TA46412 found that, when compared with no treatment, alendronate is associated with an average health gain no greater than 0.00247 QALYs, even in the highest-risk tenth of people. 12 The 2020 HTA149 introduced a little more benefit for oral bisphosphonates, but alendronate still conferred only mean gains in the range of 0.0001–0.0058 QALYs (equivalent to up to 2 quality-adjusted days of life in perfect health), depending on baseline risk.

These small benefits should not be casually dismissed, given the low cost of the intervention and the large population who could expect, on average, to gain. Nevertheless, it is clear that it would not take much to offset the expected gains. Both TA46412 and the 2020 HTA149 incorporate an estimate of QALY loss owing to adverse gastrointestinal effects for a proportion of people taking oral alendronate. However, these analyses do not reflect the fact that, even when well tolerated, oral bisphosphonates are inconvenient to take. As we have seen (see Chapter 5), even when taken weekly, oral bisphosphonates are associated with non-trivial DTD (significantly greater than for daily statins). Therefore, it is plausible that the routine downsides of taking bisphosphonates might attenuate or outweigh their benefits in preventing fractures, although this trade-off has not previously been explored.

Our aim was to explore the impact that better accounting for competing risk of death and incorporation of DTD might have on the estimated cost-effectiveness of bisphosphonates. TA46412 stipulates that prescribers should offer bisphosphonate treatment ‘when [it] is appropriate, taking into account [the person’s] risk of fracture, their risk of adverse effects from bisphosphonates, and their clinical circumstances and preferences’ (© NICE 2017 Bisphosphonates for Treating Osteoporosis. Available from www.nice.org.uk/guidance/ta464. All rights reserved. Subject to Notice of rights). Although fracture risk and adverse events are central to previous analyses of the cost-effectiveness of bisphosphonates, there has been little formal consideration of how decision-making might be affected by people’s clinical circumstances and preferences. The methods introduced in this project enable us to explore these dimensions and weigh them against the benefits and harms that existing models estimate.

Model overview

The model is a patient-level discrete-event simulation. Figure 21 illustrates the structure of the model and Table 19 summarises key design criteria. Where state-transition models define a healthcare pathway in terms of states, discrete-event models proceed by way of events. The events that are possible in the TA464 model12 are fractures (divided into four categories of hip, vertebral, proximal humerus and wrist), admission to full-time care and death (fracture related and other cause). Simulated people who experience fractures incur costs and disutility, and their risk of additional negative events (e.g. further fractures, admission to full-time care, death) increases.

FIGURE 21.

Structure of discrete-event model: events that can occur in a simulated person’s lifetime.

Our study focused exclusively on the pairwise comparison between oral alendronate and no treatment. Our DTD elicitation exercise did not distinguish between different oral bisphosphonates, and alendronate is, by some distance, the most commonly prescribed chemical in the class, consistently accounting for over 80% of prescriptions of agents classified in British National Formulary section 6.6 (drugs affecting bone metabolism). 151

General model updates

We performed fewer general updates to the bisphosphonates model than we did for the statins analysis. Producing a revised best estimate of the cost-effectiveness of bisphosphonates is beyond the scope of this project and analysis, therefore, focuses on exploring whether or not accounting for competing risks and DTD would have materially affected the historical analyses. One exception is that we have updated how the model implements persistence with oral therapy (see below), as this not only influences the effectiveness of bisphosphonates but also has potentially important implications for DTD, which is one of our central concerns.

The 2020 bisphosphonates HTA149 introduced some other modifications to the TA464 model12 that we have adopted as detailed below.

New model features specific to this project

Accounting for competing risk of non-fracture death

To estimate non-fracture death in people at risk of osteoporotic fragility fracture, we used identical methods to those we adopted for non-cardiovascular death in people at risk of CVD (see Chapter 6, Accounting for competing risk of non-cardiovascular death).

First, we fitted a model to estimate average population fracture risk as a function of age and sex (see Appendix 9, Table 57 and Figure 61). We noted that, in the CPRD data set, a relatively consistent exponential increase in risk throughout the first nine decades of life tails off somewhat in people in their nineties. We found it necessary to use a seventh-order polynomial term to capture this shape adequately. This feature is contrary to expectation, as there is no mathematical characteristic of the QFracture-2012 model that would lead to a waning of hazard. 26 Therefore, we think it is likely that selection effects in the underlying data set led to over-representation of nonagenarians with QFracture-2012 estimates that are low relative to expectation. To explore the impact of this issue, we fitted an alternative model on a restricted data set that comprised only people younger than 90 years, allowing the model to project the observed trend from people aged 30–89 years for people in their nineties, and in this case a quartic function sufficed (see Appendix 9, Table 57 and Figure 61).

As for the non-cardiovascular death model in Chapter 6, we fitted a relative survival model to estimate the multiplicative difference in time to non-fracture death between a person in the CPRD extract and someone of the same age and sex in the general population. Again, our critical covariate is Δlogit(Q), that is the difference (on a log-odds scale) between each individual’s predicted 10-year QFracture-2012 risk of major fracture and the average score for a person of the same age and sex, as estimated in the models described above.

We fitted a simple, univariable model with Δlogit(Q) as the only covariate and a more complex one with Δlogit(Q), age, polynomial terms for both and interactions between them, as per Equation 1. We tried models with Δlogit(Q) defined using each of the QFracture-2012 prediction models in Appendix 9, Table 57, and found that the relative survival model fitted better (with meaningfully lower AIC) when we used the model fitted on only people younger than 90 years. Therefore, this was the version of the relative survival model we took forward.

Table 20 shows outputs of the relative survival models for men and women.

TABLE 20 - Relative survival models: non-fracture mortality as a function of predicted fracture risk and age in men and women

Parameter	Simple model, HR (95% CI)	Polynomial model, HR (95% CI)
Men
Δlogit(Q)	0.7216 (0.7163 to 0.7270)	1.7108 (1.6808 to 1.7409)
Δlogit(Q)2		−0.2025 (−0.2215 to −0.1835)
Age		7.69 × 10−3 (6.87 × 10−3 to 8.51 × 10−3)
Age2		1.17 × 10−5 (−1.03 × 10−5 to 3.36 × 10−5)
Δlogit(Q) × age		−0.0218 (−0.0240 to −0.0195)
Δlogit(Q)2 × age		2.05 × 10−3 (6.65 × 10−4 to 3.44 × 10−3)
Δlogit(Q) × age2		−1.66 × 10−4 (−2.13 × 10−4 to −1.19 × 10−4)
Δlogit(Q)2 × age2		2.93 × 10−5 (1.89 × 10−6 to 5.66 × 10−5)
AIC	5,252,183	5,239,595
Women
Δlogit(Q)	0.8588 (0.8507 to 0.8668)	1.5458 (1.5081 to 1.5834)
Δlogit(Q)2		0.2280 (0.1987 to 0.2572)
Age		−3.07 × 10−3 (−3.97 × 10−3 to −2.16 × 10−3)
Age2		2.19 × 10−4 (1.97 × 10−4 to 2.41 × 10−4)
Δlogit(Q) × age		−1.04 × 10−3 (−3.64 × 10−3 to 1.55 × 10−3)
Δlogit(Q)2 × age		−0.0163 (−0.0184 to −0.0142)
Δlogit(Q) × age2		−5.30 × 10−4 (−5.79 × 10−4 to −4.80 × 10−4)
Δlogit(Q)2 × age2		2.07 × 10−4 (1.66 × 10−4 to 2.48 × 10−4)
AIC	5,172,844	5,164,693

Figure 22 shows estimates of non-fracture mortality generated in the way described above, compared with empirical data across different categories of cardiovascular risk and age. The empirical Kaplan–Meier curves represent observed non-fracture mortality (censored for fracture death) in the CPRD extract for people of the stated age, with baseline QFracture-2012 predictions in the specified brackets (see Figure 22). For our modelled estimates, we start from ONS 3-year lifetables for England and Wales (we use 2009–11 tables, as this is in the middle of the period covered by the CPRD data). We adjust the data to remove fracture deaths, estimated using proportions recorded under ICD-10 code M80 in the ONS’s ‘Deaths registered in England and Wales’ series. For each combination of risk and age bracket, we calculate HRs, as described above (for comparability, we fit at the mean QFracture-2012 prediction and age observed within that category in the CPRD data), and apply the HRs to the ONS curves, and this produces the adjusted curves shown in Figure 22.

FIGURE 22.

Fitted relative survival model for non-fracture mortality compared with empirical data. (a) Men; and (b) women.

The fitted models provide a less strikingly accurate fit to the observed data than we saw with non-cardiovascular mortality, but still represent a substantial improvement over the unadjusted estimates. Larger coefficients are estimated for the polynomial terms in the more complex model, meaning that it diverges from the simple model to a somewhat greater extent than we saw for non-cardiovascular mortality.

Some instances in which the fitted models appear to depart from the empirical data seem to be due to incoherencies in the observed curves. For instance, it appears that women aged 60–69 years with a 10-year risk of less than 1% and women a decade older with a 1–2% risk have worse survival than the models predict. However, if this is true, then it means that women at this low risk have worse survival than women with higher fracture risk, which is both hard to explain and inconsistent with the patterns seen elsewhere in the data, where increasing risk is invariably associated with worse survival.

We incorporated the relative survival models into the health economic decision model by calculating a fitted HR for each sampled patient and applying it when calculating their life expectancy.

Deriving results

For TA464, the modellers presented results from a single model run comprising 2 million virtual people (c. 200,000 per risk-stratified tenth), with parameter uncertainty propagated at the individual level (i.e. the model samples new parameter values from distributions reflecting uncertainty in their true value for each patient). Our preliminary investigation suggested that there is still non-trivial Monte Carlo error in the model with that many runs, but this generally appears to stabilise after around 4 million patients. Therefore, our base-case results summarise the outputs of 5 million people (c. 500,000 per tenth).

To enable us to explore results further, as in TA464, we fitted a generalised additive model (GAM) to the patient-level outputs of the discrete-event model. The response variable is incremental net monetary benefit for alendronate compared with no treatment. In its simplest version, the meta-model takes the form:

where g(.) is a link function (i.e. identity, as we assume a Gaussian distribution for INMB), f(.) is a smooth term (with a cubic regression spline basis function) that can capture non-linearities in the relationship between the predictor (i.e. baseline QFracture-2012 10-year serious fracture probability) and the response term, and Adj is a dummy variable indicating whether or not the economic model was adjusted for competing risk of non-fracture death.

We also extended the meta-modelling methods used in TA464 to provide an estimate of value for money as a function of baseline risk and age. This model takes the form:

Here, the joint effect of fracture risk and age becomes the target of the smooth term. As before, we stratify the model by adjustment condition.

As in previous iterations of the model,149,150 we found that GAM predictions become erratic above a predicted 10-year fracture risk of 30%, owing to the very small numbers of simulated patients who reach this level of risk (i.e. around 0.1% of the cohort). Therefore, we fitted the meta-model to a data set excluding the few simulated people who exceed this level of risk and present only fitted results up to this threshold.

Accounting for competing risk

Table 21 summarises expected events with alendronate compared with no treatment, first in the unadjusted model and then in the version that adjusts for competing risk of non-fracture death. Table 21 presents analogous information to Table 9 in the 2020 HTA. 149 Table 37 of the analysis underpinning TA464150 has something similar, although note that this is based on an uncorrected network meta-analysis, and the corrigendum does not include an updated version of the same information. 147

Whether adjusted for competing risk or not, the updated model estimates that alendronate prevents slightly more fractures than the previous models in all categories except wrist fracture, and this is because the longer duration of treatment (and consequent extension of offset period) gives more time for the intervention to provide benefit. The reduced efficacy for wrist fracture compared with the 2020 HTA149 is probably a result of different inputs, for example the network meta-analysis developed for the 2020 HTA149 included a wider range of comparators and estimated a greater mean effect for bisphosphonates in preventing wrist fractures.

Comparing adjusted with unadjusted models reveals some relatively subtle differences. In higher-risk people (most obviously in the highest-risk tenth), the unadjusted model slightly overstates the effectiveness of alendronate. We can see this when we take account of the fact that people with high risk of fracture also experience increased risk of competing causes of death, with the result that fewer of people live to see the benefit of future fractures prevented. With a large enough sample, we would also expect to see the unadjusted model underpredicting fractures in the people with below-average risk (because their true life expectancy is longer than the unadjusted model simulates). There may be a hint of this type of effect in the lowest-risk tenth of the population; however, we would need a colossal sample size to detect the effect reliably, as events are so rare (and intervening deaths even rarer) in the few years during which the treatment effect obtains.

Table 22 shows how these dynamics translate into costs and QALYs. In absolute terms, the adjusted model generates more QALYs than the unadjusted model for people at lowest risk of fracture, and the relationship reverses as risk rises. Comparing incremental results from the unadjusted and adjusted models shows no difference in the deciles in which oral bisphosphonates represent an effective use of NHS resources. In both versions of the model, the lower-risk seven-tenths of people gain a tiny benefit from alendronate, but this benefit is insufficient to offset the additional costs associated with the intervention. Above the seventh decile, however, QALY gains get a bit larger and (because fractures cause expense to the health and care system, as well as disutility to the person) incremental costs go down, with the result that preventative treatment provides reasonable value for money (when we value QALYs at £20,000 each). In the adjusted model, the degree of net benefit expected for people at highest risk is discernibly lower than in the unadjusted version, and this is a direct result of the smaller number of fractures prevented, as discussed above. Nevertheless, both models agree about the people for whom incremental net benefit is positive.

Figure 23 shows the output of the GAM meta-model (see Equation 4) fitted to patient-level output and it is directly analogous to the green line in Figure 3 in the corrigendum to the analysis that informed TA464147 (with the difference that we find it helpful to show the x-axis on a logarithmic scale; note that the first seven deciles of risk all fit into the first sixth of the natural scale graph). This is consistent with what we have seen in previous model outputs, that is the unadjusted model overestimates benefit in people at greatest risk as, once we adjust for competing risk of other-cause death, their capacity to benefit from reduced fracture risk diminishes.

FIGURE 23.

Cost-effectiveness of bisphosphonates (incremental net monetary benefit compared with no treatment) as a function of risk of fracture, with different assumptions about competing risk of non-fracture death. GAM fitted to model outputs comprising 5 million simulated patients. Lines show fitted model prediction, shaded areas show 95% CI, vertical dashed bars represent deciles of risk and numbered shapes show mean values for people within each tenth of the population (where these are missing, fewer than 100 of the total 5 million simulated people fell into the group). INMB, incremental net monetary benefit.

However, as seen before, the levels of risk at which this bias operates are almost exclusively above the point where – both adjusted and unadjusted models agree that – intervention provides positive net benefit. The unadjusted model suggests that the threshold at which alendronate becomes cost-effective is a 10-year MOF risk of 4.7 (95% CI 4.5 to 4.9), whereas the adjusted model estimates the same value as 5.0 (95% CI 4.8 to 5.2).

Extending the GAM meta-model to incorporate age as well as baseline risk (see Equation 5) produces outputs such as those illustrated in Figure 24 (fitted at indicative ages of 50, 60, 70 and 80 years). By and large, results replicate those in the unstratified meta-model, that is we may overestimate value for money in people at the highest risk if we do not adjust for competing risk of non-fracture death; however, this generally only affects the magnitude of expected benefit in people for whom some degree of benefit is expected. There is some indication that these expectations may diverge for the youngest people. For example, for 50-year-olds, treatment is only cost-effective in the highest-risk tenth of people once we adjust for competing risk. Note, however, that the outputs of the meta-model are much more uncertain, in this instance, as they are based on a far smaller sample size (almost all simulated people with the highest level of fracture risk are aged > 50 years).

FIGURE 24.

Cost-effectiveness of bisphosphonates (incremental net monetary benefit compared with no treatment) as a function of age and risk of fracture, with different assumptions about competing risk of non-fracture death. (a) Age 50 years; (b) age 60 years; (c) age 70 years: and (d) age 80 years. GAM fitted to model outputs comprising 5 million simulated patients. Lines show fitted model prediction, shaded areas show 95% CI, vertical dashed bars represent deciles of risk and numbered shapes show mean values for people within each tenth of the population (where these are missing, fewer than 100 of the total 5 million simulated people fell into the group). INMB, incremental net monetary benefit.

Appendix 9, Table 58, shows the threshold at which treatment becomes associated with positive net benefit according to the GAM, adjusting for age and fracture risk. Unadjusted model outputs are relatively invariant to age, whereas the adjusted model suggests that the threshold for intervention should fall as people get older (until they reach the oldest category).

As we did with the statins model, we also looked at how it would affect cost-effectiveness results if we simulated people who are fitter or less fit than average (over and above the degree that would be expected via fracture risk alone), by halving and doubling each individual’s HR for non-fracture death. Figure 25 shows results for model 2, which adjusts for fracture risk only. Figure 26 shows the same for model 2 but adjusts for age as well. Figure 27 shows a cross-categorisation of age and fracture risk, with varying assumptions about competing risk. As would be expected, we can see that the threshold at which intervention generates positive incremental net benefit is somewhat higher in people who are more likely to die of other things and somewhat lower in people whose increased life expectancy gives them every chance of surviving to realise the benefit of fractures prevented.

FIGURE 25.

Cost-effectiveness of bisphosphonates (incremental net monetary benefit compared with no treatment) as a function of risk of fracture, with different assumptions about competing risk of non-fracture death. GAM fitted to model outputs comprising 5 million simulated patients. Lines show fitted model prediction, shaded areas show 95% CI, vertical dashed bars represent deciles of risk and numbered shapes show mean values for people within each tenth of the population (where these are missing, fewer than 100 of the total 5 million simulated people fell into the group). INMB, incremental net monetary benefit.

FIGURE 26.

Cost-effectiveness of bisphosphonates (incremental net monetary benefit compared with no treatment) as a function of age and risk of fracture, with different assumptions about competing risk of non-fracture death. (a) Age 50 years; (b) age 60 years; (c) age 70 years; and (d) age 80 years. GAM fitted to model outputs comprising 5 million simulated patients. Lines show fitted model prediction, shaded areas show 95% CI, vertical dashed bars represent deciles of risk and numbered shapes show mean values for people within each tenth of the population (where these are missing, fewer than 100 of the total 5 million simulated people fell into the group). INMB, incremental net monetary benefit.

FIGURE 27.

Relationship between age, risk of fracture and cost-effectiveness of bisphosphonates (incremental net monetary benefit compared with no treatment), with different assumptions about competing risk of non-fracture death. (a) Unadjusted (as per TA464); (b) adjusted for competing risk of non-fracture death; (c) adjusted with doubled hazard of non-fracture death; and (d) adjusted with halved hazard of non-fracture death. INMB, incremental net monetary benefit.

Incorporating direct treatment disutility

We show the lifetime discounted QALY losses we would expect from DTD in Table 23. Permanent and time-limited losses are identical because, in our base case, we assume that time-limited DTD lasts for 10 years and alendronate therapy is capped at 5 years. Therefore, under both scenarios, people experience full DTD for as long as they take the drug.

To give context to our estimates of DTD, we should remember that, as shown in Table 22, once we adjust for competing risk of non-fracture death, expected QALY gains associated with alendronate compared with no treatment are small (with no more than 0.005 QALYs per person), even in people at the highest risk of fracture. Even under the assumption of diminishing DTD, we estimate that people stand to lose more than 10 times as many QALYs from DTD as they stand to gain from the fracture-preventing effect of alendronate. Any version of Figure 27 that also accounted for any degree of DTD would comprise exclusively black area, as there is no combination of patient characteristics that leads to expected QALY gains that come close to justifying estimated DTD.

In the competing risk-adjusted model, fewer than 0.2% of simulated people experience QALY gains that are greater than their expected DTD (in permanent, time-limited and diminishing DTD scenarios). Even in the 10% of people at highest risk of fracture, these percentages rise to only 0.53–0.58%, depending on DTD scenario, and this suggests that if DTD is present then at least 199 out of 200 people treated with bisphosphonates would experience more harm than benefit, even if their chances of fracture are high.

Even if we assume that people get used to taking bisphosphonates over a short period of 1 year, QALY losses associated with DTD are still at least five times greater than expected QALY gains. The conclusion is inescapable: for anybody experiencing any duration of DTD of the magnitude estimated in Chapter 5, it is impossible to anticipate long-term benefits that would offset the up-front harm.

Main findings

Although we would reiterate that our results should not be seen as a present-day best estimate of the cost-effectiveness of bisphosphonates, our results are closely comparable with the results generated in previous iterations of the same underlying model. 147,149,150 Alendronate is associated with small but positive net benefit in people with approximately the highest-third of fracture risk (i.e. anything above around 5% per year).

As with the cohort-level model explored in Chapter 6, ‘adjusting for competing risk’ in a discrete event simulation does not mean introducing a new concept that has, hitherto, been absent from such analyses. Rather, ‘adjusting for competing risk’ means taking the competing time-to-event functions the model has always simulated and ensuring that we parameterise the functions in a way that accounts for important correlations between them. We find that ‘adjusting for competing risk’ has a relatively subtle effect on model outputs, as it makes a small difference to the magnitude of expected benefit, but the people for whom adjustment makes the greatest absolute difference are people for whom adjusted and unadjusted models predict at least some degree of benefit. This is consistent with findings in our analysis of risk prediction models that miscalibration is often most obvious in people whose level of risk clearly exceeds intervention thresholds (see Chapters 2–4).

In contrast, the effects of DTD in the same decision space are anything but subtle. At the level at which we measured the effects of DTD for bisphosphonates (see Chapter 5), DTD of any duration would be enough to swamp expected benefit from fracture prevention. The low rates of persistence with oral bisphosphonates observed in practice are further evidence that people associate bisphosphonates with disutility of a magnitude that outweighs any anticipated gain, which lends some face validity to our findings. It follows that, once we factor cost into the equation, it is impossible to find any identifiable group of people for whom oral bisphosphonates represent an effective use of NHS resources if we assume population-level average DTD for everyone to whom the decision applies. In fact, there are a small number of simulated people for whom bisphosphonates would be cost-effective at conventional QALY values despite these people experiencing net harm (as this can occur when non-trivial cost savings arise from prevented fractures and associated QALY gains are marginally less than DTD-related QALY losses).

Nevertheless, in the same way we argued for statins for the primary prevention of CVD (see Chapter 6, Discussion), we would not suggest that decision-makers should advocate blanket disinvestment in bisphosphonates. On the one hand, we think it is useful information that the average person would sacrifice more quality-adjusted life expectation to avoid taking bisphosphonates than they could expect to gain from their therapeutic effect (and this is true both of people with experience of taking bisphosphonates and members of the general public). On the other hand, although our survey respondents were more eager to avoid bisphosphonates than statins, around one-sixth of the respondents still ascribed no DTD to bisphosphonates. Therefore, if adequately informed people consider the trade-off differently to the average person (and they are at sufficiently high risk of fragility fracture to justify the costs of treatment), then we would still want to offer people access to a treatment from which they could expect some benefit before any counterweighting from disadvantages they are prepared to tolerate.

Limitations

First, our epidemiological work (see Chapter 3) finds that QFracture-2012 is poorly calibrated for reasons not limited to its inability to deal with competing risks, and, on the face of it, this would seem to undermine our cost-effectiveness model that, in common with its previous iterations, predicts risk of fracture using baseline QFracture-2012 predictions. However, in the world simulated by our decision-analytic model, QFracture-2012 does not have poor predictive utility but, instead, has perfect predictive utility because the model simulates events as a direct function of each simulated person’s baseline risk (i.e. the simulated people in the model really do have the risk of fracture that QFracture-2012 would ascribe to them). Therefore, our results can be seen as assessing the cost-effectiveness of bisphosphonates in a population with accurate risk prediction. This means that our results would be valid for the more accurate risk estimates that a better risk prediction model would produce. The only areas in which QFracture-2012’s poor calibration may undermine the model are (1) the deciles of risk by which the population is subdivided bear little relation to the true risks faced in the population and (2) for TA464, the modellers derived the distributional assumptions underpinning time-to-fracture events from analysis of QFracture-2012 data (see ‘Estimating time to event from absolute fracture risk’ in Davis et al. 150). If a more accurate risk prediction model implied different functional forms and/or shape parameters, then this ought to be reflected in an updated cost-effectiveness model.

Second, as with our relative survival model for non-cardiovascular death, we acknowledge that collapsing risk to a single covariate (in this case, 10-year major fracture risk) runs the risk of lumping together people with different non-fracture survival expectations. For example, a 70-year-old woman with no long-term conditions but a family history of osteoporosis has an identical QFracture-2012 prediction to a 70-year-old woman with type 2 diabetes and CKD but no family history (as both have a 12.6% chance of major fragility fracture over 10 years). It is plausible that, although their risk of fracture is indistinguishable, these two profiles would be associated with markedly different life expectancy. Our relative survival model is blind to such dynamics, as the model uses fracture risk as a single covariate. When comparing our adjusted life expectancy estimates with observed data (see Figure 22), we acknowledge that, for some age–sex–fracture risk strata, there are larger discrepancies between modelled and observed outcomes than we saw with the cardiovascular data set. We might be able to achieve a closer fit with a more sophisticated model (i.e. instead of using a summary risk prediction measure we could estimate relative survival as a function of each of the individual covariates on which the prediction itself relies). Such an approach would theoretically be able to capture the differential effect of risk factors for fracture on life expectancy. However, this would be a very complex model to fit, and we could only validate it against ever-more stratified subsections of the empirical data, which, even when starting from a large data set, would swiftly lead to sample-size constraints. At very least, we remain confident that our simpler modelling approach results in much better estimates of life expectancy than relying on unadjusted general population data, as has been done in all previous models.

Third, our results suggest that, after one accounts for competing risk of non-fracture death, the risk threshold for intervention with bisphosphonates goes down as people get older. This is the opposite of what is recommended in NICE Quality Standard 149 (QS149),164 where intervention thresholds rise with age. The thresholds in QS149164 derive from an analysis by McCloskey et al. ,165 the predominant aim of which was to reduce the number of people in whom treatment is indicated. Our modelling suggests that a younger person with higher risk of fracture has less capacity to benefit because they are likely to be substantially less healthy than an average person of the same age and sex, with the consequence that their life expectancy leaves less room for fracture prevention. An older person with the same estimated risk will not stand out from their contemporaries in the same way and so their life expectancy will be closer to typical.

Fourth, the modelled population has a lower proportion of people with a high fracture risk than the real population because the simulation does not account for correlations between risk factors. For example, when the model generates a simuland with type 2 diabetes, then that person is no more likely than an average member of the population to have CKD or CVD, whereas, in reality, such risk factors cluster together. This has the effect that there are unrepresentatively few people with multiple risk factors and, hence, higher fracture risk in the simulated data set. This issue may have caused or compounded the problems we encountered fitting meta-models to full data sets, leading us to truncate the data to people with predicted 10-year fracture risks lower than 30% (see Deriving results).

Finally, given the computational demands of producing a probabilistic sensitivity analysis, this was not possible for this study.

Conclusion

Similar to the analysis presented in Chapter 6, this analysis has shown that including competing risk in the model-based CEA of bisphosphonates for primary prevention produced only subtle changes to the observed cost-effectiveness. However, we noted some effect and this analysis provides more evidence that competing risk should be included in CEAs of primary prevention strategies. Incorporating DTD had a dramatic effect on the cost-effectiveness of bisphosphonates. The effect was so dramatic that the observed QALY gains from taking a bisphosphonate were swamped by DTD in all scenarios. This result raises an interesting challenge to decision-makers in the context of bisphosphonates. However, rather than advise bisphosphonates should not be recommended because of the possible DTD, we advise that the impact of including DTD in sensitivity analyses on the outputs of a decision-analytic model should be an integral component of CEAs of primary prevention strategies.