case weights
@mattsecrest Please have a look at the overall idea :)
This code creates a borrowing type for case weights. If any weights were specified in the outcome object they will be silently ignored. Is that the point of #309 to create another borrowing type and remove the weights from the outcome?
Using the same same example, the results are quite different from the full or HC borrowing.
cuts = c(1, 5, 10)
cw_pem_test <- create_analysis_obj(
data_matrix = example_matrix,
outcome = outcome_surv_pem("time", "cnsr", prior_normal(0, 100000), cut_points = cuts),
borrowing = borrowing_case_weights("ext", samples = c(200, 200)),
treatment = treatment_details("trt", prior_normal(0, 100000)),
covariates = add_covariates(c("cov1", "cov2"), prior_normal(0, 100000))
)
cw_pem_res <- mcmc_sample(cw_pem_test,
iter_warmup = 200,
iter_sampling = 1000,
parallel_chains = 4
)
cw_pem_res_results2 <- cw_pem_res$summary(
variables = NULL,
mean = mean,
median = median,
exp_mean = function(x) mean(exp(x))
)
@mattsecrest I think this weight calculation could be called from the prepare_stan_data_inputs method, after the PEM cast to long step. Something like this:
setMethod(
f = "prepare_stan_data_inputs",
signature = c("OutcomeSurvPEM", "BorrowingPowerPriorCaseWeights", "ANY"),
definition = function(outcome, borrowing, analysis_obj) {
analysis_obj <- cast_mat_to_long_pem(analysis_obj)
analysis_obj <- calculate_case_weights(analysis_object)
....
# Add covariates and weights
data_in <- add_covariates_and_weights(data_in, analysis_obj, data_matrix)
return(data_in)
})
Hey @gravesti, still need some time to review as I work my way thruogh the SI text. In the meantime, is your example code working?
cuts = c(1, 5, 10)
cw_pem_test <- create_analysis_obj(
data_matrix = example_matrix,
outcome = outcome_surv_pem("time", "cnsr", prior_normal(0, 100000), cut_points = cuts),
borrowing = borrowing_case_weights("ext", samples = c(200, 200)),
treatment = treatment_details("trt", prior_normal(0, 100000)),
covariates = add_covariates(c("cov1", "cov2"), prior_normal(0, 100000))
)
cw_pem_res <- mcmc_sample(cw_pem_test,
iter_warmup = 200,
iter_sampling = 1000,
parallel_chains = 4
)
cw_pem_res_results2 <- cw_pem_res$summary(
variables = NULL,
mean = mean,
median = median,
exp_mean = function(x) mean(exp(x))
)
I'm erroring from mismatched dimensions here: https://github.com/Genentech/psborrow2/blob/d47af1c23c24b988041ea225835bcbfb902e8fe5/R/case_weight_helpers.R#L231
Browse[1]> head(sim.mm)
cov1 cov2 interval1 interval2 interval3 interval4
489 1 1 1 0 0 0
490 1 1 0 1 0 0
491 1 1 0 0 1 0
492 1 1 1 0 0 0
493 1 1 0 1 0 0
494 1 0 1 0 0 0
Browse[1]> lambda
cov1 cov2 interval1 interval2 interval3 interval4
[1,] 0.5175591 0.6154345 -2.481291 -2.521592 -2.879851 -3.683704
Browse[1]> rate_t <- exp(sim.mm[impute_these_id, ] %*% lambda)
Error in sim.mm[impute_these_id, ] %*% lambda : non-conformable arguments
Other, high-level thoughts as I keep going through the SI text:
- I have no strong feelings for the API, think this looks good. I had in mind passing the weights directly to a power_prior_case_weights() as a more flexible approach (maybe they have their own weights to specify directly) but this works as well!
# My idea
cuts <- c(1, 5, 10)
case_weights <- estimate_box_p_case_weights(
data_matrix = example_matrix,
cut_points = cuts,
samples = c(200, 200),
"ext",
"trt",
c("cov1", "cov2")
)
analysis_obj <- create_analysis_obj(
data_matrix = example_matrix,
outcome = outcome_surv_pem("time", "cnsr", prior_normal(0, 100000), cut_points = cuts),
borrowing = borrowing_power_prior_cw("ext", case_weights),
treatment = treatment_details("trt", prior_normal(0, 100000)),
covariates = add_covariates(c("cov1", "cov2"), prior_normal(0, 100000))
)
If you prefer yours, let's go with it. I only preferred this b/c it separates the two steps more cleanly.
Other, high-level thoughts as I keep going through the SI text:
- I have no strong feelings for the API, think this looks good. I had in mind passing the weights directly to a power_prior_case_weights() as a more flexible approach (maybe they have their own weights to specify directly) but this works as well!
I avoided that so we don't duplicate all the data transformation, finding the variables etc. No reason we can't have another simple power_prior borrowing method that the user can specify any weight vector they like. It's just tricky to do validation on the input because we don't know how many rows there will be until quite late in the process.