We make decisions everyday from getting up in the morning to figuring out what we want to have for dinner. In most cases, we have a good idea of what we like or want based on previous experiences. For example, you might prefer to have ice cream for dinner because you like it. But what if you also believe that having ice cream is not the healthiest option for dinner. Instead, you could have vegetables for dinner because it’s healthier (or so you believe).
These decisions are mostly determined by the most current information at our disposal. But in many instances, there is a great amount of uncertainty. Sometimes, this uncertainty can impact the decision we will make. For instance, what if we know with 100% certainty that having ice cream for dinner was going to cause a heart attack afterwards? Would you still want to have the ice cream? What if that certainty was reduced to 1%? Would you change your mind about having the ice cream?
This is a common problem in health care where different treatments have a lot of uncertainty in their effectiveness. When you go see a doctor for an illness, you will usually receive some kind of treatment after the doctor diagnoses your illnes. But what if the doctor misdiagnose your illness? Or what if the doctor correctly diagnoses the illness but the treatment isn’t effective? These questions are common in health care and difficult to predict with a high degree of uncertainty.
In the field of pharmacoeconomics, decision trees are used to predict what would happen after a decision has been made about a treatment strategy. Often times, these decisions are driven by the current evidence usually in the form of randomized controlled trials (RCTs). Sometimes the evidence comes from a meta-analysis or a clinician’s experience. But what is almost always in a decision tree is the outcomes (or payoffs). Each treatment strategy will involve a series of pathways that will ultimately lead to some kind of outcome. These outcomes can be good or bad. Some examples of outcomes common with decision tree analysis are survival, response, life years, and quality-adjusted life years (QALYs). Another element that is also evaluated is the cost of the treatment strategy. Just like the outcomes, each treatment strategy will have costs associated with each pathway. Summing these up for each pathway, a pharmacoeconomist can determine the total cost of that treatment strategy.
Pharmacoeconomists can use decision trees to evaluate the value of a treatment strategy by measuring the outcomes alongside the costs associated with each treatment pathway and comparing these to other treatment strategies (e.g., usual care or standard care). This is called a cost-effectiveness analysis, which is the systematic evaluation of the value of a treatment strategy compared to an alternative treatment strategy. In this tutorial, I will show you how to create simple decision trees to predict the cost-effectivenss of a treatment strategy compared to an existing strategy such as standard care.
A decision tee provides a visial diagram of the possible pathways for a treatment strategy. Usually a decision tree will be illustrated going from left to right and will start with a decision node (square) and end with a terminal node (triangle). Other elements of a decision tree include the chance nodes (circle), which denotes the probability of an event occurring. These symbols are common conventions for decision trees and each intersection is referred to as a node. In Figure 1, a decision is made at the decision node (square) and then there is a probability that it is effective as indicated by the chance node (circle); the decision tree will have to eventually end with an outcome at the terminal node (triangle).
Now that we know the basic elements of a decision tree, let’s take a look at an example.
Suppose that a patient sees their doctor for an acute respiratory infection (ARI). The doctor can treat the patient with two types of antibiotics (Treatment A and Treatment B). The doctor can forgoe treatment with Treatment B and instead recommend that the patient received Treatment A.
The decision tree for an ARI scenario is illustrated in Figure 2. The decision node (square) is to treat with Treatment A or Treatment B. The chance nodes (circle) indicate the part of the decision tree where the patient can have some probability of having an adverse reaction and/or a cure. The terminal node (triangle) captures the outcomes for the decision pathway. (Although unnecessary, it is prudent to have similar decision pathways for the two treatment strategies. This makes it easier to compare the results and reduces potential bias due to the structure of the decision tree.)
Each decision will have its own path. For example, if a patient is treated with Treatment A, they can experience no adverse reaction and be cured (Pathway 1). Similarly, a patient can be treated with Treatment A, experence no adverse reaction but does not get cured (Pathway 2). In this figure (Figure 3), there are a total of 8 possible pathways (4 pathways for each treatment strategy).
Suppose a patient who is treated with Treatment A has an 20% probability of not having an adverse reaction. What is the probability of having an adverse reaction? The answer is 80%. Here is how I calculated the answer:
\[\begin{equation} 80\% = 100\% - 20\% \end{equation}\]Remember that the total probability of any event happening and not happing is 100 percent. So if the probability that the patient will not experience an adverse reaction after receiving Treatment A is 80 percent, then the probability of experiencing an adverse reaction must be 20 percent.
Let’s further assume that a patient who receives Treatment A and does not experience an adverse reaction has a 90% probability of being cured. What is the expected probability of Pathway 1? The answer is 72%. Here is how I calculated the answer:
\[\begin{equation} 72\% = 80\% * 90\% \end{equation}\]Once you figure out the probability of the different treatment pathways, you can insert them in the decision tree diagram. Figure 4 depicts the ARI decision tree with the probabilities for the different treatment pathways.
The expected probabilities for each pathway will be used to estimate the expected costs and benefits for each strategy. To check to see if the calculations are correct, summarize probabilities for each treatment. These should sum up to 100% as shown in Figure 4.
We can represent these probabilities in a table (Table A). Table A lists all the parameters used in the decision tree model including the probabilities, costs, and outcomes (e.g., Life Years).
Since each decision tree pathway has a probability associated with it, we can start to figure out how much each decision path costs. For example, what is the cost of Treatment A -> No adverse reaction -? Cured? Using Table A, we can aggregate the costs:
\[\begin{equation} $100 = $100 + $0 + $0 \end{equation}\]Therefore, a patient who received Treatment A, has no adverse reaction, and is cured accrues a total cost of $100.
Similary, a patient who receives Treatment A, has no adverse reaction, but is not cured accrues a total cost of $250.
\[\begin{equation} $250 = $100 + $0 + $150 \end{equation}\]We do this for each pathway, which amounts to a total of 8 pathways. Figure 5 shows how much each pathway costs.
We do the same thing for the life years associated with eachout pathway (also called the payoffs for each pathway). For this model, we will assume that a patient will generate 20 Life Years if they experience a cure and only 10 Life Years if they are not cured. Figure 6 illustrates the payoffs in terms of life years for each pathway.
Once we have the costs and life years determined for each pathway, we can estimate the expected costs and expected life years associated with Treatment A and Treatment B. To generate the expecte costs for Pathway 1, we multiply the probability of Pathway 1 (0.72) with the cost of Pathway 1 ($100). This generates the expected costs for Pathway 1 ($72.00).
Similarly, we can estimate the expected life years for Pathway 1 by multiplying the probability of Pathway 1 (0.72) with the life years of Pathway 1 (20 life years). This generates the expected life years for Pathway 1 (0.2 life years).
Table B summarizes the expected costs and life years for the remaining pathways.
After the expected costs and life years are calcualted for each pathway, you can sum these up for each treatment strategy. Thefore, Treatment A has a total cost of $215 and Treatment B has a total cost of $260. Similarly, the total life years for Treatments A and B are 19 Life Years and 18 Life Years, respectively. Figure 7 illustrates how these are summarized.
The incremental costs and life years can be calculated using the total costs and total life years. The incremental costs between Treatment A and Treatment B are:
\[\begin{equation} Incremental\ Costs = C_{TxA} - C_{TxB} \end{equation}\] \[\begin{equation} Incremental\ Life\ Years = LY_{TxA} - LY_{TxB} \end{equation}\]Once the incremental costs and life years are estimated, you can calculate the incremental cost-effectiveness ratio (ICER), which is the main result of a cost-effectiveness analysis.
\[\begin{equation} ICER = \frac{C_{TxA} - C_{TxB}}{LY_{TxA} - LY_{TxB}} \end{equation}\]Table C summarizes the results from the decision tree.
The ICER is -$45 per Life Year gained. Since the numerator is negative and the denominator is positive, the ICER falls into the southeast quadrant of a cost-effectiveness plane. This indicates that Treatment A is a dominant strategy compared to Treatment B. Figure 8 illustrates where the ICER comparing Treatment A and Treatment B for ARI falls on a cost-effectiveness plane.
The following R code will generate a table with the Total Costs, Total Life Years, Incremental Costs, Incremental Life Years, and the Incremental Cost-Effectiveness Ratio (ICER).
################################
#### Vector of parameter inputs
################################
input <- data.frame(
p.TxA_ADR = 0.20, # Probability of adverse reaction with Treatment A
p.TxB_ADR = 0.30, # Probability of advrese reaction with Treatment B
p.TxA_Cured = 0.90, # Probability of being Cured with Treatment A
p.TxB_Cured = 0.80, # Probability of being Cured with Treatment B
ly.Cured = 20, # Life years after being cured
ly.NotCured = 10, # Life years after not being cured
c.TxA = 100, # Costs of Treatment A
c.TxB = 80, # Costs of Treatment B
c.ADR = 500, # Costs of adverse reaction
c.NoADR = 0, # Costs of not having an adverse reaction
c.NotCured = 150, # Costs of not being cured
c.Cured = 0, # Costs of being cured
wtp = 100 # Willingess to pay per life year gained
)
##################################################
#### Wrap decision tree into a function ####
##################################################
dec_tree <- function(params){
with(
as.list(params),
{
# Expected probabilities for each pathway
### Pathways for Treatment A
ep1 <- (1 - p.TxA_ADR) * p.TxA_Cured # Expected probability for Pathway 1
ep2 <- (1 - p.TxA_ADR) * (1 - p.TxA_Cured) # Expected probability for Pathway 2
ep3 <- p.TxA_ADR * p.TxA_Cured # Expected probability for Pathway 3
ep4 <- p.TxA_ADR * (1 - p.TxA_Cured) # Expected probability for Pathway 4
### Pathways for Treatment B
ep5 <- (1 - p.TxB_ADR) * p.TxB_Cured # Expected probability for Pathway 5
ep6 <- (1 - p.TxB_ADR) * (1 - p.TxB_Cured) # Expected probability for Pathway 6
ep7 <- p.TxB_ADR * p.TxB_Cured # Expected probability for Pathway 7
ep8 <- p.TxB_ADR * (1 - p.TxB_Cured) # Expected probability for Pathway 8
# Total costs for each pathway (unweighted)
### Total costs for Treatment A
tc1 <- c.TxA + c.NoADR + c.Cured # Total costs for Pathway 1
tc2 <- c.TxA + c.NoADR + c.NotCured # Total costs for Pathway 2
tc3 <- c.TxA + c.ADR + c.Cured # Total costs for Pathway 3
tc4 <- c.TxA + c.ADR + c.NotCured # Total costs for Pathway 4
### Total costs for Treatment B (unweighted)
tc5 <- c.TxB + c.NoADR + c.Cured # Total costs for Pathway 5
tc6 <- c.TxB + c.NoADR + c.NotCured # Total costs for Pathway 6
tc7 <- c.TxB + c.ADR + c.Cured # Total costs for Pathway 7
tc8 <- c.TxB + c.ADR + c.NotCured # Total costs for Pathway 8
# Expected Total Costs for each Treatment strategy is the sum of the weighted values
# (probabilities of each pathway multiplied by the total costs of each pathway)
### Expected Total Costs for Treatment A
etc.TxA <- (ep1 * tc1) + (ep2 * tc2) + (ep3 * tc3) + (ep4 * tc4)
### Expected Total Costs for Treatment B
etc.TxB <- (ep5 * tc5) + (ep6 * tc6) + (ep7 * tc7) + (ep8 * tc8)
# Expected Life Years for each Treatment strategy is the sum of the weighted values
# (probabilities of each pathway multiplied by the Life Years of each pathway)
### Expected Life Years for Treatment A
ely.TxA <- (ep1 * ly.Cured) + (ep2 * ly.NotCured) + (ep3 * ly.Cured) + (ep4 * ly.NotCured)
### Expected Life Years for Treatment B
ely.TxB <- (ep5 * ly.Cured) + (ep6 * ly.NotCured) + (ep7 * ly.Cured) + (ep8 * ly.NotCured)
# Expected total costs, expected life years, incremental costs, incremental life years , and ICER lists
C <- c(etc.TxA, etc.TxB)
LY <- c(ely.TxA, ely.TxB)
IC <- etc.TxA - etc.TxB
IE <- ely.TxA - ely.TxB
ICER <- (etc.TxA - etc.TxB)/ (ely.TxA - ely.TxB)
names(C) <- paste("C", c("TxA", "TxB"), sep = "_")
names(LY) <- paste("LY", c("TxA", "TxB"), sep = "_")
names(IC) <- paste("Incr Costs")
names(IE) <- paste("Incr Life Years")
names(ICER) <- paste("ICER")
# Generate the ouput
return(c(C, LY, IC, IE, ICER))
}
)
}
#### Now, we can use the function "dec_tree" with the inputs to estimate the ICER and its corresponding values
dec_tree(input)
## C_TxA C_TxB LY_TxA LY_TxB Incr Costs
## 215 260 19 18 -45
## Incr Life Years ICER
## 1 -45
A one-way sensitivity analysis evaluates the impact of parameter uncertainty on the decision model results. Changing one parameter’s values across a specific range allows decision modelers to investigate its impact on the output. For example, if we are interested in the impact of changing the Probability of Cured for Treatment A on the Total Costs, we can vary this probablity across a range of 0% to 100%. Since we are only changing the Probability of Cured for Treatment A, it is expected that only the Total Costs of Treatment A will change. The Total Costs of Treatment B will remain the same.
An important feature of a one-way sensitivity analysis is the break-even point or the point where the Total Costs of Treatments A and B are equal. In Figure 9, when the probability of Treatment A is 0.60, its Total Costs is equal to the Total Costs of Treatment B.
The following R code will generate the one-way sensitivity analysis from this example.
##################################################
#### One-way SA
##################################################
# In the base-case, the Probability of Cured for Treatment A is 0.90. The range is 0.00 to 1.00.
p.TxA_Cured_range <- seq(0.00, 1.00, length.out=11) # Perform 11 calculations between 0.00 and 0.100.
p.TxA_Cured_range
## [1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
p.TxB_Cured_range <- seq(0.00, 1.00, length.out=11) # Perform 11 calculations between 0.00 and 0.100.
p.TxB_Cured_range
## [1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
## Generate matrix of inputs for decision tree (Replace variable 3 which is TxA_Cured with the ranges)
m.owsa.input <- cbind(p.TxA_Cured = p.TxA_Cured_range, input[-3])
m.owsa.input
## p.TxA_Cured p.TxA_ADR p.TxB_ADR p.TxB_Cured ly.Cured ly.NotCured c.TxA c.TxB
## 1 0.0 0.2 0.3 0.8 20 10 100 80
## 2 0.1 0.2 0.3 0.8 20 10 100 80
## 3 0.2 0.2 0.3 0.8 20 10 100 80
## 4 0.3 0.2 0.3 0.8 20 10 100 80
## 5 0.4 0.2 0.3 0.8 20 10 100 80
## 6 0.5 0.2 0.3 0.8 20 10 100 80
## 7 0.6 0.2 0.3 0.8 20 10 100 80
## 8 0.7 0.2 0.3 0.8 20 10 100 80
## 9 0.8 0.2 0.3 0.8 20 10 100 80
## 10 0.9 0.2 0.3 0.8 20 10 100 80
## 11 1.0 0.2 0.3 0.8 20 10 100 80
## c.ADR c.NoADR c.NotCured c.Cured wtp
## 1 500 0 150 0 100
## 2 500 0 150 0 100
## 3 500 0 150 0 100
## 4 500 0 150 0 100
## 5 500 0 150 0 100
## 6 500 0 150 0 100
## 7 500 0 150 0 100
## 8 500 0 150 0 100
## 9 500 0 150 0 100
## 10 500 0 150 0 100
## 11 500 0 150 0 100
## Run model with captured Total Costs across the ranges 0.00 to 1.00 for probability of cured with Treatment A
outcomes_TC <- t(apply(m.owsa.input, 1, dec_tree))[ , 1:2] # t(x) is the transpose function [column 1 = C_TxA in the output]
outcomes_TC
## C_TxA C_TxB
## [1,] 350 260
## [2,] 335 260
## [3,] 320 260
## [4,] 305 260
## [5,] 290 260
## [6,] 275 260
## [7,] 260 260
## [8,] 245 260
## [9,] 230 260
## [10,] 215 260
## [11,] 200 260
#### ========== NOTE: Understanding the apply function ========== ####
apply.table1 <- (apply(m.owsa.input, 1, dec_tree)) # Note: This will give us 7 rows and 11 columns.
apply.table1
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
## C_TxA 350.00 335.00000 320 305 290.0 275 260 245 230 215
## C_TxB 260.00 260.00000 260 260 260.0 260 260 260 260 260
## LY_TxA 10.00 11.00000 12 13 14.0 15 16 17 18 19
## LY_TxB 18.00 18.00000 18 18 18.0 18 18 18 18 18
## Incr Costs 90.00 75.00000 60 45 30.0 15 0 -15 -30 -45
## Incr Life Years -8.00 -7.00000 -6 -5 -4.0 -3 -2 -1 0 1
## ICER -11.25 -10.71429 -10 -9 -7.5 -5 0 15 -Inf -45
## [,11]
## C_TxA 200
## C_TxB 260
## LY_TxA 20
## LY_TxB 18
## Incr Costs -60
## Incr Life Years 2
## ICER -30
apply.table2 <- (apply(m.owsa.input, 1, dec_tree))[ , 1] # Note: You get different columns than the transpose matrix; we don't want apply.table. We want outcomes_TC. The [ , 1] denotes the 1 column from apply.table1.
apply.table2
## C_TxA C_TxB LY_TxA LY_TxB Incr Costs
## 350.00 260.00 10.00 18.00 90.00
## Incr Life Years ICER
## -8.00 -11.25
#### ============================================================= ####
##################################################
## Plot OWSA using base R commands
##################################################
Strategies <- c("Treatment A", "Treatment B")
plot(p.TxA_Cured_range, outcomes_TC[ , 1], type="l", xlab="Prob of being Cured", ylab="Total Costs")
lines(p.TxA_Cured_range, outcomes_TC[ , 2], col="red")
legend("topright", Strategies, col=1:2, lty=c(1, 1), bty="n")
Instead of having to create multiple one-way sensitivity analyses plots, you can summarize these into a tornado diagram. A tornado diagram can rank the impact of the parameter on the model’s outcome (e.g., ICER). Rank is determined by spread or the amount of change from the base-case ICER. In our example, the base-case ICER is -$45 per Life Years gained. Any change from that base-case ICER can be captured as bars on a tornado diagram.
Figure 10 illustrates how the tornado diagram reflects the lower and upper limits of the ranges used for each parameter. For instance if we change the Probability of Adverse Reaction for Treatment B between 0.225 and 0.375, the ICER values range from -$7.5 per Life Year gained to -$82.5 per Life Year gained.
devtools::source_url("https://github.com/mbounthavong/Decision_Analysis/blob/master/tornado_diagram_code.R?raw=TRUE")
################################
#### Data inputs (parameters)
################################
p.TxA_ADR = 0.20 # Probability of adverse reaction with Treatment A
p.TxB_ADR = 0.30 # Probability of advrese reaction with Treatment B
p.TxA_Cured = 0.90 # Probability of being Cured with Treatment A
p.TxB_Cured = 0.80 # Probability of being Cured with Treatment B
ly.Cured = 20 # Life years after being cured
ly.NotCured = 10 # Life years after not being cured
c.TxA = 100 # Costs of Treatment A
c.TxB = 80 # Costs of Treatment B
c.ADR = 500 # Costs of adverse reaction
c.NoADR = 0 # Costs of not having an adverse reaction
c.NotCured = 150 # Costs of not being cured
c.Cured = 0 # Costs of being cured
wtp = 100 # Willingess to pay per life year gained
########################
#### Tornado Plot #2
########################
# Define ranges
p.TxA_ADR_range <- c(BaseCase = p.TxA_ADR, low = 0.15, high = 0.25)
p.TxB_ADR_range <- c(BaseCase = p.TxB_ADR, low = 0.225, high = 0.375)
ly.Cured_range <- c(BaseCase = ly.Cured, low = 15, high = 25)
ly.NotCured_range <- c(BaseCase = ly.NotCured, low = 7.5, high = 12.5)
c.TxA_range <- c(BaseCase = c.TxA, low = 75, high = 125)
c.TxB_range <- c(BaseCase = c.TxB, low = 60, high = 100)
c.ADR_range <- c(BaseCase = c.ADR, low = 375, high = 625)
c.NotCured_range <- c(BaseCase = c.NotCured, low = 113, high = 188)
## Parameter names
paramNames <- c( "p.TxA_ADR",
"p.TxB_ADR",
"ly.Cured",
"ly.NotCured",
"c.TxA",
"c.TxB",
"c.ADR",
"c.NotCured"
)
## List of inputs
l.tor.in <- vector("list", 8)
names(l.tor.in) <- paramNames
l.tor.in$p.TxA_ADR <- cbind(p.TxA_ADR = p.TxA_ADR_range, input[-1])
l.tor.in$p.TxB_ADR <- cbind(p.TxB_ADR = p.TxB_ADR_range, input[-2])
l.tor.in$ly.Cured <- cbind(ly.Cured = ly.Cured_range, input[-5])
l.tor.in$ly.NotCured <- cbind(ly.NotCured = ly.NotCured_range, input[-6])
l.tor.in$c.TxA <- cbind(c.TxA = c.TxA_range, input[-7])
l.tor.in$c.TxB <- cbind(c.TxB = c.TxB_range, input[-8])
l.tor.in$c.ADR <- cbind(c.ADR = c.ADR_range, input[-9])
l.tor.in$c.NotCured <- cbind(c.NotCured = c.NotCured_range, input[-11])
## List of outputs
l.tor.out <- vector("list", 8)
names(l.tor.out) <- paramNames
## Run model on different parameters
# NOTE: we select [ , 7] because that is the location of the ICER output.
for(i in 1:8){
l.tor.out[[i]] <- t(apply(l.tor.in[[i]], 1, dec_tree))[ , 7]
}
## Data structure: ymean, ymin, ymax
m.tor <- matrix(unlist(l.tor.out), nrow = 8, ncol = 3, byrow = TRUE,
dimnames = list(paramNames, c("basecase", "low", "high")))
TornadoPlot(main_title = "Tornado Plot", Parms = paramNames, Outcomes = m.tor,
outcomeName = "Incremental Cost-Effectiveness Ratio (ICER)",
xlab = "ICER",
ylab = "Parameters",
col1="#3182bd", col2="#6baed6")
Decision trees are useful in estimating the values between two or more strategies. In this example, we used a decision tree to model the potential pathways for treating ARI. Based on our decision tree, the ICER was -$45 per Life Year gained. Using R, we were able to evaluate the impact of varying the Probability of Cured of Treatment A across a reasonable range and its impact on the Total Costs. We also used R to generate a tornado diagram to rank the impact of each parameter uncertainty on the ICER. Based on our torando diagram, the Probability of Adverse Reaction for Treatment B has the largest impact on the base-case ICER.
Using decisions trees to project the consequences of treatment decisions can help us determine whether the strategy we choose will be cost-effective relative to the reference case. However, careful assessment of the parameters’ uncertainty on the decision tree’s outcomes (e.g., ICER) is necessary to determine the robustness of the model’s conclusions.
I would like to acknowledge the researchers who taught the Decision Modeling Using R workshop at the University of Washington on February 2016. They are:
Petros Pechlivanoglou PhD
Fernando Alarid-Escudero MSc, PhD(c)
Szu-Yu Zoe Kao MSc
The Hospital for Sick Children, Toronto, ON, Canada
University of Toronto, ON, Canada
University of Minnesota School of Public Health, Minneapolis, MN, USA
Their contributions to the field of decision anlaysis and R were instrumental in creating this tutorial.
I modified the codes that they developed to create this tutorial.