working

R/make_model_string_model.R

@@ -156,9 +156,9 @@ setMethod(
     has_covariates <- !is.null(analysis_obj@covariates)
 
     linear_predictor <- if (has_covariates) {
-      h_glue("lp = Z*alpha_ext' + (1-Z)*alpha' + rep_matrix(trt * beta_trt + X * beta, M);")
+      h_glue("lp = Z*alpha' + rep_matrix(trt * beta_trt + X * beta, M);")
     } else if (!has_covariates) {
-      h_glue("lp = Z*alpha_ext' + (1-Z)*alpha' + rep_matrix(trt * beta_trt, M);")
+      h_glue("lp = Z*alpha' + rep_matrix(trt * beta_trt, M);")
     }
 
     ### Add priors for relevant parameters
@@ -187,8 +187,8 @@ setMethod(
       tau ~ {{tau_prior}} ;
       real sigma;
       sigma = 1 / tau;
-      alpha_ext ~ {{alpha_2_prior}};
-      alpha ~ normal(alpha_ext, sqrt(sigma)) ;
+      alpha[, 2] ~ {{alpha_2_prior}};
+      alpha[, 1] ~ normal(alpha[,2], sqrt(sigma)) ;
       ")
 
     ### Add in likelihood function
@@ -209,50 +209,45 @@ setMethod(
 )
 
 # BorrowingNone, OutcomeSurvPiecewiseExponential -----
-#' @rdname make_model_string_model
-setMethod(
-  "make_model_string_model",
-  signature("BorrowingNone", "OutcomeSurvPiecewiseExponential", "Analysis"),
-  function(borrowing, outcome, analysis_obj) {
-    ### Treatment prior
-    beta_trt_prior <- get_prior_string(analysis_obj@treatment@trt_prior)
+make_model_string_pem_full_none <- function(borrowing, outcome, analysis_obj) {
+  ### Treatment prior
+  beta_trt_prior <- get_prior_string(analysis_obj@treatment@trt_prior)
 
-    ### Linear predictor
-    has_covariates <- !is.null(analysis_obj@covariates)
+  ### Linear predictor
+  has_covariates <- !is.null(analysis_obj@covariates)
 
-    linear_predictor <- if (has_covariates) {
-      h_glue("lp = alpha' + rep_matrix(trt * beta_trt + X * beta, M);")
-    } else if (!has_covariates) {
-      h_glue("lp = alpha' + rep_matrix(trt * beta_trt, M);")
-    }
+  linear_predictor <- if (has_covariates) {
+    h_glue("lp = rep_matrix(alpha', N) + rep_matrix(trt * beta_trt + X * beta, M);")
+  } else if (!has_covariates) {
+    h_glue("lp = rep_matrix(alpha', N) + rep_matrix(trt * beta_trt, M);")
+  }
 
-    ### Add priors for relevant parameters
-    if (NROW(analysis_obj@outcome@param_priors) > 0) {
-      names <- names(analysis_obj@outcome@param_priors)
-      values <- get_prior_string(analysis_obj@outcome@param_priors)
-      outcome_prior <- h_glue("{{names}} ~ {{values}} ;", collapse = TRUE)
-    } else {
-      outcome_prior <- ""
-    }
+  ### Add priors for relevant parameters
+  if (NROW(analysis_obj@outcome@param_priors) > 0) {
+    names <- names(analysis_obj@outcome@param_priors)
+    values <- get_prior_string(analysis_obj@outcome@param_priors)
+    outcome_prior <- h_glue("{{names}} ~ {{values}} ;", collapse = TRUE)
+  } else {
+    outcome_prior <- ""
+  }
 
-    ### Add priors on betas
-    if (has_covariates) {
-      i <- seq_along(analysis_obj@covariates@covariates)
-      value <- get_prior_string(analysis_obj@covariates@priors)
-      index <- if (test_named(value)) get_vars(analysis_obj@covariates) else rep(1, length(i))
-      covariate_prior <- h_glue("beta[{{i}}] ~ {{value[index]}} ;", collapse = TRUE)
-    } else {
-      covariate_prior <- ""
-    }
+  ### Add priors on betas
+  if (has_covariates) {
+    i <- seq_along(analysis_obj@covariates@covariates)
+    value <- get_prior_string(analysis_obj@covariates@priors)
+    index <- if (test_named(value)) get_vars(analysis_obj@covariates) else rep(1, length(i))
+    covariate_prior <- h_glue("beta[{{i}}] ~ {{value[index]}} ;", collapse = TRUE)
+  } else {
+    covariate_prior <- ""
+  }
 
-    tau_prior <- get_prior_string(analysis_obj@borrowing@tau_prior)
-    alpha_prior <- get_prior_string(analysis_obj@outcome@baseline_prior)
-    borrowing_string <- h_glue("alpha ~ {{alpha_prior}};")
+  alpha_prior <- get_prior_string(analysis_obj@outcome@baseline_prior)
+  borrowing_string <- h_glue("alpha ~ {{alpha_prior}};")
 
-    ### Add in likelihood function
-    likelihood_string <- "target += sum((lp .* D) - (exp(lp) .* T));"
+  ### Add in likelihood function
+  likelihood_string <- "target += sum((lp .* D) - (exp(lp) .* T));"
 
-    h_glue("
+  h_glue("
   model {
     matrix[N,M] lp;
     beta_trt ~ {{beta_trt_prior}};
@@ -262,5 +257,18 @@ setMethod(
     {{borrowing_string}}
     {{likelihood_string}}
   }")
-  }
+}
+
+#' @rdname make_model_string_model
+setMethod(
+  "make_model_string_model",
+  signature("BorrowingNone", "OutcomeSurvPiecewiseExponential", "Analysis"),
+  make_model_string_pem_full_none
+)
+
+#' @rdname make_model_string_model
+setMethod(
+  "make_model_string_model",
+  signature("BorrowingFull", "OutcomeSurvPiecewiseExponential", "Analysis"),
+  make_model_string_pem_full_none
 )

R/make_model_string_parameters.R

@@ -14,14 +14,15 @@
 #'   borrowing = borrowing_full("ext"),
 #'   treatment = treatment_details("trt", prior_normal(0, 100))
 #' )
-#' 
+#'
 #' make_model_string_parameters(anls_obj)
 #' @noRd
 make_model_string_parameters <- function(analysis_obj) {
   ## Parameters string
   trt_string <- h_glue("real{{eval_constraints(analysis_obj@treatment@trt_prior)}} beta_trt;")
 
   is_bdb <- isTRUE(is(analysis_obj@borrowing, "BorrowingHierarchicalCommensurate"))
+  is_pem <- isTRUE(is(analysis_obj@outcome, "OutcomeSurvPiecewiseExponential"))
   ### Set tau
   borrowing_string <- if (is_bdb) h_glue("real{{eval_constraints(analysis_obj@borrowing@tau_prior)}} tau;") else ""
 
@@ -32,6 +33,7 @@ make_model_string_parameters <- function(analysis_obj) {
     constraint = eval_constraints(analysis_obj@outcome@baseline_prior),
     n = if (is_bdb) "[2]" else ""
   )
+  if (is_pem) intercept_string <- ""
 
   ### Add outcome specific parameters
   if (NROW(analysis_obj@outcome@param_priors) > 0) {
@@ -41,6 +43,11 @@ make_model_string_parameters <- function(analysis_obj) {
     outcome_string <- analysis_obj@outcome@param_stan_code
   }
 
+  if (is_pem && is_bdb) {
+    intercept_string <- "matrix[M, 2] alpha;"
+    outcome_string <- ""
+  }
+
   ### Add in vector of coefficients if covariates are provided
   covariate_string <- if (!is.null(analysis_obj@covariates)) "vector<lower=L_beta, upper=U_beta>[K] beta;" else ""
 

R/outcome_surv_piecewise_exponential.R

@@ -74,8 +74,6 @@ outcome_surv_piecewise_exp <- function(time_var, cens_var, start_times, baseline
       "
       vector[N] time;
       vector[N] cens;
-      int<lower = 1> M = {{length(start_times)}};
-      vector[M] starts = [{{toString(start_times)}}]';
       {{weight}}",
       weight = if (has_weight) "vector[N] weight;" else ""
     ),
@@ -86,6 +84,8 @@ outcome_surv_piecewise_exp <- function(time_var, cens_var, start_times, baseline
     }"),
     param_stan_code = "vector[M] alpha;",
     transformed_data_stan_code = h_glue("
+      int<lower = 1> M = {{length(start_times)}};
+      vector[M] starts = [{{toString(start_times)}}]';
       vector[M] durations;
       matrix[N,M] T;
       matrix[N,M] D;