Personalizing Medical Treatment Decisions: Integrating Meta-analytic Treatment Comparisons with Patient-Specific Risks and Preferences

Overview

Using data from an NMA, the modeling framework simulates 10,000 identical versions of a patient on each treatment alternative and computes treatment scores that account for both estimated values and uncertainty in the probability of disease outcomes and patient preferences. To account for patient risks, the framework simulates events and then creates scale ranges to facilitate comparison of changes in different disease outcomes. To account for patient preferences, the framework includes steps for preference elicitation, creation of a utility function that translates preferences and outcomes onto a consistent numeric scale, incorporation of the preferences and outcomes into a personalized cost-effectiveness ratio, and risk adjustment (Figure 1).

OPEN IN VIEWER

Figure 1

Overview of modeling framework. We simulate 10,000 identical versions of a patient on each treatment and compute discounted numbers of events (e.g., disease events and occurrence of severe side effects) in each case. From this simulation, for each outcome, we compute a scale range. We then elicit patient preferences and compute utility scores for each version of the patient on each treatment. Finally, we use these utilities to compute QALYs and ICERs, and then risk adjust.

Event Simulation

Using discrete event simulation, 10,000 identical versions of a patient on each treatment alternative are simulated using data from an NMA over the time period for which data are available. At each time step, the framework computes baseline (pretreatment) absolute risks of each disease outcome, adjusts for the RR of each outcome for each treatment alternative based on the NMA data, and probabilistically estimates the occurrence of each disease outcome over the simulated period. Machine learning models are used to relate patient features to the pretreatment risk of each outcome. Machine learning is a type of artificial intelligence that uses computers to automatically learn from data without being explicitly programmed. Machine-learning models have been trained and validated to estimate the probability of numerous disease outcomes using a wide variety of approaches.^21–23 If not already in the literature, models can be trained using data from the placebo/comparator arms of the NMA. We demonstrate this latter approach in our schizophrenia example, where we train random forest models and use cross-validation to estimate out-of-sample performance.

NMAs typically report results as odds ratios (ORs) or standardized mean differences (SMDs). Importantly, NMAs report relative treatment effects, which must be integrated with an individual patient’s baseline risk of each outcome to calculate absolute treatment effects. For example, a patient with a low baseline risk of a condition may not experience a substantial absolute risk reduction from a given RR reduction, whereas a patient with a high baseline risk of a condition may experience a substantial absolute risk reduction from that same RR reduction.

To adjust baseline absolute risks for each treatment alternative, the framework converts ORs of disease outcomes from the NMA of randomized trials to RRs using baseline prevalence data from published studies. To estimate each outcome probabilistically, a binomial function with probability equal to the adjusted absolute risk is used. Finally, discounted numbers of events for each version of the patient on each treatment are computed using a 3% discount rate. Thus, for each treatment, the discounted number of events is calculated for each of the 10,000 versions of the patient.

If the NMA reports SMDs of treatments versus placebo for an outcome (e.g., an outcome such as body weight, rather than a binomial outcome that could be codified as an event), the SMDs are converted to mean differences using pooled standard deviations calculated from published studies. The corresponding patient baseline risk factor values are then modified by the mean differences to simulate the effects of the different treatments.

Scale Ranges

For each outcome, 95% credible intervals (crIs) are computed for the identical versions of the patient on each treatment. The scale range is taken as the interval hull of these 95% crIs for each outcome. The interval hull of m intervals is the smallest interval that includes all m intervals. For example, if the intervals were [1, 3], [2, 5], and [4, 6], the interval hull would be [1, 6].

Preference Elicitation

There are many methods for preference elicitation. ²⁴ In general, it is unclear which methods are best in which situations. ²⁴ Any method can be used in the framework, with no additional modification, provided it can generate preferences for outcomes on a scale from 0 to 1 that sum to 1. Even if a method generates preferences on a different scale, it can still be used if the preferences can be standardized onto this scale. Preference standardization is discussed in the Utility Scores subsection.

Two common methods for preference elicitation include swing weighting and discrete choice experiments (DCE).^25,26 For swing weighting, patient preferences can be elicited using sliders, a method that has been shown to be visually intuitive and easy for patients with varied numeracy skills to understand.^27–29 Using sliders, a patient allocates weights to swings across the scale ranges such that the most important is assigned 1 and the rest are relative on a scale of 0 to 1. For example, in a case with 2 outcomes, nausea and having to inject a medication, the patient might be asked the relative importance of having no nausea as opposed to daily nausea, versus no injection as opposed to daily injections. The patient might input the weights 1 and 0.5 for the swings, respectively, in response. The sliders are marked with ticks every 0.1 increment, and the framework assumes 95% crIs from the next lowest to the next highest tick. For DCE, a patient repeatedly selects which set of outcomes he or she would prefer from 2 different hypothetical options. This information can then be used to compute preference weights and associated crIs. For a comparison of these 2 approaches and discussion regarding when they should be used, see Tervonen et al. ²⁵ In the context of a single patient, swing weighting might be preferable, as DCE might require the patient to state preferences over a large number of hypothetical alternatives.

Utility Scores

The framework assumes linear partial value functions and a linear, additive utility function over the scale ranges, although other utility functions can also be used. A partial value function gives the value over a scale range for a single outcome. These forms are the most common in the literature^6,9 and are based on 2 key assumptions. First, it is assumed that there are no significant nonlinearities in value over the scale ranges: for example, the value associated with a 1-unit change at the lower end of a scale range is approximately equal to that associated with a 1-unit change at the upper end of a scale range. Second, it is assumed that the disutility associated with 2 conditions equals the sum of their individual disutility values. The validity of these assumptions must be considered for each application, possibly through solicitation of expert opinion of patients and physicians. Of note, only linearity of the partial value functions over the scale ranges is required, which reflect probable ranges of outcomes. This is a much weaker assumption than assuming linearity over all values of all outcomes. In our schizophrenia example, we make these assumptions based on the expert opinion of the senior author.

If these assumptions are violated, an alternative utility function can be used in the framework, with no additional modification. The appropriate functional form depends on the specific application. For example, for the scale range for toxin level, a piecewise partial value function may be appropriate in which the function is relatively constant below some physiologically determined threshold and then rapidly decreases as the toxin reaches unsafe levels. If the appropriate functional form is not known, a technique such as conjoint analysis can be used. ³⁰

To compute utilities from the partial value functions, preference weights are drawn from distributions for each of the 10,000 identical versions of the patient. These values are then standardized so that the weights w i sum to 1 over all outcomes i = 1, …, n. The partial value functions u i ( c i ) , where c i is the value of outcome i, are also standardized such that, for undesirable outcomes, the value of being at the upper end of a scale range is 0 and the value of being at the lower end of a scale range is 1. Equation 1 is used to compute utilities for each of the 10,000 versions of the patient on each of the treatments:

Here, u is utility, c is the vector of outcomes, w is the vector of preference weights for swings over scale ranges, and u i ( c i ) are the partial value functions. Thus, for example, for undesirable outcomes, a patient at the upper end of all scale ranges would have a utility of 0, while a patient at the lower end of all scale ranges would have a utility of 1.

Personalized Cost-Effectiveness

The quality-adjusted life-years (QALYs) associated with moving across all scale ranges from the upper to the lower ends for each of the n outcomes, Q UL , i , are summed using equation 2:

The number Q UL is then multiplied by the utility u for a patient on a treatment (from equation 1) to obtain QALYs relative to the patient at the upper end of all the scale ranges. However, this reference point is problematic because the patient at the upper end of all the scale ranges would have zero utility and thus zero QALYs but may not be dead over the model period. Therefore, the offset between QALYs for the dead patient and the patient at the upper end of all the scale ranges, which we denote by Q DU , is added to Q UL u :

The quantity Q is QALYs experienced by the patient over the model period. To calculate the offset Q DU , the framework calculates QALYs for the patient if completely healthy over the time period under consideration versus dead ( Q DH ) and subtracts QALYs associated with moving across all scale ranges ( Q UL ) and QALYs associated with moving from the lower end of all scale ranges to completely healthy ( Q LH ) :

Figure 2 illustrates these QALY components.

OPEN IN VIEWER

Figure 2

QALY components. The total QALYs for a patient from death to full health, Q_DH, can be computed as the sum of three components: (1) Q_DU, the QALYs from death to the upper ends of the scale ranges, (2) Q_UL, the QALYs from the upper to the lower ends of the scale ranges, and (3) Q_LH, the QALYs from the lower ends of the scale ranges to full health. The scale ranges are assumed to be undesirable in that being at the upper ends results in worse health; upper and lower can be interchanged if this is not the case.

To calculate a cost-effectiveness ratio, the increase in cost is divided by QALYs gained, which is calculated as the difference in Q between 2 treatment alternatives, both estimated by equation 3. This enables the framework to convert personalized measures of utility to measures of cost-effectiveness, enabling comparison of disparate interventions. If cost-effectiveness of a population-level program is desired, incremental costs and QALYs for the population, respectively, can be summed before dividing.

Risk Adjustment

Past methods for personalizing treatment decisions have generally not accounted for patient risk preferences.^6,9 However, consider a treatment that guarantees 10 years of life and one that offers a 50/50 chance of 0 or 20 years. Even though these treatments yield the same expected number of years of life, different patients may prefer one or the other.

Numerous methods can be used to characterize patient risk preference. The framework assumes an exponential utility function based on previous proposals ³¹ :

Here, u R is risk-adjusted utility, Q is QALYs, r is risk odds (computed from a single gamble), and a and b are standardizing constants. Other utility functions can be used, as appropriate. One advantage of this approach is that a patient’s utility function can be fully defined based on the result of a single gamble. See Supplementary Text SA1 for additional details.

A simple, generic illustration of the modeling framework for the case of 3 drugs and 2 outcomes (disease events and serious adverse events) is presented in Supplementary Text SA2. This example is intended to make the modeling framework easier to understand; however, it is not necessary to read the example to understand the present study.

Outcomes

We considered 9 outcomes, including the 6 outcomes from the NMA (efficacy, weight gain, extrapyramidal side effects, prolactin increase, QTc prolongation, and sedation) and the additional outcomes of agranulocytosis (a possible adverse immunological effect of clozapine), prescription drug cost, and years of life lost. The incidence of agranulocytosis in each of the first 2 years after initiation of clozapine was obtained from the literature ³³ and assumed to be zero thereafter. ³⁴ We obtained costs for schizophrenia drugs from various sources (Supplementary Table S2) and costs for conditions associated with schizophrenia (Supplementary Table S3) from a past study. ³⁵ Years of life lost for a patient were calculated as compared with expected years of life for a healthy person of the same age and sex in the United States. ³⁶ All outcomes except agranulocytosis were assumed to be irreversible based on the senior author’s experiences treating patients with schizophrenia and current knowledge of the natural history of the disease. ³⁷ If patients develop agranulocytosis, it was assumed to affect them only in the year that they develop it as the condition is treatable.

Population Sample

We generated a set of 1000 patients demographically representative of the adult population of the United States with schizophrenia, each with their own risk factors (age, sex, race, Positive and Negative Syndrome Scale [PANSS] score, prior psychiatric hospitalizations, weight, prolactin level, and QTc interval), based on data from the Collaborative Psychiatric Epidemiology Surveys (CPES). ³⁸ From the CPES, we obtained means, standard deviations, and a correlation matrix for patient features. We randomly generated patient features using a multivariate random sample from the correlated distributions using copulas.

Personalized Risk Estimates

We generated baseline mortality rates for hypothetical versions of the patients of the same age and sex in the United States without schizophrenia based on data from the World Health Organization. ³⁶ We calculated mortality rates for these hypothetical patients with moderate schizophrenia using a mortality multiplier from the literature that accounts for common non–therapy-related causes of heightened mortality among patients with schizophrenia (e.g., pneumonia related to homelessness, interpersonal violence related to poverty, etc.). ³⁹ We then calculated mortality rates for the simulated patients by assigning both of these health states PANSS scores, based on the literature, ⁴⁰ and assuming that mortality scales linearly with PANSS score, based on the expert opinion of the senior author. We also assumed that prolonged QTc interval and agranulocytosis increase mortality, based on prior literature.^33,41

We developed models relating patient features to associated baseline absolute risks for each potential outcome using data from the placebo arms of trials cited by the NMA. ¹ We used random forest machine learning models ⁴² and leave-one-out cross-validation to estimate out-of-sample performance (Supplementary Table S4) measured by root mean square error, mean absolute error, and coefficient of determination (R²). Supplementary Figures S5 and S6 show the variable importance metrics and partial dependence plots for the random forest models.

Patient Preferences

In practice, preferences of individual patients would be elicited as described in the Preference Elicitation subsection, above. Here, we generated preferences for the simulated patients by repeatedly sampling from distributions based on a previous DCE ⁴³ and assuming independence between preferences. Specifically, we constructed gamma distributions for preferences based on reported preference means and uncertainty estimates. We used gamma distributions as these preferences were on a scale from 0 to positive infinity. After drawing sets of preferences, we standardized them so that they summed to 1. For discrete variables, no further work was required. For continuous variables, we scaled the preferences by the outcome scale ranges for the patient divided by the changes in the outcomes over which preferences were solicited in the DCE. This assumes that these preferences scale linearly with the outcomes. Of note, the changes in outcomes over which these preferences were solicited in the DCE were relatively broad, so we rarely had to extrapolate outside these ranges. Thus, for each patient, we sampled preferences from distributions for each outcome and then standardized them so they summed to 1.

Selection of Antipsychotic Drugs for Schizophrenia

We used these variables and our modeling framework to personalize the selection of antipsychotic drugs for schizophrenia. Briefly, we simulated 10,000 versions of each patient in our population sample on each treatment alternative to project potential outcomes and scale ranges. The event simulation was run for 5 years to reflect clinical decision timelines and duration of studies reported in the NMA.

We used patient preferences from a previous DCE and equation 1 to compute utilities for each of the 10,000 versions of each patient on each treatment alternative. Utilities were converted to QALYs using equation 3. We next computed mean QALYs, which were used to determine the best treatment alternative for each patient. We also calculated total US societal costs associated with each treatment plan. Outcomes, utilities, QALYs, and costs were discounted over the duration of each patient’s life at a rate of 3% annually, following current cost-effectiveness guidelines. ⁴⁴ In this example, we did not risk adjust QALYs in the base case analysis, as we did not have risk preference information specific to schizophrenia. In supplemental analysis, we considered the effects of risk adjustment.

We then calculated the incremental cost-effectiveness ratio (ICER) of personalized treatment determined using our modeling framework versus prescribing the best treatment for the population as a whole to all patients, which was also determined using our modeling framework. We adopted a US societal perspective for both costs and QALYs. We used our previously calculated QALYs to select drugs but did not use them to calculate the ICER as they were influenced by patient preferences. Instead, to calculate the ICER, we used separately calculated QALYs from simulation outcomes using societal QALY weights. Specifically, we summed incremental simulated costs and summed the incremental simulated QALYs for the sample population before dividing the 2 quantities to compute the ICER.