Draw Time-to-Event with Time-Dependent Covariates and Luck Adjustment

Simulate a time-to-event (TTE) from a parametric distribution with parameters varying over time. User provides parameter functions and distribution name. The function uses internal survival and conditional quantile functions, plus luck adjustment to simulate the event time.

Usage

qtimecov(
  luck,
  a_fun,
  b_fun = function(t) NA,
  dist,
  dt = 0.1,
  max_time = 100,
  return_luck = FALSE,
  start_time = 0
)

Arguments

luck: Numeric between 0 and 1. Initial random quantile (luck).
a_fun: Function of time t returning the first distribution parameter (e.g., rate, shape, meanlog).
b_fun: Function of time t returning the second distribution parameter (e.g., scale, sdlog). Defaults to a function returning NA.
dist: Character string specifying the distribution. Supported: "exp", "gamma", "lnorm", "norm", "weibull", "llogis", "gompertz".
dt: Numeric. Time step increment to update parameters and survival. Default 0.1.
max_time: Numeric. Max allowed event time to prevent infinite loops. Default 100.
return_luck: Boolean. If TRUE, returns a list with tte and luck (useful if max_time caps TTE)
start_time: Numeric. Time to use as a starting point of reference (e.g., curtime).

Value

Numeric. Simulated time-to-event.

Details

The objective of this function is to avoid the user to have cycle events with the only scope of updating some variables that depend on time and re-evaluate a TTE. The idea is that this function should only be called at start and when an event impacts a variable (e.g., stroke event impacting death TTE), in which case it would need to be called again at that point. In that case, the user would need to call e.g., a <- qtimecov with max_time = curtime, return_luck = TRUE arguments, and then call it again with no max_time, return_luck = FALSE, and luck = a$luck, start_time=a$tte (so there is no need to add curtime to the resulting time).

It's recommended to play with dt argument to balance running time and precision of the estimates. For example, if we know we only update the equation annually (not continuously), then we could just set dt = 1, which would make computations faster.

Examples


param_fun_factory <- function(p0, p1, p2, p3) {
  function(t) p0 + p1*t + p2*t^2 + p3*(floor(t) + 1)
}

set.seed(42)

# 1. Exponential Example
rate_exp <- param_fun_factory(0.1, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = rate_exp,
  dist = "exp"
)
#> [1] 24.52825


# 2. Gamma Example
shape_gamma <- param_fun_factory(2, 0, 0, 0)
rate_gamma <- param_fun_factory(0.2, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = shape_gamma,
  b_fun = rate_gamma,
  dist = "gamma"
)
#> [1] 22.22001


# 3. Lognormal Example
meanlog_lnorm <- param_fun_factory(log(10) - 0.5*0.5^2, 0, 0, 0)
sdlog_lnorm <- param_fun_factory(0.5, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = meanlog_lnorm,
  b_fun = sdlog_lnorm,
  dist = "lnorm"
)
#> [1] 6.555032


# 4. Normal Example
mean_norm <- param_fun_factory(10, 0, 0, 0)
sd_norm <- param_fun_factory(2, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = mean_norm,
  b_fun = sd_norm,
  dist = "norm"
)
#> [1] 11.8122


# 5. Weibull Example
shape_weibull <- param_fun_factory(2, 0, 0, 0)
scale_weibull <- param_fun_factory(10, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = shape_weibull,
  b_fun = scale_weibull,
  dist = "weibull"
)
#> [1] 10.03201


# 6. Loglogistic Example
shape_llogis <- param_fun_factory(2.5, 0, 0, 0)
scale_llogis <- param_fun_factory(7.6, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = shape_llogis,
  b_fun = scale_llogis,
  dist = "llogis"
)
#> [1] 7.736007


# 7. Gompertz Example
shape_gomp <- param_fun_factory(0.01, 0, 0, 0)
rate_gomp <- param_fun_factory(0.091, 0, 0, 0)
qtimecov(
  luck = runif(1),
  a_fun = shape_gomp,
  b_fun = rate_gomp,
  dist = "gompertz"
)
#> [1] 13.57996

#Time varying example, with change at time 8
rate_exp <- function(t) 0.1 + 0.01*t * 0.00001*t^2
rate_exp2 <- function(t) 0.2 + 0.02*t
time_change <- 8
init_luck <- 0.95

a <- qtimecov(luck = init_luck,a_fun = rate_exp,dist = "exp", dt = 0.005,
                      max_time = time_change, return_luck = TRUE)
qtimecov(luck = a$luck,a_fun = rate_exp2,dist = "exp", dt = 0.005, start_time=a$tte)
#> [1] 11.69083


#An example of how it would work in the model, this would also work with time varying covariates!
rate_exp <- function(t) 0.1
rate_exp2 <- function(t) 0.2
time_change <- 10
init_luck <- 0.95
#at start, we would just draw TTE
qtimecov(luck = init_luck,a_fun = rate_exp,dist = "exp", dt = 0.005)
#> [1] 29.95232

#at event in which rate changes (at time 10) we need to do this:
a <- qtimecov(luck = init_luck,a_fun = rate_exp,dist = "exp", dt = 0.005,
                      max_time = time_change, return_luck = TRUE)
qtimecov(luck = a$luck,a_fun = rate_exp2,dist = "exp", dt = 0.005, start_time=a$tte)
#> [1] 19.97366