# NOT RUN {
# head-and-neck cancer study data. See Efron, B. (1988) doi:10.2307/2288857
is_censored <- c(
6, 27, 34, 36, 42, 46, 48:51, 51 + c(15, 30:28, 33, 35:37, 39, 40, 42:45))
head_neck_cancer <- data.frame(
id = 1:96,
stop = c(
1, 2, 2, rep(3, 6), 4, 4, rep(5, 8),
rep(6, 7), 7, 8, 8, 8, 9, 9, 10, 10, 10, 11, 14, 14, 14, 15, 18, 18, 20,
20, 37, 37, 38, 41, 45, 47, 47,
2, 2, 3, rep(4, 4), rep(5, 5), rep(6, 5),
7, 7, 7, 9, 10, 11, 12, 15, 16, 18, 18, 18, 21,
21, 24, 25, 27, 36, 41, 44, 52, 54, 59, 59, 63, 67, 71, 76),
event = !(1:96 %in% is_censored),
group = factor(c(rep(1, 45 + 6), rep(2, 45))))
# fit model
set.seed(61364778)
ctrl <- PF_control(
N_fw_n_bw = 500, N_smooth = 2500, N_first = 2000,
n_max = 1, # set to one as an example
n_threads = max(parallel::detectCores(logical = FALSE), 1),
eps = .001, Q_tilde = as.matrix(.3^2), est_a_0 = FALSE)
pf_fit <- suppressWarnings(
PF_EM(
survival::Surv(stop, event) ~ ddFixed(group),
data = head_neck_cancer, by = 1, Q_0 = 1, Q = 0.1^2, control = ctrl,
max_T = 30))
# the log-likelihood in the final iteration
(end_log_like <- logLik(pf_fit))
# gives the same
fw_ps <- PF_forward_filter(
survival::Surv(stop, event) ~ ddFixed(group), N_fw = 500, N_first = 2000,
data = head_neck_cancer, by = 1, Q_0 = 1, Q = 0.1^2,
a_0 = pf_fit$a_0, fixed_effects = -0.5370051,
control = ctrl, max_T = 30, seed = pf_fit$seed)
all.equal(c(end_log_like), c(logLik(fw_ps)))
# will differ since we use different number of particles
fw_ps <- PF_forward_filter(
survival::Surv(stop, event) ~ ddFixed(group), N_fw = 1000, N_first = 3000,
data = head_neck_cancer, by = 1, Q_0 = 1, Q = 0.1^2,
a_0 = pf_fit$a_0, fixed_effects = -0.5370051,
control = ctrl, max_T = 30, seed = pf_fit$seed)
all.equal(c(end_log_like), c(logLik(fw_ps)))
# will differ since we use the final estimates
fw_ps <- PF_forward_filter(pf_fit, N_fw = 500, N_first = 2000)
all.equal(c(end_log_like), c(logLik(fw_ps)))
# }
Run the code above in your browser using DataLab