Nothing
## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>"
)
# Numeric notation in this document
# With thanks to https://stackoverflow.com/questions/18965637/set-global-thousand-separator-on-knitr/18967590#18967590
knitr::knit_hooks$set(
inline = function(x) {
if(!is.numeric(x)){
x
} else if(x<100) {
prettyNum(round(x, 3), big.mark=",")
} else {
prettyNum(round(x, 0), big.mark=",")
}
}
)
## ----setup, echo=TRUE, message=FALSE------------------------------------------
library(dplyr)
library(lubridate)
library(heemod)
library(dynamicpv)
## ----constants, include=FALSE-------------------------------------------------
# Time constants
days_in_year <- 365.25
days_in_week <- 7
cycle_years <- days_in_week / days_in_year # Duration of a one week cycle in years
# Time horizon (years) and number of cycles
thoz <- 20
Ncycles <- ceiling(thoz/cycle_years)
# Real discounting
disc_year <- 0.03 # Per year
disc_cycle <- (1+disc_year)^(cycle_years) - 1 # Per cycle
# Price inflation
infl_year <- 0.025 # Per year
infl_cycle <- (1+infl_year)^(cycle_years) - 1 # Per cycle
# Nominal discounting
nomdisc_year <- (1+disc_year)*(1+infl_year) - 1
nomdisc_cycle <- (1+nomdisc_year)^(cycle_years) - 1 # Per cycle
## ----heemod, include=FALSE----------------------------------------------------
# State names
state_names = c(
progression_free = "PF",
progression = "PD",
death = "Death"
)
# PFS distribution for SoC with Exp() distribution and mean of 50 weeks
surv_pfs_soc <- heemod::define_surv_dist(
distribution = "exp",
rate = 1/50
)
# OS distribution for SoC with Lognorm() distribution, meanlog = 4.5, sdlog = 1
# This implies a mean of exp(4 + 0.5 * 1^2) = exp(4.5) = 90 weeks
surv_os_soc <- heemod::define_surv_dist(
distribution = "lnorm",
meanlog = 4,
sdlog = 1
)
# PFS and OS distributions for new
surv_pfs_new <- heemod::apply_hr(surv_pfs_soc, hr=0.5)
surv_os_new <- heemod::apply_hr(surv_os_soc, hr=0.6)
# Define partitioned survival model, soc
psm_soc <- heemod::define_part_surv(
pfs = surv_pfs_soc,
os = surv_os_soc,
terminal_state = FALSE,
state_names = state_names
)
# Define partitioned survival model, soc
psm_new <- heemod::define_part_surv(
pfs = surv_pfs_new,
os = surv_os_new,
terminal_state = FALSE,
state_names = state_names
)
# Parameters
params <- heemod::define_parameters(
# Discount rate
disc = disc_cycle,
# Disease management costs
cman_pf = 80,
cman_pd = 20,
# Drug acquisition costs - the SoC regime only uses SoC drug, the New regime only uses New drug
cdaq_soc = dispatch_strategy(
soc = 400,
new = 0
),
cdaq_new = dispatch_strategy(
soc = 0,
new = 1500
),
# Drug administration costs
cadmin = dispatch_strategy(
soc = 50,
new = 75
),
# Adverse event risks
risk_ae = dispatch_strategy(
soc = 0.08,
new = 0.1
),
# Adverse event average costs
uc_ae = 2000,
# Subsequent treatments
csubs = dispatch_strategy(
soc = 1200,
new = 300
),
# Health state utilities
u_pf = 0.8,
u_pd = 0.6,
)
# Define PF states
state_PF <- heemod::define_state(
# Costs for the state
cost_daq_soc = discount(cdaq_soc, disc_cycle),
cost_daq_new = discount(cdaq_new, disc_cycle),
cost_dadmin = discount(cadmin, disc_cycle),
cost_dman = discount(cman_pf, disc_cycle),
cost_ae = risk_ae * uc_ae,
cost_subs = 0,
cost_total = cost_daq_soc + cost_daq_new + cost_dadmin + cost_dman + cost_ae + cost_subs,
# Health utility, QALYs and life years
pf_year = discount(cycle_years, disc_cycle),
life_year = discount(cycle_years, disc_cycle),
qaly = discount(cycle_years * u_pf, disc_cycle)
)
# Define PD states
state_PD <- heemod::define_state(
# Costs for the state
cost_daq_soc = 0,
cost_daq_new = 0,
cost_dadmin = 0,
cost_dman = discount(cman_pd, disc_cycle),
cost_ae = 0,
cost_subs = discount(csubs, disc_cycle),
cost_total = cost_daq_soc + cost_daq_new + cost_dadmin + cost_dman + cost_ae + cost_subs,
# Health utility, QALYs and life years
pf_year = 0,
life_year = heemod::discount(cycle_years, disc_cycle),
qaly = heemod::discount(cycle_years * u_pd, disc_cycle)
)
# Define Death state
state_Death <- heemod::define_state(
# Costs are zero
cost_daq_soc = 0,
cost_daq_new = 0,
cost_dadmin = 0,
cost_dman = 0,
cost_ae = 0,
cost_subs = 0,
cost_total = cost_daq_soc + cost_daq_new + cost_dadmin + cost_dman + cost_ae + cost_subs,
# Health outcomes are zero
pf_year = 0,
life_year = 0,
qaly = 0,
)
# Define strategy for SoC
strat_soc <- heemod::define_strategy(
transition = psm_soc,
"PF" = state_PF,
"PD" = state_PD,
"Death" = state_Death
)
# Define strategy for new
strat_new <- heemod::define_strategy(
transition = psm_new,
"PF" = state_PF,
"PD" = state_PD,
"Death" = state_Death
)
# Create heemod model
heemodel <- heemod::run_model(
soc = strat_soc,
new = strat_new,
parameters = params,
cycles = Ncycles,
cost = cost_total,
effect = qaly,
init = c(1, 0, 0),
method = "life-table"
)
## ----dynpricing, include=FALSE------------------------------------------------
# Dates
# Date of calculation = 1 September 2025
doc <- lubridate::ymd("20250901")
# Date of LOE for SoC = 1 January 2028
loe_soc_start <- lubridate::ymd("20280101")
# Maturation of SoC prices by LOE + 1 year, i.e. = 1 January 2029
loe_soc_end <- lubridate::ymd("20290101")
# Date of LOE for new treatment = 1 January 2031
loe_new_start <- lubridate::ymd("20310101")
# Maturation of new treatment prices by LOE + 1 year, i.e. = 1 January 2032
loe_new_end <- lubridate::ymd("20320101")
# Effect of LoEs on prices once mature
loe_effect_soc <- 0.7
loe_effect_new <- 0.5
# Calculation of weeks since DoC for LoEs and price maturities
wk_start_soc <- floor((loe_soc_start-doc) / lubridate::dweeks(1))
wk_end_soc <- floor((loe_soc_end-doc) / lubridate::dweeks(1))
wk_start_new <- floor((loe_new_start-doc) / lubridate::dweeks(1))
wk_end_new <- floor((loe_new_end-doc) / lubridate::dweeks(1))
# Price maturity times
wk_mature_soc <- wk_end_soc - wk_start_soc
wk_mature_new <- wk_end_new - wk_start_new
# Create a tibble of price indices of length 2T, then pull out columns as needed
# We only need of length T for now, but need of length 2T for future calculations later
pricetib <- dplyr::tibble(
model_time = 1:(2*Ncycles),
model_year = model_time * cycle_years,
static = 1,
geninfl = (1 + infl_cycle)^(model_time - 1),
loef_soc = pmin(pmax(model_time - wk_start_soc, 0), wk_mature_soc) / wk_mature_soc,
loef_new = pmin(pmax(model_time - wk_start_new, 0), wk_mature_new) / wk_mature_new,
dyn_soc = geninfl * (1 - loe_effect_soc * loef_soc),
dyn_new = geninfl * (1 - loe_effect_new * loef_new)
)
# Price indices required for calculations
prices_oth <- pricetib$geninfl
prices_static <- pricetib$static
prices_dyn_soc <- pricetib$dyn_soc
prices_dyn_new <- pricetib$dyn_new
## ----dynuptake, include=FALSE-------------------------------------------------
# Time for uptake to occur
uptake_years <- 2
# Uptake vector for non-dynamic uptake
uptake_single <- c(1, rep(0, Ncycles-1))
# Uptake vector for dynamic uptake
uptake_weeks <- round(uptake_years / cycle_years)
share_multi <- c((1:uptake_weeks)/uptake_weeks, rep(1, Ncycles-uptake_weeks))
uptake_multi <- rep(1, Ncycles) * share_multi
## ----static2, include=FALSE---------------------------------------------------
# Pull out the payoffs of interest from oncpsm
payoffs <- get_dynfields(
heemodel = heemodel,
payoffs = c("cost_daq_new", "cost_daq_soc", "cost_total", "qaly", "life_year"),
discount = "disc"
) |>
dplyr::mutate(
model_years = model_time * cycle_years,
# Derive costs other than drug acquisition, as at time zero
cost_nondaq = cost_total - cost_daq_new - cost_daq_soc,
# ... and at the start of each timestep
cost_nondaq_rup = cost_total_rup - cost_daq_new_rup - cost_daq_soc_rup
)
# Create and view dataset for SoC
hemout_soc <- payoffs |> dplyr::filter(int=="soc")
head(hemout_soc)
# Create and view dataset for new intervention
hemout_new <- payoffs |> dplyr::filter(int=="new")
head(hemout_new)
## ----bimshare-----------------------------------------------------------------
# BIM settings
bi_horizon_yrs <- 5
bi_horizon_wks <- round(bi_horizon_yrs / cycle_years)
bi_discount <- 0
# Newly eligible patients
newly_eligible <- rep(1, Ncycles)
# Time for uptake to occur
uptake_years <- 2
uptake_weeks <- round(uptake_years / cycle_years)
# Market share of new intervention
share_multi <- c((1:uptake_weeks)/uptake_weeks, rep(1, Ncycles-uptake_weeks))
# Newly eligible patients receiving each intervention, "world with"
uptake_new <- newly_eligible * share_multi
uptake_soc <- newly_eligible - uptake_new
## ----bim_wout1----------------------------------------------------------------
# World without new intervention
# SoC, drug acquisition costs
wout1_soc_daqcost <- dynpv(
uptakes = newly_eligible,
payoffs = hemout_soc$cost_daq_soc_rup,
horizon = bi_horizon_wks,
prices = prices_static,
discrate = bi_discount
)
# SoC, other costs
wout1_soc_othcost <- dynpv(
uptakes = newly_eligible,
payoffs = hemout_soc$cost_nondaq_rup,
horizon = bi_horizon_wks,
prices = prices_static,
discrate = bi_discount
)
# Total budgetary costs
budget_wout1 <- budget_wout1_soc <- wout1_soc_daqcost + wout1_soc_othcost
## ----bim_with1----------------------------------------------------------------
# World with
# SoC, drug acquisition costs
with1_soc_daqcost <- dynpv(
uptakes = uptake_soc,
payoffs = hemout_soc$cost_daq_soc_rup,
horizon = bi_horizon_wks,
prices = prices_static,
discrate = bi_discount
)
# SoC, other costs
with1_soc_othcost <- dynpv(
uptakes = uptake_soc,
payoffs = hemout_soc$cost_nondaq_rup,
horizon = bi_horizon_wks,
prices = prices_static,
discrate = bi_discount
)
# New intervention, drug acquisition costs
with1_new_daqcost <- dynpv(
uptakes = uptake_new,
payoffs = hemout_new$cost_daq_new_rup,
horizon = bi_horizon_wks,
prices = prices_static,
discrate = bi_discount
)
# New intervention, other costs
with1_new_othcost <- dynpv(
uptakes = uptake_new,
payoffs = hemout_new$cost_nondaq_rup,
horizon = bi_horizon_wks,
prices = prices_static,
discrate = bi_discount
)
# Total
budget_with1_soc <- with1_soc_daqcost + with1_soc_othcost
budget_with1_new <- with1_new_daqcost + with1_new_othcost
budget_with1 <- budget_with1_soc + budget_with1_new
## ----warnlength---------------------------------------------------------------
# The uptake vector for the new intervention is long
length(trim_vec(uptake_new))
# The uptake vector for the SoC is short, once trimmed of excess zeros
length(trim_vec(uptake_soc))
## ----bim_iwith1---------------------------------------------------------------
# Budget impact
bi1_soc <- budget_with1_soc - budget_wout1_soc
bi1_new <- budget_with1_new
bi1 <- budget_with1 - budget_wout1
summary(bi1)
## ----bim_wout2----------------------------------------------------------------
# World without new intervention
# SoC, drug acquisition costs
wout2_soc_daqcost <- dynpv(
uptakes = newly_eligible,
payoffs = hemout_soc$cost_daq_soc_rup,
horizon = bi_horizon_wks,
prices = prices_dyn_soc,
discrate = bi_discount
)
# SoC, other costs - unchanged from static calculations
wout2_soc_othcost <- wout1_soc_othcost
# Total budgetary costs
budget_wout2 <- budget_wout2_soc <- wout2_soc_daqcost + wout2_soc_othcost
## ----bim_with2----------------------------------------------------------------
# World with
# SoC, drug acquisition costs
with2_soc_daqcost <- dynpv(
uptakes = uptake_soc,
payoffs = hemout_soc$cost_daq_soc_rup,
horizon = bi_horizon_wks,
prices = prices_dyn_soc,
discrate = bi_discount
)
# SoC, other costs
with2_soc_othcost <- with1_soc_othcost
# New intervention, drug acquisition costs
with2_new_daqcost <- dynpv(
uptakes = uptake_new,
payoffs = hemout_new$cost_daq_new_rup,
horizon = bi_horizon_wks,
prices = prices_dyn_new,
discrate = bi_discount
)
# New intervention, other costs
with2_new_othcost <- with1_new_othcost
# Total
budget_with2_soc <- with2_soc_daqcost + with2_soc_othcost
budget_with2_new <- with2_new_daqcost + with2_new_othcost
budget_with2 <- budget_with2_soc + budget_with2_new
## ----bim_iwith2---------------------------------------------------------------
# Budget impact
bi2_soc <- budget_with2_soc - budget_wout2_soc
bi2_new <- budget_with2_new
bi2 <- budget_with2 - budget_wout2
summary(bi2)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.