R/01_simulate_process.R

In vsritter/sim.etm: What the Package Does (One Line, Title Case)

#' Title
#'
#' @param N
#' @param gen
#' @param uni_cens
#' @param init_states
#' @param init_times
#' @param max_iter
#'
#' @return
#' @export
#'
#' @examples
simulate_process <- function(N, gen,
                             uni_cens = NULL, init_states = NULL,
                             init_times = NULL, max_iter = 10) {
  t <- Matrix::Matrix(0, nrow = N, ncol = max_iter)
  if (is.null(init_times)) {
    init_times <- rep(0, N)
  }
  t[, 1] <- init_times

  states <- Matrix::Matrix(0, nrow = N, ncol = max_iter)
  if (is.null(init_states)) {
    init_states <- rep(1, N)
  }
  states[, 1] <- init_states

  trans_mat <- apply(gen, 1:2, function(x) sum(!is.na(x)) > 0)

  non_absorb_states <- which(rowSums(trans_mat) > 0)

  iter = 1
  while (iter < max_iter & sum(states[,iter] %in% non_absorb_states) > 0) {
    for (k in non_absorb_states) {
      ids <- which(states[, iter] == k)

      if (length(ids) > 0) {
        iter_out <- iterate_process(t[ids, iter], gen[k,,], uni_cens)

        states[ids, iter + 1] <- iter_out[[1]]
        t[ids, iter + 1] <- iter_out[[2]] # t[ids, iter] + iter_out[[2]]
      }
    }
    iter = iter + 1
    if (iter == max_iter)
      warning("Maximum number of iterations reached.")
  }

  sim_out <- data.frame(id = states@i + 1,
                        state = states@x,
                        time = c(init_times, t@x))

  return(sim_out)
}