Simulate health state transitions in a cohort discrete time state transition model.
An R6::R6Class object.
create_CohortDtstmTrans()
creates a CohortDtstmTrans
object from either
a fitted statistical model or a parameter object. A complete economic model can be implemented
with the CohortDtstm
class.
params
Parameters for simulating health state transitions.
Supports objects of class tparams_transprobs
or params_mlogit_list
.
input_data
An object of class input_mats
.
cycle_length
The length of a model cycle in terms of years.
The default is 1
meaning that model cycles are 1 year long.
absorbing
A numeric vector denoting the states that are
absorbing states; i.e., states that cannot be transitioned from.
Each element should correspond to a state_id
,
which should, in turn, be the index of the health state.
start_stateprobs
A non-negative vector with length equal to the number of
health states containing the probability that the cohort is in each health
state at the start of the simulation. For example,
if there were three states and the cohort began the simulation in state 1,
then start_stateprobs = c(1, 0, 0)
. Automatically normalized to sum to 1.
If NULL
, then a vector with the first element equal to 1 and
all remaining elements equal to 0.
trans_mat
A transition matrix describing the states and transitions
in a discrete-time multi-state model. Only required if the model is
parameterized using multinomial logistic regression. The (i,j)
element
represents a transition from state i
to state j
. Each possible transition
from row i
should be based on a separate multinomial logistic regression
and ordered from 0
to K - 1
where K
is the number of
possible transitions. Transitions that are not possible should be NA
.
and the reference category for each row should be 0
.
new()
Create a new CohortDtstmTrans
object.
CohortDtstmTrans$new(
params,
input_data = NULL,
trans_mat = NULL,
start_stateprobs = NULL,
cycle_length = 1,
absorbing = NULL
)
params
The params
field.
input_data
The input_data
field.
trans_mat
The trans_mat
field.
start_stateprobs
The start_stateprobs
field.
cycle_length
The cycle_length
field.
absorbing
The absorbing
field. If NULL
, then the constructor
will determine which states are absorbing automatically; non NULL
values
will override this behavior.
sim_stateprobs()
Simulate probability of being in each health state during each model cycle.
An object of class stateprobs
.
library("msm")
library("data.table")
set.seed(101)
# We consider two examples that have the same treatment strategies and patients.
# One model is parameterized by fitting a multi-state model with the "msm"
# package; in the second model, the parameters are entered "manually" with
# a "params_mlogit_list" object.
# MODEL SETUP
strategies <- data.table(
strategy_id = c(1, 2, 3),
strategy_name = c("SOC", "New 1", "New 2")
)
patients <- data.table(patient_id = 1:2)
hesim_dat <- hesim_data(
strategies = strategies,
patients = patients
)
# EXAMPLE #1: msm
## Fit multi-state model with panel data via msm
qinit <- rbind(
c(0, 0.28163, 0.01239),
c(0, 0, 0.10204),
c(0, 0, 0)
)
fit <- msm(state_id ~ time, subject = patient_id,
data = onc3p[patient_id %in% sample(patient_id, 100)],
covariates = list("1-2" =~ strategy_name),
qmatrix = qinit)
## Simulation model
transmod_data <- expand(hesim_dat)
transmod <- create_CohortDtstmTrans(fit,
input_data = transmod_data,
cycle_length = 1/2,
fixedpars = 2,
n = 2)
transmod$sim_stateprobs(n_cycles = 2)
#> sample strategy_id patient_id grp_id state_id t prob
#> <num> <int> <int> <int> <num> <num> <num>
#> 1: 1 1 1 1 1 0.0 1.000000000
#> 2: 1 1 1 1 1 0.5 0.792695102
#> 3: 1 1 1 1 1 1.0 0.628365524
#> 4: 1 1 1 1 2 0.0 0.000000000
#> 5: 1 1 1 1 2 0.5 0.202046511
#> ---
#> 104: 2 3 2 1 2 0.5 0.184425447
#> 105: 2 3 2 1 2 1.0 0.325475890
#> 106: 2 3 2 1 3 0.0 0.000000000
#> 107: 2 3 2 1 3 0.5 0.004553194
#> 108: 2 3 2 1 3 1.0 0.016768466
# EXAMPLE #2: params_mlogit_list
## Input data
transmod_data[, intercept := 1]
#> strategy_id patient_id strategy_name intercept
#> <num> <int> <char> <num>
#> 1: 1 1 SOC 1
#> 2: 1 2 SOC 1
#> 3: 2 1 New 1 1
#> 4: 2 2 New 1 1
#> 5: 3 1 New 2 1
#> 6: 3 2 New 2 1
transmod_data[, new1 := ifelse(strategy_name == "New 1", 1, 0)]
#> strategy_id patient_id strategy_name intercept new1
#> <num> <int> <char> <num> <num>
#> 1: 1 1 SOC 1 0
#> 2: 1 2 SOC 1 0
#> 3: 2 1 New 1 1 1
#> 4: 2 2 New 1 1 1
#> 5: 3 1 New 2 1 0
#> 6: 3 2 New 2 1 0
transmod_data[, new2 := ifelse(strategy_name == "New 2", 1, 0)]
#> strategy_id patient_id strategy_name intercept new1 new2
#> <num> <int> <char> <num> <num> <num>
#> 1: 1 1 SOC 1 0 0
#> 2: 1 2 SOC 1 0 0
#> 3: 2 1 New 1 1 1 0
#> 4: 2 2 New 1 1 1 0
#> 5: 3 1 New 2 1 0 1
#> 6: 3 2 New 2 1 0 1
## Parameters
n <- 10
transmod_params <- params_mlogit_list(
## Transitions from stable state (stable -> progression, stable -> death)
stable = params_mlogit(
coefs = list(
progression = data.frame(
intercept = rnorm(n, -0.65, .1),
new1 = rnorm(n, log(.8), .02),
new2 = rnorm(n, log(.7, .02))
),
death = data.frame(
intercept = rnorm(n, -3.75, .1),
new1 = rep(0, n),
new2 = rep(0, n)
)
)
),
## Transition from progression state (progression -> death)
progression = params_mlogit(
coefs = list(
death = data.frame(
intercept = rnorm(n, 2.45, .1),
new1 = rep(0, n),
new2 = rep(0, n)
)
)
)
)
transmod_params
#> A "params_mlogit_list" object
#>
#> Summary of coefficients:
#> from to term mean sd 2.5%
#> <char> <char> <char> <num> <num> <num>
#> 1: stable progression intercept -0.6634794 0.05869755 -0.7388932
#> 2: stable progression new1 -0.2152838 0.01251555 -0.2372537
#> 3: stable progression new2 -0.2360914 1.24487406 -2.0402193
#> 4: stable death intercept -3.7448951 0.09512508 -3.8635279
#> 5: stable death new1 0.0000000 0.00000000 0.0000000
#> 6: stable death new2 0.0000000 0.00000000 0.0000000
#> 7: progression death intercept 2.4654098 0.07608317 2.3275534
#> 8: progression death new1 0.0000000 0.00000000 0.0000000
#> 9: progression death new2 0.0000000 0.00000000 0.0000000
#> 97.5%
#> <num>
#> 1: -0.5587877
#> 2: -0.1993740
#> 3: 1.9666138
#> 4: -3.6104667
#> 5: 0.0000000
#> 6: 0.0000000
#> 7: 2.5418337
#> 8: 0.0000000
#> 9: 0.0000000
#>
#> Number of parameter samples: 10
#> Number of starting (non-absorbing) states: 2
#> Number of transitions by starting state: 2 1
## Simulation model
tmat <- rbind(c(0, 1, 2),
c(NA, 0, 1),
c(NA, NA, NA))
transmod <- create_CohortDtstmTrans(transmod_params,
input_data = transmod_data,
trans_mat = tmat, cycle_length = 1)
transmod$sim_stateprobs(n_cycles = 2)
#> sample strategy_id patient_id grp_id state_id t prob
#> <num> <int> <int> <int> <num> <num> <num>
#> 1: 1 1 1 1 1 0 1.000000000
#> 2: 1 1 1 1 1 1 0.644542208
#> 3: 1 1 1 1 1 2 0.415434658
#> 4: 1 1 1 1 2 0 0.000000000
#> 5: 1 1 1 1 2 1 0.338329178
#> ---
#> 536: 10 3 2 1 2 1 0.842673723
#> 537: 10 3 2 1 2 2 0.191956639
#> 538: 10 3 2 1 3 0 0.000000000
#> 539: 10 3 2 1 3 1 0.003998375
#> 540: 10 3 2 1 3 2 0.784533916
# \dontshow{
pb <- expmat(coef(fit)$baseline)[, , 1]
## From stable
b1 <- log(pb[1, 2]/(1 - pb[1, 2] - pb[1, 3]))
b2 <- log(pb[1, 3]/(1 - pb[1, 2] - pb[1, 3]))
exp(b1)/(1 + exp(b1) + exp(b2))
#> [1] 0.3379348
exp(b2)/(1 + exp(b1) + exp(b2))
#> [1] 0.01522282
### From progression
b <- qlogis(pb[2, 2])
# }