ddhazard: Fitting dynamic hazard models

Description

Function to fit dynamic hazard models using state space models.

Usage

ddhazard(formula, data, model = "logit", by, max_T, id, a_0, Q_0, Q = Q_0,
  order = 1, weights, control = list(), verbose = F)

Arguments

formula

coxph like formula with Surv(tstart, tstop, event) on the left hand site of ~.

data

data.frame or environment containing the outcome and co-variates.

model

"logit" or "exponential" for the logistic link function in the first case or for the continuous time model in the latter case.

interval length of the bins in which parameters are fixed.

max_T

end of the last interval interval.

vector of ids for each row of the in the design matrix.

a_0

vector \(a_0\) for the initial coefficient vector for the first iteration (optional). Default is estimates from static model (see static_glm).

Q_0

covariance matrix for the prior distribution.

initial covariance matrix for the state equation.

order

order of the random walk.

weights

weights to use if e.g. a skewed sample is used.

control

list of control variables (see the control section below).

verbose

TRUE if you want status messages during execution.

Value

A list with class ddhazard. The list contains

formula: the passed formula.
call: the matched call.
state_vecs: 2D matrix with the estimated state vectors (regression parameters) in each bin.
state_vars: 3D array with smoothed variance estimates for each state vector.
lag_one_cov: 3D array with lagged correlation matrix for each for each change in the state vector. Only present when the model is logit and the method is EKF.
n_risk: the number of observations in each interval.
times: the interval borders.
risk_set: the object from get_risk_obj if saved.
data: the data argument if saved.
weights: weights used in estimation if saved.
id: ids used to match rows in data to individuals.
order: order of the random walk.
F_: matrix which map from one state vector to the next.
method: method used in the E-step.
est_Q_0: TRUE if Q_0 was estimated in the EM-algorithm.
family: Rcpp Module with C++ functions used for estimation given the model argument.
discrete_hazard_func: the hazard function corresponding to the model argument.
terms: the terms object used.
has_fixed_intercept: TRUE if the model has a time-invariant intercept.
xlev: a record of the levels of the factors used in fitting.

Control

The control argument allows you to pass a list to select additional parameters. See vignette("ddhazard", "dynamichazard") for more information. Unspecified elements of the list will yield default values.

method: set to the method to use in the E-step. Either "EKF" for the Extended Kalman Filter, "UKF" for the Unscented Kalman Filter, "SMA" for the sequential posterior mode approximation method or "GMA" for the global mode approximation method. "EKF" is the default.
LR: learning rate.
NR_eps: tolerance for the Extended Kalman filter. Default is NULL which means that no extra iteration is made in the correction step.
alpha: hyper parameter \(\alpha\) in the unscented Kalman Filter.
beta: hyper parameter \(\beta\) in the unscented Kalman Filter.
kappa: hyper parameter \(\kappa\) in the unscented Kalman Filter.
n_max: maximum number of iteration in the EM-algorithm.
eps: tolerance parameter for the EM-algorithm.
est_Q_0: TRUE if you want the EM-algorithm to estimate Q_0. Default is FALSE.
save_risk_set: TRUE if you want to save the list from get_risk_obj used to estimate the model. It may be needed for later calls to e.g., residuals, plot and logLike.
save_data: TRUE if you want to keep the data argument. It may be needed for later calls to e.g., residuals, plot and logLike.
denom_term: term added to denominators in either the EKF or UKF.
n_threads: maximum number of threads to use.
fixed_parems_start: starting value for fixed terms.
fixed_terms_method: the method used to estimate the fixed effects. Either 'M_step' or 'E_step' for estimation in the M-step or E-step respectively.
Q_0_term_for_fixed_E_step: the diagonal value of the initial covariance matrix, Q_0, for the fixed effects if fixed effects are estimated in the E-step.
eps_fixed_parems: tolerance used in the M-step of the Fisher's scoring algorithm for the fixed effects.
permu: TRUE if the risk sets should be permutated before computation. This is TRUE by default for posterior mode approximation method and FALSE for all other methods.
posterior_version: the implementation version of the posterior approximation method. Either "woodbury" or "cholesky".
GMA_max_rep: maximum number of iterations in the correction step if method = 'GMA'.
GMA_NR_eps: tolerance for the convergence criteria for the relative change in the norm of the coefficients in the correction step if method = 'GMA'.

Details

This function can be used to estimate survival models where the regression parameters follows a given order random walk. The order is specified by the order argument. 1. and 2. order random walks is implemented. The regression parameters are updated at time by, 2by, ..., max_T. See the vignette("ddhazard", "dynamichazard") for details.

All filter methods needs a state covariance matrix Q_0 and state vector a_0. An estimate from a time-invariant model is used for a_0 if it is not supplied (the same model you would get from static_glm). A diagonal matrix with large entries is recommended for Q_0. What is large dependents on the data set and model. Further, a covariance matrix for the first iteration Q is needed. The Q and a_0 are estimated with an EM-algorithm.

The model is specified through the model argument. The logistic model is where outcomes are binned into the intervals. Be aware that there can be "loss" of information due to binning. It is key for the logit model that the id argument is provided if individuals in the data set have time-varying co-variates. The the exponential model use an exponential model for the arrival times where there is no "loss" information due to binning.

It is recommended to see the Shiny app demo for this function by calling ddhazard_app().

References

Fahrmeir, Ludwig. Dynamic modelling and penalized likelihood estimation for discrete time survival data. Biometrika 81.2 (1994): 317-330.

Durbin, James, and Siem Jan Koopman. Time series analysis by state space methods. No. 38. Oxford University Press, 2012.

Examples

Run this code

# NOT RUN {
# example with first order model
library(dynamichazard)
fit <- ddhazard(
 Surv(time, status == 2) ~ log(bili), pbc, id = pbc$id, max_T = 3600,
 Q_0 = diag(1, 2), Q = diag(1e-4, 2), by = 50,
 control = list(method = "GMA"))
plot(fit)

# example with second order model
fit <- ddhazard(
 Surv(time, status == 2) ~ log(bili), pbc, id = pbc$id, max_T = 3600,
 Q_0 = diag(1, 4), Q = diag(1e-4, 2), by = 50,
 control = list(method = "GMA"),
 order = 2)
plot(fit)

# }

Run the code above in your browser using DataLab