diff --git a/R/DataJoint.R b/R/DataJoint.R index 640c754c..943932fa 100755 --- a/R/DataJoint.R +++ b/R/DataJoint.R @@ -186,8 +186,8 @@ as.list.DataJoint <- function(x, ...) { subset.DataJoint <- function(x, patients, ...) { data <- as.list(x) dat <- data.frame( - time = data[["Times"]], - event = as.numeric(seq_along(data[["Times"]]) %in% data[["dead_ind_index"]]), + time = data[["event_times"]], + event = as.numeric(seq_along(data[["event_times"]]) %in% data[["subject_event_index"]]), patient = names(data[["pt_to_ind"]]) ) subset_and_add_grouping(dat, patients) @@ -246,9 +246,9 @@ extract_observed_values <- function(object) { assert_class(object, "DataJoint") data <- as.list(object) x <- data.frame( - subject = names(data$pt_to_ind)[data$ind_index], - time = data$Tobs, - Yob = data$Yobs + subject = names(data$pt_to_ind)[data$subject_tumour_index], + time = data$tumour_time, + Yob = data$tumour_value ) row.names(x) <- NULL x diff --git a/R/DataLongitudinal.R b/R/DataLongitudinal.R index 13d79913..a4e16ce1 100644 --- a/R/DataLongitudinal.R +++ b/R/DataLongitudinal.R @@ -188,20 +188,20 @@ as_stan_list.DataLongitudinal <- function(object, subject_var, ...) { sparse_mat_inds_cens_y <- rstan::extract_sparse_parts(mat_sld_index[, index_cen]) model_data <- list( - Nta_total = nrow(df), + n_tumour_all = nrow(df), # Number of individuals and tumour assessments. - Nta_obs_y = length(index_obs), - Nta_cens_y = length(index_cen), + n_tumour_obs = length(index_obs), + n_tumour_cens = length(index_cen), # Index vectors - ind_index = as.numeric(df[[subject_var]]), - obs_y_index = index_obs, - cens_y_index = index_cen, + subject_tumour_index = as.numeric(df[[subject_var]]), + subject_tumour_index_obs = index_obs, + subject_tumour_index_cens = index_cen, - Yobs = df[[vars$outcome]], - Tobs = df[[vars$time]], - Ythreshold = adj_threshold, + tumour_value = df[[vars$outcome]], + tumour_time = df[[vars$time]], + tumour_value_lloq = adj_threshold, # Sparse matrix parameters # Matrix of individuals x observed tumour assessments. diff --git a/R/DataSubject.R b/R/DataSubject.R index 3405c35e..3d83df4c 100644 --- a/R/DataSubject.R +++ b/R/DataSubject.R @@ -106,11 +106,11 @@ as_stan_list.DataSubject <- function(object, ...) { df <- as.data.frame(harmonise(object)) vars <- extractVariableNames(object) list( - Nind = nrow(df), + n_subjects = nrow(df), n_studies = length(unique(df[[vars$study]])), n_arms = length(unique(df[[vars$arm]])), - pt_study_index = as.numeric(df[[vars$study]]), - pt_arm_index = as.numeric(df[[vars$arm]]), + subject_study_index = as.numeric(df[[vars$study]]), + subject_arm_index = as.numeric(df[[vars$arm]]), pt_to_ind = stats::setNames( seq_len(nlevels(df[[vars$subject]])), levels(df[[vars$subject]]) diff --git a/R/DataSurvival.R b/R/DataSurvival.R index bc09a26c..82850c44 100644 --- a/R/DataSurvival.R +++ b/R/DataSurvival.R @@ -131,9 +131,9 @@ as_stan_list.DataSurvival <- function(object, ...) { ) model_data <- list( - Nind_dead = sum(df[[vars$event]]), - dead_ind_index = which(df[[vars$event]] == 1), - Times = df[[vars$time]], + n_subject_event = sum(df[[vars$event]]), + subject_event_index = which(df[[vars$event]] == 1), + event_times = df[[vars$time]], p_os_cov_design = ncol(design_mat), os_cov_design = design_mat, n_nodes = length(gh_parameters$nodes), diff --git a/R/LongitudinalGSF.R b/R/LongitudinalGSF.R index 7bd672d2..336f74be 100755 --- a/R/LongitudinalGSF.R +++ b/R/LongitudinalGSF.R @@ -76,7 +76,7 @@ LongitudinalGSF <- function( Parameter( name = "lm_gsf_psi_phi", prior = prior_init_only(prior_beta(a_phi@init, b_phi@init)), - size = "Nind" + size = "n_subjects" ), Parameter(name = "lm_gsf_sigma", prior = sigma, size = 1) @@ -88,24 +88,24 @@ LongitudinalGSF <- function( Parameter( name = "lm_gsf_psi_bsld", prior = prior_init_only(prior_lognormal(mu_bsld@init, omega_bsld@init)), - size = "Nind" + size = "n_subjects" ), Parameter( name = "lm_gsf_psi_ks", prior = prior_init_only(prior_lognormal(mu_ks@init, omega_ks@init)), - size = "Nind" + size = "n_subjects" ), Parameter( name = "lm_gsf_psi_kg", prior = prior_init_only(prior_lognormal(mu_kg@init, omega_kg@init)), - size = "Nind" + size = "n_subjects" ) ) } else { list( - Parameter(name = "lm_gsf_eta_tilde_bsld", prior = prior_std_normal(), size = "Nind"), - Parameter(name = "lm_gsf_eta_tilde_ks", prior = prior_std_normal(), size = "Nind"), - Parameter(name = "lm_gsf_eta_tilde_kg", prior = prior_std_normal(), size = "Nind") + Parameter(name = "lm_gsf_eta_tilde_bsld", prior = prior_std_normal(), size = "n_subjects"), + Parameter(name = "lm_gsf_eta_tilde_ks", prior = prior_std_normal(), size = "n_subjects"), + Parameter(name = "lm_gsf_eta_tilde_kg", prior = prior_std_normal(), size = "n_subjects") ) } parameters <- append(parameters, parameters_extra) diff --git a/R/LongitudinalQuantities.R b/R/LongitudinalQuantities.R index 3bb8fe74..73ab7ab6 100644 --- a/R/LongitudinalQuantities.R +++ b/R/LongitudinalQuantities.R @@ -61,7 +61,7 @@ LongitudinalQuantities <- function( assert_character(groups, null.ok = TRUE) data <- as.list(object@data) patients <- decompose_patients(groups, names(data$pt_to_ind)) - time_grid <- expand_time_grid(time_grid, max(data[["Tobs"]])) + time_grid <- expand_time_grid(time_grid, max(data[["tumour_time"]])) gq <- generateQuantities( object, diff --git a/R/LongitudinalRandomSlope.R b/R/LongitudinalRandomSlope.R index 8a500ec1..3a883815 100755 --- a/R/LongitudinalRandomSlope.R +++ b/R/LongitudinalRandomSlope.R @@ -49,7 +49,7 @@ LongitudinalRandomSlope <- function( Parameter( name = "lm_rs_ind_rnd_slope", prior = prior_init_only(prior_normal(slope_mu@init, slope_sigma@init)), - size = "Nind" + size = "n_subjects" ) ) ) diff --git a/R/LongitudinalSteinFojo.R b/R/LongitudinalSteinFojo.R index 516a2212..9cf44fa5 100755 --- a/R/LongitudinalSteinFojo.R +++ b/R/LongitudinalSteinFojo.R @@ -72,24 +72,24 @@ LongitudinalSteinFojo <- function( Parameter( name = "lm_sf_psi_bsld", prior = prior_init_only(prior_lognormal(mu_bsld@init, omega_bsld@init)), - size = "Nind" + size = "n_subjects" ), Parameter( name = "lm_sf_psi_ks", prior = prior_init_only(prior_lognormal(mu_ks@init, omega_ks@init)), - size = "Nind" + size = "n_subjects" ), Parameter( name = "lm_sf_psi_kg", prior = prior_init_only(prior_lognormal(mu_kg@init, omega_kg@init)), - size = "Nind" + size = "n_subjects" ) ) } else { list( - Parameter(name = "lm_sf_eta_tilde_bsld", prior = prior_std_normal(), size = "Nind"), - Parameter(name = "lm_sf_eta_tilde_ks", prior = prior_std_normal(), size = "Nind"), - Parameter(name = "lm_sf_eta_tilde_kg", prior = prior_std_normal(), size = "Nind") + Parameter(name = "lm_sf_eta_tilde_bsld", prior = prior_std_normal(), size = "n_subjects"), + Parameter(name = "lm_sf_eta_tilde_ks", prior = prior_std_normal(), size = "n_subjects"), + Parameter(name = "lm_sf_eta_tilde_kg", prior = prior_std_normal(), size = "n_subjects") ) } parameters <- append(parameters, parameters_extra) diff --git a/R/SurvivalQuantities.R b/R/SurvivalQuantities.R index ac43b706..49fc83fe 100644 --- a/R/SurvivalQuantities.R +++ b/R/SurvivalQuantities.R @@ -78,7 +78,7 @@ SurvivalQuantities <- function( data <- as.list(object@data) patients <- decompose_patients(groups, names(data$pt_to_ind)) - time_grid <- expand_time_grid(time_grid, max(data[["Times"]])) + time_grid <- expand_time_grid(time_grid, max(data[["event_times"]])) assert_that( all(time_grid >= 0), diff --git a/inst/WORDLIST b/inst/WORDLIST index 56bbf2f6..dfb3f543 100644 --- a/inst/WORDLIST +++ b/inst/WORDLIST @@ -8,7 +8,7 @@ Kaplan Kerioui LogLogistic Modelling -Nind +n_subjects Overal Pharm postprocessing @@ -136,3 +136,4 @@ u SimLongitudinal SimSurvival geq +LLoQ diff --git a/inst/stan/base/base.stan b/inst/stan/base/base.stan index 6d6be06e..54eacf11 100755 --- a/inst/stan/base/base.stan +++ b/inst/stan/base/base.stan @@ -11,11 +11,11 @@ data{ // // Source - base/base.stan // - int Nind; // Number of individuals. + int n_subjects; // Number of individuals. int n_studies; // Number of studies. int n_arms; // Number of treatment arms. - array[Nind] int pt_study_index; // Index of study per pt (PT index sorted) - array[Nind] int pt_arm_index; // Index of treatment arm per pt (PT index sorted) + array[n_subjects] int subject_study_index; // Index of study per pt (PT index sorted) + array[n_subjects] int subject_arm_index; // Index of treatment arm per pt (PT index sorted) {{ survival.data }} {{ link.data }} @@ -48,7 +48,7 @@ transformed parameters{ // // Log-likelihood values for using the loo package. - vector[Nind] log_lik = rep_vector(0.0, Nind); + vector[n_subjects] log_lik = rep_vector(0.0, n_subjects); {{ longitudinal.transformed_parameters }} {{ link.transformed_parameters }} diff --git a/inst/stan/base/generated_quantities_data.stan b/inst/stan/base/generated_quantities_data.stan index 00febed0..0cb5aac9 100644 --- a/inst/stan/base/generated_quantities_data.stan +++ b/inst/stan/base/generated_quantities_data.stan @@ -8,5 +8,5 @@ data { int n_sm_time_grid; // Number of time points in the grid (os) vector[n_sm_time_grid] sm_time_grid; // Time points grid. int n_pt_select_index; - array[n_pt_select_index] int pt_select_index; // Index vector of which patients we want to return + array[n_pt_select_index] int pt_select_index; // Index vector of which patients we want to return } diff --git a/inst/stan/base/link_none.stan b/inst/stan/base/link_none.stan index c373301c..33e286bd 100644 --- a/inst/stan/base/link_none.stan +++ b/inst/stan/base/link_none.stan @@ -21,7 +21,7 @@ transformed data { // // If user has requested link_none then provide a dummy pars_lm object // that contains nothing - matrix[Nind, 0] link_function_inputs = rep_matrix(0, Nind, 0); + matrix[n_subjects, 0] link_function_inputs = rep_matrix(0, n_subjects, 0); vector[0] link_coefficients; } diff --git a/inst/stan/base/longitudinal.stan b/inst/stan/base/longitudinal.stan index 6ce86c87..814dbc59 100755 --- a/inst/stan/base/longitudinal.stan +++ b/inst/stan/base/longitudinal.stan @@ -10,34 +10,34 @@ data{ // // Longitudinal data - int Nta_total; // Total number of tumour assessments. - int Nta_obs_y; // Number of observed tumour assessments (not censored). - int Nta_cens_y; // Number of censored tumour assessments (below threshold). + int n_tumour_all; // Total number of tumour assessments. + int n_tumour_obs; // Number of observed tumour assessments (not censored). + int n_tumour_cens; // Number of censored tumour assessments (below threshold). - array[Nta_total] int ind_index; // Index of individuals for each tumour assessment. - array[Nta_obs_y] int obs_y_index; // Index of observed tumour assessments (not censored). - array[Nta_cens_y] int cens_y_index; // Index of censored tumour assessments. + array[n_tumour_all] int subject_tumour_index; // Index of individuals for each tumour assessment. + array[n_tumour_obs] int subject_tumour_index_obs; // Index of observed tumour assessments (not censored). + array[n_tumour_cens] int subject_tumour_index_cens; // Index of censored tumour assessments. - vector[Nta_total] Yobs; // Array of individual responses. - vector[Nta_total] Tobs; // Individual timepoints. - real Ythreshold; // Censoring threshold. + vector[n_tumour_all] tumour_value; // Array of individual responses. + vector[n_tumour_all] tumour_time; // Individual timepoints. + real tumour_value_lloq; // Censoring threshold. // Matrix of individuals x observed tumour assessments (sparse matrix of 0s and 1s), - // so the dimension is Nind x Nta_obs_y. + // so the dimension is n_subjects x n_tumour_obs. array [3] int n_mat_inds_obs_y; vector[n_mat_inds_obs_y[1]] w_mat_inds_obs_y; array[n_mat_inds_obs_y[2]] int v_mat_inds_obs_y; array[n_mat_inds_obs_y[3]] int u_mat_inds_obs_y; // Matrix of individuals x censored tumour assessments (sparse matrix of 0s and 1s). - // so the dimension is Nind x Nta_cens_y. + // so the dimension is n_subjects x n_tumour_cens. array [3] int n_mat_inds_cens_y; vector[n_mat_inds_cens_y[1]] w_mat_inds_cens_y; array[n_mat_inds_cens_y[2]] int v_mat_inds_cens_y; array[n_mat_inds_cens_y[3]] int u_mat_inds_cens_y; // Matrix of all individuals x tumour assessments (sparse matrix of 0s and 1s). - // so the dimension is Nind x Nta_total. + // so the dimension is n_subjects x n_tumour_all. array [3] int n_mat_inds_all_y; vector[n_mat_inds_all_y[1]] w_mat_inds_all_y; array[n_mat_inds_all_y[2]] int v_mat_inds_all_y; @@ -58,7 +58,7 @@ transformed parameters { // // Source - base/longitudinal.stan // - vector[Nta_total] Ypred_log_lik = rep_vector(0, Nta_total); + vector[n_tumour_all] Ypred_log_lik = rep_vector(0, n_tumour_all); {{ stan.transformed_parameters }} @@ -66,8 +66,8 @@ transformed parameters { // Source - base/longitudinal.stan // log_lik += csr_matrix_times_vector( - Nind, - Nta_total, + n_subjects, + n_tumour_all, w_mat_inds_all_y, v_mat_inds_all_y, u_mat_inds_all_y, diff --git a/inst/stan/base/survival.stan b/inst/stan/base/survival.stan index 35d3c3c2..737c79f1 100755 --- a/inst/stan/base/survival.stan +++ b/inst/stan/base/survival.stan @@ -152,11 +152,11 @@ data{ // Source - base/survival.stan // - int Nind_dead; // Number of dead individuals (observed survival time). - array[Nind_dead] int dead_ind_index; // Index of dead individuals (observed survival time). - vector[Nind] Times; + int n_subject_event; // Number of dead individuals (observed survival time). + array[n_subject_event] int subject_event_index; // Index of dead individuals (observed survival time). + vector[n_subjects] event_times; int p_os_cov_design; - matrix[Nind, p_os_cov_design] os_cov_design; + matrix[n_subjects, p_os_cov_design] os_cov_design; // Integration parameters ---- // These are the x positions and weights required to evaluate a polynomial function @@ -168,9 +168,9 @@ data{ transformed data { - array[rows(Times)] int time_positive = is_positive(Times); - int n_positive = sum(time_positive); - array[n_positive] int time_positive_index = which(time_positive); + array[rows(event_times)] int time_positive_flag = is_positive(event_times); + int n_times_positive = sum(time_positive_flag); + array[n_times_positive] int time_positive_index = which(time_positive_flag); {{ stan.transformed_data }} } @@ -199,7 +199,7 @@ transformed parameters { // // Calculate covariate contributions to log hazard function - vector[Nind] os_cov_contribution = get_os_cov_contribution( + vector[n_subjects] os_cov_contribution = get_os_cov_contribution( os_cov_design, beta_os_cov ); @@ -211,10 +211,10 @@ transformed parameters { // // Log of survival function at the observed time points. - vector[Nind] log_surv_fit_at_obs_times; - log_surv_fit_at_obs_times = rep_vector(0.0, Nind); + vector[n_subjects] log_surv_fit_at_obs_times; + log_surv_fit_at_obs_times = rep_vector(0.0, n_subjects); log_surv_fit_at_obs_times[time_positive_index] += log_survival( - Times[time_positive_index], + event_times[time_positive_index], pars_os, link_function_inputs[time_positive_index], link_coefficients, @@ -227,13 +227,13 @@ transformed parameters { log_lik += log_surv_fit_at_obs_times; // In case of death we add the log-hazard on top. - log_lik[dead_ind_index] += to_vector( + log_lik[subject_event_index] += to_vector( log_hazard( - to_matrix(Times[dead_ind_index]), + to_matrix(event_times[subject_event_index]), pars_os, - link_function_inputs[dead_ind_index], + link_function_inputs[subject_event_index], link_coefficients, - os_cov_contribution[dead_ind_index] + os_cov_contribution[subject_event_index] ) ); } diff --git a/inst/stan/lm-gsf/link.stan b/inst/stan/lm-gsf/link.stan index d70de732..acd916bf 100755 --- a/inst/stan/lm-gsf/link.stan +++ b/inst/stan/lm-gsf/link.stan @@ -5,7 +5,7 @@ transformed parameters { // // Source - lm-gsf/link.stan // - matrix[Nind, 4] link_function_inputs; + matrix[n_subjects, 4] link_function_inputs; link_function_inputs[,1] = lm_gsf_psi_bsld; link_function_inputs[,2] = lm_gsf_psi_ks; link_function_inputs[,3] = lm_gsf_psi_kg; diff --git a/inst/stan/lm-gsf/model.stan b/inst/stan/lm-gsf/model.stan index 1d0e8bc6..37d0626d 100755 --- a/inst/stan/lm-gsf/model.stan +++ b/inst/stan/lm-gsf/model.stan @@ -16,17 +16,17 @@ parameters{ real lm_gsf_omega_kg; {% if centred -%} - vector[Nind] lm_gsf_psi_bsld; - vector[Nind] lm_gsf_psi_ks; - vector[Nind] lm_gsf_psi_kg; + vector[n_subjects] lm_gsf_psi_bsld; + vector[n_subjects] lm_gsf_psi_ks; + vector[n_subjects] lm_gsf_psi_kg; {% else -%} - vector[Nind] lm_gsf_eta_tilde_bsld; - vector[Nind] lm_gsf_eta_tilde_ks; - vector[Nind] lm_gsf_eta_tilde_kg; + vector[n_subjects] lm_gsf_eta_tilde_bsld; + vector[n_subjects] lm_gsf_eta_tilde_ks; + vector[n_subjects] lm_gsf_eta_tilde_kg; {%- endif -%} // Phi Parameters - vector[Nind] lm_gsf_psi_phi; + vector[n_subjects] lm_gsf_psi_phi; vector[n_arms] lm_gsf_a_phi; vector[n_arms] lm_gsf_b_phi; @@ -45,38 +45,38 @@ transformed parameters{ // {% if not centred -%} - vector[Nind] lm_gsf_psi_bsld = exp( - lm_gsf_mu_bsld[pt_study_index] + (lm_gsf_eta_tilde_bsld * lm_gsf_omega_bsld) + vector[n_subjects] lm_gsf_psi_bsld = exp( + lm_gsf_mu_bsld[subject_study_index] + (lm_gsf_eta_tilde_bsld * lm_gsf_omega_bsld) ); - vector[Nind] lm_gsf_psi_ks = exp( - lm_gsf_mu_ks[pt_arm_index] + (lm_gsf_eta_tilde_ks * lm_gsf_omega_ks) + vector[n_subjects] lm_gsf_psi_ks = exp( + lm_gsf_mu_ks[subject_arm_index] + (lm_gsf_eta_tilde_ks * lm_gsf_omega_ks) ); - vector[Nind] lm_gsf_psi_kg = exp( - lm_gsf_mu_kg[pt_arm_index] + (lm_gsf_eta_tilde_kg * lm_gsf_omega_kg) + vector[n_subjects] lm_gsf_psi_kg = exp( + lm_gsf_mu_kg[subject_arm_index] + (lm_gsf_eta_tilde_kg * lm_gsf_omega_kg) ); {%- endif -%} - vector[Nta_total] Ypred; + vector[n_tumour_all] Ypred; Ypred = sld( - Tobs, - lm_gsf_psi_bsld[ind_index], - lm_gsf_psi_ks[ind_index], - lm_gsf_psi_kg[ind_index], - lm_gsf_psi_phi[ind_index] + tumour_time, + lm_gsf_psi_bsld[subject_tumour_index], + lm_gsf_psi_ks[subject_tumour_index], + lm_gsf_psi_kg[subject_tumour_index], + lm_gsf_psi_phi[subject_tumour_index] ); - Ypred_log_lik[obs_y_index] = vect_normal_log_dens( - Yobs[obs_y_index], - Ypred[obs_y_index], - Ypred[obs_y_index] * lm_gsf_sigma + Ypred_log_lik[subject_tumour_index_obs] = vect_normal_log_dens( + tumour_value[subject_tumour_index_obs], + Ypred[subject_tumour_index_obs], + Ypred[subject_tumour_index_obs] * lm_gsf_sigma ); - if (Nta_cens_y > 0 ) { - Ypred_log_lik[cens_y_index] = vect_normal_log_cum( - Ythreshold, - Ypred[cens_y_index], - Ypred[cens_y_index] * lm_gsf_sigma + if (n_tumour_cens > 0 ) { + Ypred_log_lik[subject_tumour_index_cens] = vect_normal_log_cum( + tumour_value_lloq, + Ypred[subject_tumour_index_cens], + Ypred[subject_tumour_index_cens] * lm_gsf_sigma ); } } @@ -87,11 +87,11 @@ model { // Source - lm-gsf/model.stan // {% if centred %} - lm_gsf_psi_bsld ~ lognormal(lm_gsf_mu_bsld[pt_study_index], lm_gsf_omega_bsld); - lm_gsf_psi_ks ~ lognormal(lm_gsf_mu_ks[pt_arm_index], lm_gsf_omega_ks); - lm_gsf_psi_kg ~ lognormal(lm_gsf_mu_kg[pt_arm_index], lm_gsf_omega_kg); + lm_gsf_psi_bsld ~ lognormal(lm_gsf_mu_bsld[subject_study_index], lm_gsf_omega_bsld); + lm_gsf_psi_ks ~ lognormal(lm_gsf_mu_ks[subject_arm_index], lm_gsf_omega_ks); + lm_gsf_psi_kg ~ lognormal(lm_gsf_mu_kg[subject_arm_index], lm_gsf_omega_kg); {%- endif -%} - lm_gsf_psi_phi ~ beta(lm_gsf_a_phi[pt_arm_index], lm_gsf_b_phi[pt_arm_index]); + lm_gsf_psi_phi ~ beta(lm_gsf_a_phi[subject_arm_index], lm_gsf_b_phi[subject_arm_index]); } @@ -99,7 +99,7 @@ generated quantities { // // Source - lm-gsf/model.stan // - matrix[Nind, 4] long_gq_parameters; + matrix[n_subjects, 4] long_gq_parameters; long_gq_parameters[, 1] = lm_gsf_psi_bsld; long_gq_parameters[, 2] = lm_gsf_psi_ks; long_gq_parameters[, 3] = lm_gsf_psi_kg; diff --git a/inst/stan/lm-random-slope/link.stan b/inst/stan/lm-random-slope/link.stan index 1c94e43f..e770230e 100644 --- a/inst/stan/lm-random-slope/link.stan +++ b/inst/stan/lm-random-slope/link.stan @@ -4,7 +4,7 @@ transformed parameters { // Source - lm-random-slope/link.stan // - matrix[Nind, 2] link_function_inputs; + matrix[n_subjects, 2] link_function_inputs; link_function_inputs[, 1] = lm_rs_ind_intercept; link_function_inputs[, 2] = lm_rs_ind_rnd_slope; } diff --git a/inst/stan/lm-random-slope/model.stan b/inst/stan/lm-random-slope/model.stan index 7508553c..0b3a93a6 100755 --- a/inst/stan/lm-random-slope/model.stan +++ b/inst/stan/lm-random-slope/model.stan @@ -22,7 +22,7 @@ parameters { array [n_arms] real lm_rs_slope_mu; real lm_rs_slope_sigma; real lm_rs_sigma; - vector[Nind] lm_rs_ind_rnd_slope; + vector[n_subjects] lm_rs_ind_rnd_slope; } @@ -30,21 +30,21 @@ transformed parameters { // // Source - lm-random-slope/model.stan // - vector[Nind] lm_rs_ind_intercept = to_vector(lm_rs_intercept[pt_study_index]); - vector[Nta_total] lm_rs_rslope_ind = to_vector(lm_rs_ind_rnd_slope[ind_index]); + vector[n_subjects] lm_rs_ind_intercept = to_vector(lm_rs_intercept[subject_study_index]); + vector[n_tumour_all] lm_rs_rslope_ind = to_vector(lm_rs_ind_rnd_slope[subject_tumour_index]); - vector[Nta_total] Ypred = lm_rs_ind_intercept[ind_index] + lm_rs_rslope_ind .* Tobs; + vector[n_tumour_all] Ypred = lm_rs_ind_intercept[subject_tumour_index] + lm_rs_rslope_ind .* tumour_time; - Ypred_log_lik[obs_y_index] = vect_normal_log_dens( - Yobs[obs_y_index], - Ypred[obs_y_index], - rep_vector(lm_rs_sigma, Nta_obs_y) + Ypred_log_lik[subject_tumour_index_obs] = vect_normal_log_dens( + tumour_value[subject_tumour_index_obs], + Ypred[subject_tumour_index_obs], + rep_vector(lm_rs_sigma, n_tumour_obs) ); - if (Nta_cens_y > 0 ) { - Ypred_log_lik[cens_y_index] = vect_normal_log_cum( - Ythreshold, - Ypred[cens_y_index], - rep_vector(lm_rs_sigma, Nta_cens_y) + if (n_tumour_cens > 0 ) { + Ypred_log_lik[subject_tumour_index_cens] = vect_normal_log_cum( + tumour_value_lloq, + Ypred[subject_tumour_index_cens], + rep_vector(lm_rs_sigma, n_tumour_cens) ); } } @@ -56,7 +56,7 @@ model { // lm_rs_ind_rnd_slope ~ normal( - lm_rs_slope_mu[pt_arm_index], + lm_rs_slope_mu[subject_arm_index], lm_rs_slope_sigma ); } @@ -66,7 +66,7 @@ generated quantities { // // Source - lm-random-slope/model.stan // - matrix[Nind, 2] long_gq_parameters; + matrix[n_subjects, 2] long_gq_parameters; long_gq_parameters[, 1] = lm_rs_ind_intercept; long_gq_parameters[, 2] = lm_rs_ind_rnd_slope; } diff --git a/inst/stan/lm-stein-fojo/link.stan b/inst/stan/lm-stein-fojo/link.stan index 698edbe2..7dab3335 100755 --- a/inst/stan/lm-stein-fojo/link.stan +++ b/inst/stan/lm-stein-fojo/link.stan @@ -5,7 +5,7 @@ transformed parameters { // // Source - lm-stein-fojo/link.stan // - matrix[Nind, 3] link_function_inputs; + matrix[n_subjects, 3] link_function_inputs; link_function_inputs[,1] = lm_sf_psi_bsld; link_function_inputs[,2] = lm_sf_psi_ks; link_function_inputs[,3] = lm_sf_psi_kg; diff --git a/inst/stan/lm-stein-fojo/model.stan b/inst/stan/lm-stein-fojo/model.stan index da8b11a9..ff0a7ce9 100755 --- a/inst/stan/lm-stein-fojo/model.stan +++ b/inst/stan/lm-stein-fojo/model.stan @@ -16,13 +16,13 @@ parameters{ real lm_sf_omega_kg; {% if centred -%} - vector[Nind] lm_sf_psi_bsld; - vector[Nind] lm_sf_psi_ks; - vector[Nind] lm_sf_psi_kg; + vector[n_subjects] lm_sf_psi_bsld; + vector[n_subjects] lm_sf_psi_ks; + vector[n_subjects] lm_sf_psi_kg; {% else -%} - vector[Nind] lm_sf_eta_tilde_bsld; - vector[Nind] lm_sf_eta_tilde_ks; - vector[Nind] lm_sf_eta_tilde_kg; + vector[n_subjects] lm_sf_eta_tilde_bsld; + vector[n_subjects] lm_sf_eta_tilde_ks; + vector[n_subjects] lm_sf_eta_tilde_kg; {%- endif -%} // Standard deviation of the error term @@ -40,36 +40,36 @@ transformed parameters{ // {% if not centred -%} - vector[Nind] lm_sf_psi_bsld = exp( - lm_sf_mu_bsld[pt_study_index] + (lm_sf_eta_tilde_bsld * lm_sf_omega_bsld) + vector[n_subjects] lm_sf_psi_bsld = exp( + lm_sf_mu_bsld[subject_study_index] + (lm_sf_eta_tilde_bsld * lm_sf_omega_bsld) ); - vector[Nind] lm_sf_psi_ks = exp( - lm_sf_mu_ks[pt_arm_index] + (lm_sf_eta_tilde_ks * lm_sf_omega_ks) + vector[n_subjects] lm_sf_psi_ks = exp( + lm_sf_mu_ks[subject_arm_index] + (lm_sf_eta_tilde_ks * lm_sf_omega_ks) ); - vector[Nind] lm_sf_psi_kg = exp( - lm_sf_mu_kg[pt_arm_index] + (lm_sf_eta_tilde_kg * lm_sf_omega_kg) + vector[n_subjects] lm_sf_psi_kg = exp( + lm_sf_mu_kg[subject_arm_index] + (lm_sf_eta_tilde_kg * lm_sf_omega_kg) ); {%- endif -%} - vector[Nta_total] Ypred; + vector[n_tumour_all] Ypred; Ypred = sld( - Tobs, - lm_sf_psi_bsld[ind_index], - lm_sf_psi_ks[ind_index], - lm_sf_psi_kg[ind_index] + tumour_time, + lm_sf_psi_bsld[subject_tumour_index], + lm_sf_psi_ks[subject_tumour_index], + lm_sf_psi_kg[subject_tumour_index] ); - Ypred_log_lik[obs_y_index] = vect_normal_log_dens( - Yobs[obs_y_index], - Ypred[obs_y_index], - Ypred[obs_y_index] * lm_sf_sigma + Ypred_log_lik[subject_tumour_index_obs] = vect_normal_log_dens( + tumour_value[subject_tumour_index_obs], + Ypred[subject_tumour_index_obs], + Ypred[subject_tumour_index_obs] * lm_sf_sigma ); - if (Nta_cens_y > 0 ) { - Ypred_log_lik[cens_y_index] = vect_normal_log_cum( - Ythreshold, - Ypred[cens_y_index], - Ypred[cens_y_index] * lm_sf_sigma + if (n_tumour_cens > 0 ) { + Ypred_log_lik[subject_tumour_index_cens] = vect_normal_log_cum( + tumour_value_lloq, + Ypred[subject_tumour_index_cens], + Ypred[subject_tumour_index_cens] * lm_sf_sigma ); } } @@ -80,9 +80,9 @@ model { // Source - lm-stein-fojo/model.stan // {% if centred %} - lm_sf_psi_bsld ~ lognormal(lm_sf_mu_bsld[pt_study_index], lm_sf_omega_bsld); - lm_sf_psi_ks ~ lognormal(lm_sf_mu_ks[pt_arm_index], lm_sf_omega_ks); - lm_sf_psi_kg ~ lognormal(lm_sf_mu_kg[pt_arm_index], lm_sf_omega_kg); + lm_sf_psi_bsld ~ lognormal(lm_sf_mu_bsld[subject_study_index], lm_sf_omega_bsld); + lm_sf_psi_ks ~ lognormal(lm_sf_mu_ks[subject_arm_index], lm_sf_omega_ks); + lm_sf_psi_kg ~ lognormal(lm_sf_mu_kg[subject_arm_index], lm_sf_omega_kg); {%- endif -%} } @@ -90,7 +90,7 @@ generated quantities { // // Source - lm-stein-fojo/model.stan // - matrix[Nind, 3] long_gq_parameters; + matrix[n_subjects, 3] long_gq_parameters; long_gq_parameters[, 1] = lm_gsf_psi_bsld; long_gq_parameters[, 2] = lm_gsf_psi_ks; long_gq_parameters[, 3] = lm_gsf_psi_kg; diff --git a/tests/testthat/test-DataJoint.R b/tests/testthat/test-DataJoint.R index 297f0075..8312d87d 100644 --- a/tests/testthat/test-DataJoint.R +++ b/tests/testthat/test-DataJoint.R @@ -36,21 +36,21 @@ test_that("DataJoint basic usage", { ) li <- as.list(d_joint) - expect_equal(li$Nind, 3) - expect_equal(li$Nta_total, 6) + expect_equal(li$n_subjects, 3) + expect_equal(li$n_tumour_all, 6) expect_equal(li$pt_to_ind, c("C" = 1, "B" = 2, "A" = 3)) expect_equal(li$n_arms, 3) - expect_equal(li$pt_study_index, c(2, 1, 1)) + expect_equal(li$subject_study_index, c(2, 1, 1)) expect_equal(li$n_studies, 2) - expect_equal(li$pt_arm_index, c(3, 2, 1)) - expect_equal(li$Times, c(150, 200, 100)) - expect_equal(li$ind_index, c(1, 2, 2, 2, 3, 3)) - expect_equal(li$obs_y_index, c(1:6)) - expect_equal(li$cens_y_index, integer(0)) - expect_equal(li$Yobs, c(5, 3, 4, 6, 1, 2)) - expect_equal(li$dead_ind_index, c(1, 2)) + expect_equal(li$subject_arm_index, c(3, 2, 1)) + expect_equal(li$event_times, c(150, 200, 100)) + expect_equal(li$subject_tumour_index, c(1, 2, 2, 2, 3, 3)) + expect_equal(li$subject_tumour_index_obs, c(1:6)) + expect_equal(li$subject_tumour_index_cens, integer(0)) + expect_equal(li$tumour_value, c(5, 3, 4, 6, 1, 2)) + expect_equal(li$subject_event_index, c(1, 2)) expect_equal(li$os_cov_design, matrix(c(4, 2, 5), ncol = 1), ignore_attr = TRUE) - expect_equal(li$Ythreshold, -999999) + expect_equal(li$tumour_value_lloq, -999999) }) diff --git a/tests/testthat/test-DataLongitudinal.R b/tests/testthat/test-DataLongitudinal.R index e9167e4e..f68f3135 100644 --- a/tests/testthat/test-DataLongitudinal.R +++ b/tests/testthat/test-DataLongitudinal.R @@ -1,8 +1,8 @@ test_that("DataLongitudinal being rendered to list is as expected for simple inputs", { expected_vars <- c( - "Nta_total", "Nta_obs_y", "Nta_cens_y", "ind_index", "obs_y_index", - "cens_y_index", "Yobs", "Tobs", "Ythreshold", + "n_tumour_all", "n_tumour_obs", "n_tumour_cens", "subject_tumour_index", "subject_tumour_index_obs", + "subject_tumour_index_cens", "tumour_value", "tumour_time", "tumour_value_lloq", "n_mat_inds_obs_y", "w_mat_inds_obs_y", "v_mat_inds_obs_y", "u_mat_inds_obs_y", "n_mat_inds_cens_y", "w_mat_inds_cens_y", "v_mat_inds_cens_y", "u_mat_inds_cens_y", "n_mat_inds_all_y", "w_mat_inds_all_y", "v_mat_inds_all_y", "u_mat_inds_all_y" @@ -26,15 +26,15 @@ test_that("DataLongitudinal being rendered to list is as expected for simple inp ) expect_equal(names(li), expected_vars) - expect_equal(li$Nta_total, 6) - expect_equal(li$Nta_obs_y, 3) - expect_equal(li$Nta_cens_y, 3) - expect_equal(li$ind_index, c(1, 3, 3, 1, 2, 3)) - expect_equal(li$obs_y_index, c(3, 4, 5)) - expect_equal(li$cens_y_index, c(1, 2, 6)) - expect_equal(li$Yobs, x$voutcome) - expect_equal(li$Tobs, x$vtime) - expect_equal(li$Ythreshold, 3) + expect_equal(li$n_tumour_all, 6) + expect_equal(li$n_tumour_obs, 3) + expect_equal(li$n_tumour_cens, 3) + expect_equal(li$subject_tumour_index, c(1, 3, 3, 1, 2, 3)) + expect_equal(li$subject_tumour_index_obs, c(3, 4, 5)) + expect_equal(li$subject_tumour_index_cens, c(1, 2, 6)) + expect_equal(li$tumour_value, x$voutcome) + expect_equal(li$tumour_time, x$vtime) + expect_equal(li$tumour_value_lloq, 3) }) diff --git a/tests/testthat/test-DataSubject.R b/tests/testthat/test-DataSubject.R index 1bd99c5b..6b0f9413 100644 --- a/tests/testthat/test-DataSubject.R +++ b/tests/testthat/test-DataSubject.R @@ -19,19 +19,19 @@ test_that("DataSubject works as expected", { expected_variables <- c( - "Nind", "n_studies", "n_arms", "pt_study_index", - "pt_arm_index", "pt_to_ind" + "n_subjects", "n_studies", "n_arms", "subject_study_index", + "subject_arm_index", "pt_to_ind" ) li <- as_stan_list(obj) expect_equal(names(li), expected_variables) - expect_equal(li$Nind, 4) + expect_equal(li$n_subjects, 4) expect_equal(li$n_studies, 2) expect_equal(li$n_arms, 3) - expect_equal(li$pt_study_index, c(2, 1, 1, 2)) - expect_equal(li$pt_arm_index, c(3, 2, 1, 3)) + expect_equal(li$subject_study_index, c(2, 1, 1, 2)) + expect_equal(li$subject_arm_index, c(3, 2, 1, 3)) expect_equal(li$pt_to_ind, c("C" = 1, "B" = 2, "A" = 3, "D" = 4)) }) diff --git a/tests/testthat/test-DataSurvival.R b/tests/testthat/test-DataSurvival.R index 592ba31b..c9a3c252 100644 --- a/tests/testthat/test-DataSurvival.R +++ b/tests/testthat/test-DataSurvival.R @@ -22,16 +22,16 @@ test_that("DataSurvival being rendered to list is as expected for simple inputs" expect_equal( c( - "Nind_dead", "dead_ind_index", "Times", "p_os_cov_design", + "n_subject_event", "subject_event_index", "event_times", "p_os_cov_design", "os_cov_design", "n_nodes", "nodes", "weights" ), names(res) ) - expect_equal(res$Nind_dead, 3) + expect_equal(res$n_subject_event, 3) expect_equal(res$p_os_cov_design, 3) expect_equal(res$os_cov_design, covmat) - expect_equal(res$dead_ind_index, c(1, 2, 4)) - expect_equal(res$Times, c(10, 20, 30, 25, 15)) + expect_equal(res$subject_event_index, c(1, 2, 4)) + expect_equal(res$event_times, c(10, 20, 30, 25, 15)) ## Dropped rows works as expected @@ -61,16 +61,16 @@ test_that("DataSurvival being rendered to list is as expected for simple inputs" expect_equal( c( - "Nind_dead", "dead_ind_index", "Times", "p_os_cov_design", + "n_subject_event", "subject_event_index", "event_times", "p_os_cov_design", "os_cov_design", "n_nodes", "nodes", "weights" ), names(res) ) - expect_equal(res$Nind_dead, 2) + expect_equal(res$n_subject_event, 2) expect_equal(res$p_os_cov_design, 3) expect_equal(res$os_cov_design, covmat) - expect_equal(res$dead_ind_index, c(1, 2)) - expect_equal(res$Times, c(10, 20, 30, 15)) + expect_equal(res$subject_event_index, c(1, 2)) + expect_equal(res$event_times, c(10, 20, 30, 15)) }) diff --git a/tests/testthat/test-Link.R b/tests/testthat/test-Link.R index 88bdd33b..4e631f87 100644 --- a/tests/testthat/test-Link.R +++ b/tests/testthat/test-Link.R @@ -29,7 +29,7 @@ test_that("link_none() works as expected", { ) expect_true( grepl( - "matrix[Nind, 0] link_function_inputs = rep_matrix(0, Nind, 0);", + "matrix[n_subjects, 0] link_function_inputs = rep_matrix(0, n_subjects, 0);", as.character(x), fixed = TRUE ) diff --git a/vignettes/extending-jmpost.Rmd b/vignettes/extending-jmpost.Rmd index c4153684..dfaa6aaa 100644 --- a/vignettes/extending-jmpost.Rmd +++ b/vignettes/extending-jmpost.Rmd @@ -63,7 +63,7 @@ For reference the following is roughly the implementation for the `LongitudinalG enableLink.LongitudinalGSF <- function(object, ...) { stan <- StanModule(" transformed parameters { - matrix[Nind, 4] link_function_inputs; + matrix[n_subjects, 4] link_function_inputs; link_function_inputs[,1] = lm_gsf_psi_bsld; link_function_inputs[,2] = lm_gsf_psi_ks; link_function_inputs[,3] = lm_gsf_psi_kg; @@ -196,7 +196,7 @@ respectively. Under the hood this library works by merging multiple Stan programs together into a single program. It does this by parsing the program to extract out each block independently. Unfortunately, the formal Stan parser (`stanc`) provided by the Stan team -only works with complete programs whereas most of the programs within jmpost are +only works with complete programs whereas most of the programs within `jmpost` are incomplete fragments. This package has therefore implemented its own simple parser; as a result, in order to not have to traverse the full abstract syntax tree (AST), a few addition constraints are made on how Stan programs can be formatted. @@ -232,3 +232,188 @@ parameter { real mu; // non-whitespace after opening `{` x ~ normal(mu, sigma); } // some comment // non-whitespace after closing `}` ``` + +## Stan Data Objects + +When writing your own Stan code to extend `jmpost` it is important to note that +many different data objects have already been defined in the `data` block of the +base Stan template. This section outlines the different data objects that are +made available to user. Note that some objects are only made available if +the corresponding model is used; for example death times are only available if +the user specifies a `SurvivalModel` object to `JointModel()`. + +### Global Data Objects + +**Number of unique subjects** +```stan +int n_subjects; +``` + +**Number of unique studies** +```stan +int n_studies; +``` + +**Number of unique treatment arms** +```stan +int n_arms; +``` + +**Study index for each subject** +```stan +array[n_subjects] int subject_study_index; +``` +- Note that this is sorted based upon the subject's factor level within the R `data.frame`. +For example lets say subject `"A"` has a factor level of 15 and that their corresponding +study value has a factor level of 2 then `subject_study_index[15]` will be 2. + +**Treatment arm index for each subject** +```stan + array[n_subjects] int subject_arm_index; +``` +- A mirror of `subject_study_index` but for the treatment arm. + + +### Survival Data Objects + + +**Number of events** +```stan +int n_subject_event; +``` + +**Event/Censor Times** +```stan +vector[n_subjects] event_times; +``` +- Ordered by patient factor level + +**Event Index** +```stan +array[n_subject_event] int subject_event_index; +``` +- Is the index into `event_times` to identify which times are an event. The rest are censored. + + +**Number of covariates for the survival model** +```stan +int p_os_cov_design; +``` +- Note that this does not include an intercept term, which would conflict with the baseline distribution parameters. + +**Covariate design matrix for the survival model** +```stan +matrix[n_subjects, p_os_cov_design] os_cov_design; +``` +- Note that this does not include an intercept term, which would conflict with the baseline distribution parameters. + + +**Time >0 flag** +```stan +array[rows(event_times)] int time_positive_flag; +``` + +**Number of times >0** +```stan +int n_times_positive; +``` + +**Positive time index** +```stan +array[n_times_positive] int time_positive_index +``` + + +**Gaussian Quadrature Integration Parameters** +```stan +int n_nodes; +vector[n_nodes] nodes; +vector[n_nodes] weights; +``` +- These are the nodes and weights for the Gaussian quadrature integration. + + +### Longitudinal Data Objects + +**Total number of tumour assessments** +```stan +int n_tumour_all; +``` + +**Number of tumour assessments above LLoQ (Lower Limit of Quantification)** +```stan +int n_tumour_obs; +``` + +**Number of tumour assessments below LLoQ (Lower Limit of Quantification)** +```stan +int n_tumour_cens; +``` + +**Tumour assessments values** +```stan +vector[n_tumour_all] tumour_value; +``` + +**Tumour assessments time points** +```stan +vector[n_tumour_all] tumour_time; +``` + +**LLoQ threshold** +```stan +real tumour_value_lloq; +``` + +**Individual tumour assessment index** +```stan +array[n_tumour_all] int subject_tumour_index; +``` +- That is if tumour assessment 1 belongs to the patient with factor level 3 then +`subject_tumour_index[1]` will be 3. + + +**Tumour assessment index for observations above LLoQ (Lower Limit of Quantification)** +```stan +array[n_tumour_obs] int subject_tumour_index_obs; +``` +- For example if only tumour assessments 3 and 5 were above the LLoQ then +`subject_tumour_index_obs` will be `[3, 5]`. + +**Tumour assessment index for observations below LLoQ (Lower Limit of Quantification)** +```stan +array[n_tumour_cens] int subject_tumour_index_cens; +``` +- For example if only tumour assessments 1 and 2 were below the LLoQ then +`subject_tumour_index_obs` will be `[1, 2]`. + + +**Sparse matrix components for patient indexes of the tumour assessments** +```stan +array [3] int n_mat_inds_all_y; +vector[n_mat_inds_all_y[1]] w_mat_inds_all_y; +array[n_mat_inds_all_y[2]] int v_mat_inds_all_y; +array[n_mat_inds_all_y[3]] int u_mat_inds_all_y; +``` +- This is the sparse matrix representation of the binary matrix that has one row per patient and one column per tumour assessment. That is, for row 3 of this matrix all columns that have an entry of 1 indicate that the corresponding entry in `tumour_value` belongs to the patient with factor level 3. This matrix is primarily used to calculate the sum +of log-likelihood for all tumour assessments per patient in an efficient way. +- See [the Stan CSR documentation](https://mc-stan.org/docs/functions-reference/sparse_matrix_operations.html#CSR) for more information on the sparse matrix representation. + + +**Sparse matrix components for patient indexes of the tumour assessments above the LLoQ** +```stan +array [3] int n_mat_inds_obs_y; +vector[n_mat_inds_obs_y[1]] w_mat_inds_obs_y; +array[n_mat_inds_obs_y[2]] int v_mat_inds_obs_y; +array[n_mat_inds_obs_y[3]] int u_mat_inds_obs_y; +``` +- Same as above but only for the tumour assessments above the LLoQ. + +**Sparse matrix components for patient indexes of the tumour assessments below the LLoQ** +```stan +array [3] int n_mat_inds_cens_y; +vector[n_mat_inds_cens_y[1]] w_mat_inds_cens_y; +array[n_mat_inds_cens_y[2]] int v_mat_inds_cens_y; +array[n_mat_inds_cens_y[3]] int u_mat_inds_cens_y; +``` +- Same as above but only for the tumour assessments below the LLoQ.