Simulate expected values as a function of simulated state occupancy probabilities, with simulation of costs and quality-adjusted life-years (QALYs) as particular use cases.

# S3 method for stateprobs
sim_ev(
object,
models = NULL,
dr = 0.03,
integrate_method = c("trapz", "riemann_left", "riemann_right"),
value_name = "value",
outcome_name = "outcome",
...
)

sim_qalys(
object,
model,
dr = 0.03,
integrate_method = c("trapz", "riemann_left", "riemann_right"),
lys = TRUE
)

sim_costs(
object,
models,
dr = 0.03,
integrate_method = c("trapz", "riemann_left", "riemann_right")
)

## Arguments

object

A stateprobs object.

dr

Discount rate.

integrate_method

Method used to integrate state values when computing costs or QALYs. Options are trapz (the default) for the trapezoid rule, riemann_left left for a left Riemann sum, and riemann_right right for a right Riemann sum.

value_name

Name of the column containing values of the outcome. Default is "value".

outcome_name

Name of the column indicating the outcome corresponding to each model. Only used if models is a list. Default is "outcome".

...

Currently unused.

model, models

An object or list of objects of class StateVals used to model state values. When using sim_qalys(), this should be a single model for utility. With sim_costs(), a list of models should be used with one model for each cost category. Finally, with sim_ev(), this may either be a single model or a list of models. May also be NULL, in which case length of stay is computed based on the state probabilities contained in object.

lys

If TRUE, then life-years are simulated in addition to QALYs.

## Value

sim_ev() returns a data.table with the following columns:

sample

A random sample from the PSA.

strategy_id

The treatment strategy ID.

patient_id

The patient ID.

grp_id

The subgroup ID.

state_id

The health state ID.

dr

The rate used to discount costs.

outcome

The outcome corresponding to each model in models. Only included if models is a list.

value

The expected value.

The names of the outcome and value columns may be changed with the value_name and outcome_name arguments. sim_costs() and sim_qalys()return similar objects, that are of class costs and qalys, respectively.

## Details

Expected values in cohort models (i.e., those implemented with the CohortDtstm and Psm classes) are mean outcomes for patients comprising the cohort. The method used to simulate expected values depends on the $method field in the StateVals object(s) stored in model(s). If $method = "starting", then state values represent a one-time value that occurs at time 0.

The more common use case is $method = "wlos", or a "weighted length of stay". That is, expected values for each health state can be thought of as state values weighted by the time a patient spends in each state (and discounted by a discount factor that depends on the discount rate dr). The precise computation proceeds in four steps. In the first step, the probability of being in each health state at each discrete time point is simulated (this is the output contained in the stateprobs object). Second, a StateVals model is used to predict state values at each time point. Third an expected value at each time point is computed by multiplying the state probability, the state value, and the discount factor. Fourth, the expected values at each time point are summed across all time points. The summation in the fourth step can be thought of as a discrete approximation of an integral. In particular, the limits of integration can be partitioned into time intervals, with each interval containing a start and an end. The integrate_method argument determines the approach used for this approximation: 1. A left Riemann sum (integrate_method = "riemann_left") uses expected values at the start of each time interval. 2. A right Riemann sum (integrate_method = "riemann_right") uses expected values at the end of each time interval. 3. The trapezoid rule (integrate_method = "trapz") averages expected values at the start and end of each time interval. (This will generally be the most accurate and is recommended.) Mathematical details are provided in the reference within the "References" section below. ## Note The ID variables in the state value models in models must be consistent with the ID variables contained in object. In particular, the models should predict state values for each non-absorbing health state in object; that is, the number of health states modeled with the models should equal the number of health states in object less the number of absorbing states. The absorbing states are saved as an attribute named absorbing to stateprobs objects. When simulating state probabilities with a CohortDtstmTrans object, the absorbing state is determined by the absorbing field in the class; in a Psm (or with sim_stateprobs.survival()), the absorbing state is always equal to the final health state. ## References Incerti and Jansen (2021). See Section 2.1 for mathematical details. ## See also State probabilities can be simulated using the $sim_stateprobs() methods from either the CohortDtstmTrans (or CohortDtstm) or Psm classes. State probabilities can also be computed directly from survival curves with the generic method sim_stateprobs.survival().

Costs and QALYs are typically computed within the R6 model classes using the $sim_costs() and $sim_qalys() methods. For instance, see the documentation and examples for the CohortDtstm and Psm classes. The sim_qalys() and sim_costs() functions are exported to give users additional flexibility when creating their own modeling pipelines. sim_ev() may be useful for computing outcomes other than costs or QALYs.

costs and qalys objects can be passed to summarize_ce() to create a cost-effectiveness object for performing a cost-effectiveness analysis with cea(). Although note that typically the $summarize() method belonging to the CohortDtstm or Psm classes would be used instead. Use the IndivCtstm class to simulate costs and QALYs with an individual continuous-time state transition model. ## Examples # We need (i) a state probability object and (ii) a model for state values ## We should start by setting up our decision problem hesim_dat <- hesim_data(strategies = data.frame(strategy_id = 1:2), patients = data.frame(patient_id = 1:3), states = data.frame(state_id = 1)) input_data <- expand(hesim_dat, by = c("strategies", "patients")) ## (i) Simulate a state probability object tpmat_id <- tpmatrix_id(input_data, n_samples = 2) p_12 <- ifelse(tpmat_id$strategy_id == 1, .15, .1)
tpmat <- tpmatrix(
C, p_12,
0, 1
)
transmod <- CohortDtstmTrans$new(params = tparams_transprobs(tpmat, tpmat_id)) stprobs <- transmod$sim_stateprobs(n_cycles = 3)

## Construct model for state values
outcome_tbl <- stateval_tbl(
data.frame(
state_id = 1,
est = 5000
),
dist = "fixed"
)
outmod <- create_StateVals(outcome_tbl, n = 2, hesim_data = hesim_dat)

# We can then simulate expected values
## The generic expected values function
sim_ev(stprobs, models = outmod)
#>     sample strategy_id patient_id grp_id state_id   dr    value
#>  1:      1           1          1      1        1 0.03 11429.69
#>  2:      1           1          2      1        1 0.03 11429.69
#>  3:      1           1          3      1        1 0.03 11429.69
#>  4:      1           2          1      1        1 0.03 12346.79
#>  5:      1           2          2      1        1 0.03 12346.79
#>  6:      1           2          3      1        1 0.03 12346.79
#>  7:      2           1          1      1        1 0.03 11429.69
#>  8:      2           1          2      1        1 0.03 11429.69
#>  9:      2           1          3      1        1 0.03 11429.69
#> 10:      2           2          1      1        1 0.03 12346.79
#> 11:      2           2          2      1        1 0.03 12346.79
#> 12:      2           2          3      1        1 0.03 12346.79

## We can also pass a list of models
sim_ev(stprobs, models = list(Outcome 1 = outmod))
#>     sample strategy_id patient_id grp_id state_id   dr   outcome    value
#>  1:      1           1          1      1        1 0.03 Outcome 1 11429.69
#>  2:      1           1          2      1        1 0.03 Outcome 1 11429.69
#>  3:      1           1          3      1        1 0.03 Outcome 1 11429.69
#>  4:      1           2          1      1        1 0.03 Outcome 1 12346.79
#>  5:      1           2          2      1        1 0.03 Outcome 1 12346.79
#>  6:      1           2          3      1        1 0.03 Outcome 1 12346.79
#>  7:      2           1          1      1        1 0.03 Outcome 1 11429.69
#>  8:      2           1          2      1        1 0.03 Outcome 1 11429.69
#>  9:      2           1          3      1        1 0.03 Outcome 1 11429.69
#> 10:      2           2          1      1        1 0.03 Outcome 1 12346.79
#> 11:      2           2          2      1        1 0.03 Outcome 1 12346.79
#> 12:      2           2          3      1        1 0.03 Outcome 1 12346.79

## Suppose the outcome were a cost category. Then we might
## prefer the following:
sim_costs(stprobs, models = list(drug = outmod))
#>     sample strategy_id patient_id grp_id state_id   dr category    costs
#>  1:      1           1          1      1        1 0.03     drug 11429.69
#>  2:      1           1          2      1        1 0.03     drug 11429.69
#>  3:      1           1          3      1        1 0.03     drug 11429.69
#>  4:      1           2          1      1        1 0.03     drug 12346.79
#>  5:      1           2          2      1        1 0.03     drug 12346.79
#>  6:      1           2          3      1        1 0.03     drug 12346.79
#>  7:      2           1          1      1        1 0.03     drug 11429.69
#>  8:      2           1          2      1        1 0.03     drug 11429.69
#>  9:      2           1          3      1        1 0.03     drug 11429.69
#> 10:      2           2          1      1        1 0.03     drug 12346.79
#> 11:      2           2          2      1        1 0.03     drug 12346.79
#> 12:      2           2          3      1        1 0.03     drug 12346.79

## Length of stay is computed if there is no state value model
sim_ev(stprobs)
#>     sample strategy_id patient_id grp_id state_id    value
#>  1:      1           1          1      1        1 2.285938
#>  2:      1           1          2      1        1 2.285938
#>  3:      1           1          3      1        1 2.285938
#>  4:      1           2          1      1        1 2.469358
#>  5:      1           2          2      1        1 2.469358
#>  6:      1           2          3      1        1 2.469358
#>  7:      2           1          1      1        1 2.285938
#>  8:      2           1          2      1        1 2.285938
#>  9:      2           1          3      1        1 2.285938
#> 10:      2           2          1      1        1 2.469358
#> 11:      2           2          2      1        1 2.469358
#> 12:      2           2          3      1        1 2.469358