From c86dff9f85273705b20add4e53aaafa0cd89ff45 Mon Sep 17 00:00:00 2001 From: gowerc Date: Thu, 28 Mar 2024 17:08:09 +0000 Subject: [PATCH 1/5] added data object documentation --- vignettes/extending-jmpost.Rmd | 185 +++++++++++++++++++++++++++++++++ 1 file changed, 185 insertions(+) diff --git a/vignettes/extending-jmpost.Rmd b/vignettes/extending-jmpost.Rmd index c4153684..4c91ec46 100644 --- a/vignettes/extending-jmpost.Rmd +++ b/vignettes/extending-jmpost.Rmd @@ -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 Nind; +``` + +**Number of unique studies** +```stan +int n_studies; +``` + +**Number of unique treatment arms** +```stan +int n_arms; +``` + +**Study index for each patient** +```stan +array[Nind] int pt_study_index; +``` +- Note that this is sorted based upon the patients factor level within the R `data.frame`. +For example lets say patient `"A"` has a factor level of 15 and that their corresponding +study value has a factor level of 2 then `pt_study_index[15]` will be 2. + +**Treatment arm index for each patient** +```stan + array[Nind] int pt_arm_index; +``` +- A mirror of `pt_study_index` but for the treatment arm. + + +### Survival Data Objects + + +**Number of events** +```stan +int Nind_dead; +``` + +**Event/Censor Times** +```stan +vector[Nind] Times; +``` +- Ordered by patient factor level + +**Event Index** +```stan +array[Nind_dead] int dead_ind_index; +``` +- Is the index into `Times` to identiy 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 the intercept term. + +**Covariate design matrix for the survival model** +```stan +matrix[Nind, p_os_cov_design] os_cov_design; +``` +- Note that this does not include the intercept term. + + +**Time >0 flag** +```stan +array[rows(Times)] int time_positive; +``` + +**Number of times >0** +```stan +int n_positive; +``` + +**Positive time index** +```stan +array[n_positive] int time_positive_index +``` + + +**Gausian 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 Nta_total; +``` + +**Number of tumour assessments above LLoQ** +```stan +int Nta_obs_y; +``` + +**Number of tumour assessments below LLoQ** +```stan +int Nta_cens_y; +``` + +**Tumour assessments values** +```stan +vector[Nta_total] Yobs; +``` + +**Tumour assessments time points** +```stan +vector[Nta_total] Tobs; +``` + +**LLoQ threshold** +```stan +real Ythreshold; +``` + +**Individual tumour assessment index** +```stan +array[Nta_total] int ind_index; +``` +- That is if tumour assessment 1 belongs to the patient with factor level 3 then +`ind_index[1]` will be 3. + + +**Tumour assessment index for observations above LLoQ** +```stan +array[Nta_obs_y] int obs_y_index; +``` +- For example if only tumour assessments 3 and 5 were above the LLoQ then +`obs_y_index` will be `[3, 5]`. + +**Tumour assessment index for observations below LLoQ** +```stan +array[Nta_cens_y] int cens_y_index; +``` +- For example if only tumour assessments 1 and 2 were below the LLoQ then +`obs_y_index` 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 `Yobs` 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. From 7da69f1c0183ec9ab329bbc34747e2d476fff85c Mon Sep 17 00:00:00 2001 From: gowerc Date: Thu, 28 Mar 2024 17:31:19 +0000 Subject: [PATCH 2/5] Fix spellings --- inst/WORDLIST | 1 + vignettes/extending-jmpost.Rmd | 4 ++-- 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/inst/WORDLIST b/inst/WORDLIST index 56bbf2f6..7c1a1948 100644 --- a/inst/WORDLIST +++ b/inst/WORDLIST @@ -136,3 +136,4 @@ u SimLongitudinal SimSurvival geq +LLoQ diff --git a/vignettes/extending-jmpost.Rmd b/vignettes/extending-jmpost.Rmd index 4c91ec46..199b3ac1 100644 --- a/vignettes/extending-jmpost.Rmd +++ b/vignettes/extending-jmpost.Rmd @@ -292,7 +292,7 @@ vector[Nind] Times; ```stan array[Nind_dead] int dead_ind_index; ``` -- Is the index into `Times` to identiy which times are an event. The rest are censored. +- Is the index into `Times` to identify which times are an event. The rest are censored. **Number of covariates for the survival model** @@ -324,7 +324,7 @@ array[n_positive] int time_positive_index ``` -**Gausian Quadrature Integration Parameters** +**Gaussian Quadrature Integration Parameters** ```stan int n_nodes; vector[n_nodes] nodes; From 4ff36d6dd2fc1a0d0fe81a28a3039ea4839d0187 Mon Sep 17 00:00:00 2001 From: gowerc Date: Tue, 2 Apr 2024 17:23:31 +0100 Subject: [PATCH 3/5] backticks + n_subjects --- inst/WORDLIST | 2 +- inst/stan/base/base.stan | 8 ++++---- inst/stan/base/generated_quantities_data.stan | 2 +- inst/stan/base/link_none.stan | 2 +- inst/stan/base/longitudinal.stan | 8 ++++---- inst/stan/base/survival.stan | 8 ++++---- inst/stan/lm-gsf/model.stan | 20 +++++++++---------- inst/stan/lm-random-slope/model.stan | 4 ++-- inst/stan/lm-stein-fojo/model.stan | 18 ++++++++--------- vignettes/extending-jmpost.Rmd | 12 +++++------ 10 files changed, 42 insertions(+), 42 deletions(-) diff --git a/inst/WORDLIST b/inst/WORDLIST index 7c1a1948..dfb3f543 100644 --- a/inst/WORDLIST +++ b/inst/WORDLIST @@ -8,7 +8,7 @@ Kaplan Kerioui LogLogistic Modelling -Nind +n_subjects Overal Pharm postprocessing diff --git a/inst/stan/base/base.stan b/inst/stan/base/base.stan index 6d6be06e..4b711700 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 pt_study_index; // Index of study per pt (PT index sorted) + array[n_subjects] int pt_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..30445738 100755 --- a/inst/stan/base/longitudinal.stan +++ b/inst/stan/base/longitudinal.stan @@ -23,21 +23,21 @@ data{ real Ythreshold; // 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 Nta_obs_y. 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 Nta_cens_y. 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 Nta_total. 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; @@ -66,7 +66,7 @@ transformed parameters { // Source - base/longitudinal.stan // log_lik += csr_matrix_times_vector( - Nind, + n_subjects, Nta_total, w_mat_inds_all_y, v_mat_inds_all_y, diff --git a/inst/stan/base/survival.stan b/inst/stan/base/survival.stan index 35d3c3c2..d5e503a3 100755 --- a/inst/stan/base/survival.stan +++ b/inst/stan/base/survival.stan @@ -154,9 +154,9 @@ data{ 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; + vector[n_subjects] 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 @@ -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,7 +211,7 @@ transformed parameters { // // Log of survival function at the observed time points. - vector[Nind] log_surv_fit_at_obs_times; + vector[n_subjects] log_surv_fit_at_obs_times; log_surv_fit_at_obs_times = rep_vector(0.0, Nind); log_surv_fit_at_obs_times[time_positive_index] += log_survival( Times[time_positive_index], diff --git a/inst/stan/lm-gsf/model.stan b/inst/stan/lm-gsf/model.stan index 1d0e8bc6..9a147168 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,13 +45,13 @@ transformed parameters{ // {% if not centred -%} - vector[Nind] lm_gsf_psi_bsld = exp( + vector[n_subjects] lm_gsf_psi_bsld = exp( lm_gsf_mu_bsld[pt_study_index] + (lm_gsf_eta_tilde_bsld * lm_gsf_omega_bsld) ); - vector[Nind] lm_gsf_psi_ks = exp( + vector[n_subjects] lm_gsf_psi_ks = exp( lm_gsf_mu_ks[pt_arm_index] + (lm_gsf_eta_tilde_ks * lm_gsf_omega_ks) ); - vector[Nind] lm_gsf_psi_kg = exp( + vector[n_subjects] lm_gsf_psi_kg = exp( lm_gsf_mu_kg[pt_arm_index] + (lm_gsf_eta_tilde_kg * lm_gsf_omega_kg) ); {%- endif -%} diff --git a/inst/stan/lm-random-slope/model.stan b/inst/stan/lm-random-slope/model.stan index 7508553c..f768de3b 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,7 +30,7 @@ transformed parameters { // // Source - lm-random-slope/model.stan // - vector[Nind] lm_rs_ind_intercept = to_vector(lm_rs_intercept[pt_study_index]); + vector[n_subjects] 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[Nta_total] Ypred = lm_rs_ind_intercept[ind_index] + lm_rs_rslope_ind .* Tobs; diff --git a/inst/stan/lm-stein-fojo/model.stan b/inst/stan/lm-stein-fojo/model.stan index da8b11a9..1634ffff 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,13 +40,13 @@ transformed parameters{ // {% if not centred -%} - vector[Nind] lm_sf_psi_bsld = exp( + vector[n_subjects] lm_sf_psi_bsld = exp( lm_sf_mu_bsld[pt_study_index] + (lm_sf_eta_tilde_bsld * lm_sf_omega_bsld) ); - vector[Nind] lm_sf_psi_ks = exp( + vector[n_subjects] lm_sf_psi_ks = exp( lm_sf_mu_ks[pt_arm_index] + (lm_sf_eta_tilde_ks * lm_sf_omega_ks) ); - vector[Nind] lm_sf_psi_kg = exp( + vector[n_subjects] lm_sf_psi_kg = exp( lm_sf_mu_kg[pt_arm_index] + (lm_sf_eta_tilde_kg * lm_sf_omega_kg) ); {%- endif -%} diff --git a/vignettes/extending-jmpost.Rmd b/vignettes/extending-jmpost.Rmd index 199b3ac1..c0deec1c 100644 --- a/vignettes/extending-jmpost.Rmd +++ b/vignettes/extending-jmpost.Rmd @@ -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. @@ -235,7 +235,7 @@ parameter { real mu; // non-whitespace after opening `{` ## Stan Data Objects -When writing your own Stan code to extend jmpost it is important to note that +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 @@ -246,7 +246,7 @@ the user specifies a `SurvivalModel` object to `JointModel()`. **Number of unique subjects** ```stan -int Nind; +int n_subjects; ``` **Number of unique studies** @@ -261,7 +261,7 @@ int n_arms; **Study index for each patient** ```stan -array[Nind] int pt_study_index; +array[n_subjects] int pt_study_index; ``` - Note that this is sorted based upon the patients factor level within the R `data.frame`. For example lets say patient `"A"` has a factor level of 15 and that their corresponding @@ -269,7 +269,7 @@ study value has a factor level of 2 then `pt_study_index[15]` will be 2. **Treatment arm index for each patient** ```stan - array[Nind] int pt_arm_index; + array[n_subjects] int pt_arm_index; ``` - A mirror of `pt_study_index` but for the treatment arm. @@ -284,7 +284,7 @@ int Nind_dead; **Event/Censor Times** ```stan -vector[Nind] Times; +vector[n_subjects] Times; ``` - Ordered by patient factor level From ba073a4cb08d2ce1908dfb192d030534fb60d3cc Mon Sep 17 00:00:00 2001 From: gowerc Date: Tue, 2 Apr 2024 17:43:06 +0100 Subject: [PATCH 4/5] mass update of variable names --- R/DataJoint.R | 10 ++-- R/DataLongitudinal.R | 18 +++---- R/DataSubject.R | 6 +-- R/DataSurvival.R | 6 +-- R/LongitudinalGSF.R | 14 +++--- R/LongitudinalQuantities.R | 2 +- R/LongitudinalRandomSlope.R | 2 +- R/LongitudinalSteinFojo.R | 12 ++--- R/SurvivalQuantities.R | 2 +- inst/stan/base/base.stan | 4 +- inst/stan/base/longitudinal.stan | 28 +++++------ inst/stan/base/survival.stan | 24 ++++----- inst/stan/lm-gsf/link.stan | 2 +- inst/stan/lm-gsf/model.stan | 46 ++++++++--------- inst/stan/lm-random-slope/link.stan | 2 +- inst/stan/lm-random-slope/model.stan | 28 +++++------ inst/stan/lm-stein-fojo/link.stan | 2 +- inst/stan/lm-stein-fojo/model.stan | 42 ++++++++-------- tests/testthat/test-DataJoint.R | 22 ++++----- tests/testthat/test-DataLongitudinal.R | 22 ++++----- tests/testthat/test-DataSubject.R | 10 ++-- tests/testthat/test-DataSurvival.R | 16 +++--- tests/testthat/test-Link.R | 2 +- vignettes/extending-jmpost.Rmd | 68 +++++++++++++------------- 24 files changed, 195 insertions(+), 195 deletions(-) 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/stan/base/base.stan b/inst/stan/base/base.stan index 4b711700..54eacf11 100755 --- a/inst/stan/base/base.stan +++ b/inst/stan/base/base.stan @@ -14,8 +14,8 @@ data{ int n_subjects; // Number of individuals. int n_studies; // Number of studies. int n_arms; // Number of treatment arms. - array[n_subjects] int pt_study_index; // Index of study per pt (PT index sorted) - array[n_subjects] 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 }} diff --git a/inst/stan/base/longitudinal.stan b/inst/stan/base/longitudinal.stan index 30445738..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 n_subjects 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 n_subjects 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 n_subjects 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 }} @@ -67,7 +67,7 @@ transformed parameters { // log_lik += csr_matrix_times_vector( n_subjects, - Nta_total, + 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 d5e503a3..737c79f1 100755 --- a/inst/stan/base/survival.stan +++ b/inst/stan/base/survival.stan @@ -152,9 +152,9 @@ 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[n_subjects] 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[n_subjects, p_os_cov_design] os_cov_design; @@ -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 }} } @@ -212,9 +212,9 @@ transformed parameters { // Log of survival function at the observed time points. vector[n_subjects] log_surv_fit_at_obs_times; - log_surv_fit_at_obs_times = rep_vector(0.0, Nind); + 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 9a147168..37d0626d 100755 --- a/inst/stan/lm-gsf/model.stan +++ b/inst/stan/lm-gsf/model.stan @@ -46,37 +46,37 @@ transformed parameters{ {% if not centred -%} vector[n_subjects] lm_gsf_psi_bsld = exp( - lm_gsf_mu_bsld[pt_study_index] + (lm_gsf_eta_tilde_bsld * lm_gsf_omega_bsld) + lm_gsf_mu_bsld[subject_study_index] + (lm_gsf_eta_tilde_bsld * lm_gsf_omega_bsld) ); vector[n_subjects] lm_gsf_psi_ks = exp( - lm_gsf_mu_ks[pt_arm_index] + (lm_gsf_eta_tilde_ks * lm_gsf_omega_ks) + lm_gsf_mu_ks[subject_arm_index] + (lm_gsf_eta_tilde_ks * lm_gsf_omega_ks) ); vector[n_subjects] lm_gsf_psi_kg = exp( - lm_gsf_mu_kg[pt_arm_index] + (lm_gsf_eta_tilde_kg * lm_gsf_omega_kg) + 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 f768de3b..0b3a93a6 100755 --- a/inst/stan/lm-random-slope/model.stan +++ b/inst/stan/lm-random-slope/model.stan @@ -30,21 +30,21 @@ transformed parameters { // // Source - lm-random-slope/model.stan // - vector[n_subjects] 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 1634ffff..ff0a7ce9 100755 --- a/inst/stan/lm-stein-fojo/model.stan +++ b/inst/stan/lm-stein-fojo/model.stan @@ -41,35 +41,35 @@ transformed parameters{ {% if not centred -%} vector[n_subjects] lm_sf_psi_bsld = exp( - lm_sf_mu_bsld[pt_study_index] + (lm_sf_eta_tilde_bsld * lm_sf_omega_bsld) + lm_sf_mu_bsld[subject_study_index] + (lm_sf_eta_tilde_bsld * lm_sf_omega_bsld) ); vector[n_subjects] lm_sf_psi_ks = exp( - lm_sf_mu_ks[pt_arm_index] + (lm_sf_eta_tilde_ks * lm_sf_omega_ks) + lm_sf_mu_ks[subject_arm_index] + (lm_sf_eta_tilde_ks * lm_sf_omega_ks) ); vector[n_subjects] lm_sf_psi_kg = exp( - lm_sf_mu_kg[pt_arm_index] + (lm_sf_eta_tilde_kg * lm_sf_omega_kg) + 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 c0deec1c..b3bc4194 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; @@ -259,19 +259,19 @@ int n_studies; int n_arms; ``` -**Study index for each patient** +**Study index for each subject** ```stan -array[n_subjects] int pt_study_index; +array[n_subjects] int subject_study_index; ``` -- Note that this is sorted based upon the patients factor level within the R `data.frame`. -For example lets say patient `"A"` has a factor level of 15 and that their corresponding -study value has a factor level of 2 then `pt_study_index[15]` will be 2. +- 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 patient** +**Treatment arm index for each subject** ```stan - array[n_subjects] int pt_arm_index; + array[n_subjects] int subject_arm_index; ``` -- A mirror of `pt_study_index` but for the treatment arm. +- A mirror of `subject_study_index` but for the treatment arm. ### Survival Data Objects @@ -279,20 +279,20 @@ study value has a factor level of 2 then `pt_study_index[15]` will be 2. **Number of events** ```stan -int Nind_dead; +int n_subject_event; ``` **Event/Censor Times** ```stan -vector[n_subjects] Times; +vector[n_subjects] event_times; ``` - Ordered by patient factor level **Event Index** ```stan -array[Nind_dead] int dead_ind_index; +array[n_subject_event] int subject_event_index; ``` -- Is the index into `Times` to identify which times are an event. The rest are censored. +- Is the index into `event_times` to identify which times are an event. The rest are censored. **Number of covariates for the survival model** @@ -303,24 +303,24 @@ int p_os_cov_design; **Covariate design matrix for the survival model** ```stan -matrix[Nind, p_os_cov_design] os_cov_design; +matrix[n_subjects, p_os_cov_design] os_cov_design; ``` - Note that this does not include the intercept term. **Time >0 flag** ```stan -array[rows(Times)] int time_positive; +array[rows(event_times)] int time_positive_flag; ``` **Number of times >0** ```stan -int n_positive; +int n_times_positive; ``` **Positive time index** ```stan -array[n_positive] int time_positive_index +array[n_times_positive] int time_positive_index ``` @@ -337,55 +337,55 @@ vector[n_nodes] weights; **Total number of tumour assessments** ```stan -int Nta_total; +int n_tumour_all; ``` -**Number of tumour assessments above LLoQ** +**Number of tumour assessments above LLoQ (Lower Limit of Quantification)** ```stan -int Nta_obs_y; +int n_tumour_obs; ``` -**Number of tumour assessments below LLoQ** +**Number of tumour assessments below LLoQ (Lower Limit of Quantification)** ```stan -int Nta_cens_y; +int n_tumour_cens; ``` **Tumour assessments values** ```stan -vector[Nta_total] Yobs; +vector[n_tumour_all] tumour_value; ``` **Tumour assessments time points** ```stan -vector[Nta_total] Tobs; +vector[n_tumour_all] tumour_time; ``` **LLoQ threshold** ```stan -real Ythreshold; +real tumour_value_lloq; ``` **Individual tumour assessment index** ```stan -array[Nta_total] int ind_index; +array[n_tumour_all] int subject_tumour_index; ``` - That is if tumour assessment 1 belongs to the patient with factor level 3 then -`ind_index[1]` will be 3. +`subject_tumour_index[1]` will be 3. -**Tumour assessment index for observations above LLoQ** +**Tumour assessment index for observations above LLoQ (Lower Limit of Quantification)** ```stan -array[Nta_obs_y] int obs_y_index; +array[n_tumour_obs] int subject_tumour_index_obs; ``` - For example if only tumour assessments 3 and 5 were above the LLoQ then -`obs_y_index` will be `[3, 5]`. +`subject_tumour_index_obs` will be `[3, 5]`. -**Tumour assessment index for observations below LLoQ** +**Tumour assessment index for observations below LLoQ (Lower Limit of Quantification)** ```stan -array[Nta_cens_y] int cens_y_index; +array[n_tumour_cens] int subject_tumour_index_cens; ``` - For example if only tumour assessments 1 and 2 were below the LLoQ then -`obs_y_index` will be `[1, 2]`. +`subject_tumour_index_obs` will be `[1, 2]`. **Sparse matrix components for patient indexes of the tumour assessments** @@ -395,7 +395,7 @@ 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 `Yobs` belongs to the patient with factor level 3. This matrix is primarily used to calculate the sum +- 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. From 398b39393f0f6072a841468edcf7a58337add665 Mon Sep 17 00:00:00 2001 From: gowerc Date: Wed, 3 Apr 2024 14:41:17 +0100 Subject: [PATCH 5/5] Updated design matrix description --- vignettes/extending-jmpost.Rmd | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/vignettes/extending-jmpost.Rmd b/vignettes/extending-jmpost.Rmd index b3bc4194..dfaa6aaa 100644 --- a/vignettes/extending-jmpost.Rmd +++ b/vignettes/extending-jmpost.Rmd @@ -299,13 +299,13 @@ array[n_subject_event] int subject_event_index; ```stan int p_os_cov_design; ``` -- Note that this does not include the intercept term. +- 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 the intercept term. +- Note that this does not include an intercept term, which would conflict with the baseline distribution parameters. **Time >0 flag**