From 1d8ad2c4aa414b09f3e9c05be0ee8bb09139957b Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Mon, 6 Jan 2025 15:05:29 +1100 Subject: [PATCH 1/7] Removed redundant global error parameter from affine model. --- inst/stan/affine_single_ind.stan | 3 --- 1 file changed, 3 deletions(-) diff --git a/inst/stan/affine_single_ind.stan b/inst/stan/affine_single_ind.stan index 04425fc..e16a7a0 100644 --- a/inst/stan/affine_single_ind.stan +++ b/inst/stan/affine_single_ind.stan @@ -76,9 +76,6 @@ parameters { real ind_y_0; real ind_const; real ind_beta_1; - - //Global level - real global_error_sigma; } // The model to be estimated. From 6f7d5c7125fbe92b286d54e9ae9419668ffcfd01 Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Mon, 6 Jan 2025 23:22:46 +1100 Subject: [PATCH 2/7] Simplified affine DE. --- R/hmde_models.R | 1 - inst/stan/affine_single_ind.stan | 46 ++++++++++-------------- tests/testthat/test-hmde_models_affine.R | 2 +- 3 files changed, 19 insertions(+), 30 deletions(-) diff --git a/R/hmde_models.R b/R/hmde_models.R index 80f5c1c..43ed631 100644 --- a/R/hmde_models.R +++ b/R/hmde_models.R @@ -116,7 +116,6 @@ hmde_affine_single_ind <- function(){ y_obs = NULL, obs_index = NULL, time = NULL, - y_bar = NULL, int_method = NULL, prior_means = c(1,1), prior_sds = c(2,2), diff --git a/inst/stan/affine_single_ind.stan b/inst/stan/affine_single_ind.stan index e16a7a0..c06822f 100644 --- a/inst/stan/affine_single_ind.stan +++ b/inst/stan/affine_single_ind.stan @@ -1,9 +1,9 @@ //Growth function functions{ //Growth function for use with Runge-Kutta method - //pars = (beta_0, beta_1, y_bar) + //pars = (beta_0, beta_1) real DE_rk4(real y, array[] real pars){ //change number of pars - return pars[1] - (pars[2] * (y-pars[3])); //growth function + return pars[1] - pars[2] * y; //growth function } real rk4_step(real y, array[] real pars, real interval){ @@ -46,14 +46,13 @@ functions{ return y_hat; } - vector DE_RK45(real t, vector y, real ind_const, real beta_1, real y_bar){ - vector[size(y)] dydt = ind_const - (beta_1 * (y-y_bar)); + vector DE_RK45(real t, vector y, real beta_0, real beta_1){ + vector[size(y)] dydt = beta_0 - beta_1 * y; return dydt; } - real analytic_solution(real t, real y_0, real ind_const, real beta_1, real y_bar){ - real y_0_translate = y_0 - y_bar; - return ind_const/beta_1 + (y_0_translate - (ind_const/beta_1)) * exp(-beta_1 * t) + y_bar; + real analytic_solution(real t, real y_0, real beta_0, real beta_1){ + return beta_0/beta_1 + (y_0 - (beta_0/beta_1)) * exp(-beta_1 * t); } } @@ -64,7 +63,6 @@ data { real y_obs[n_obs]; int obs_index[n_obs]; real time[n_obs]; - real y_bar; int int_method; real prior_means[2]; #vector of means for beta parameter priors real prior_sds[2]; #Vector of SDs for beta parameter priors @@ -74,19 +72,18 @@ data { parameters { //Individual level real ind_y_0; - real ind_const; + real ind_beta_0; real ind_beta_1; } // The model to be estimated. model { real y_hat[n_obs]; - array[3] real pars; + array[2] real pars; vector[1] y_temp; - pars[1] = ind_const; + pars[1] = ind_beta_0; pars[2] = ind_beta_1; - pars[3] = y_bar; for(i in 1:n_obs){ @@ -104,15 +101,14 @@ model { y_temp[1] = y_hat[i]; y_hat[i+1] = ode_rk45(DE_RK45, y_temp, time[i], {time[i+1]}, - ind_const, ind_beta_1, y_bar)[1][1]; + ind_beta_0, ind_beta_1)[1][1]; } if(int_method == 3){ //Analytic solution y_hat[i+1] = analytic_solution((time[i+1] - time[1]), ind_y_0, - ind_const, - ind_beta_1, - y_bar); + ind_beta_0, + ind_beta_1); } } } @@ -122,22 +118,17 @@ model { //Priors //Individual level - ind_const ~lognormal(log(prior_means[1]), prior_sds[1]); + ind_beta_0 ~lognormal(log(prior_means[1]), prior_sds[1]); ind_beta_1 ~lognormal(log(prior_means[2]), prior_sds[2]); } generated quantities{ real y_hat[n_obs]; - array[3] real pars; - real ind_beta_0; + array[2] real pars; vector[1] y_temp; - int vers = 1; - ind_beta_0 = ind_const + ind_beta_1*y_bar; - - pars[1] = ind_const; + pars[1] = ind_beta_0; pars[2] = ind_beta_1; - pars[3] = y_bar; for(i in 1:n_obs){ @@ -155,15 +146,14 @@ generated quantities{ y_temp[1] = y_hat[i]; y_hat[i+1] = ode_rk45(DE_RK45, y_temp, time[i], {time[i+1]}, - ind_const, ind_beta_1, y_bar)[1][1]; + ind_beta_0, ind_beta_1)[1][1]; } if(int_method == 3){ //Analytic solution y_hat[i+1] = analytic_solution((time[i+1] - time[1]), ind_y_0, - ind_const, - ind_beta_1, - y_bar); + ind_beta_0, + ind_beta_1); } } } diff --git a/tests/testthat/test-hmde_models_affine.R b/tests/testthat/test-hmde_models_affine.R index b9014f6..72683b4 100644 --- a/tests/testthat/test-hmde_models_affine.R +++ b/tests/testthat/test-hmde_models_affine.R @@ -3,7 +3,7 @@ test_that("Model structures: affine", { # Single individual single_model <- hmde_model("affine_single_ind") expect_named(single_model, c("step_size", "n_obs", "y_obs", "obs_index", - "time", "y_bar", "int_method", "prior_means", + "time", "int_method", "prior_means", "prior_sds", "model")) expect_type(single_model, "list") expect_visible(single_model) From 0112570fb260b624d5378e2027483681bf8f9423 Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Tue, 7 Jan 2025 14:07:11 +1100 Subject: [PATCH 3/7] Revert "Removed redundant global error parameter from affine model." This reverts commit 1d8ad2c4aa414b09f3e9c05be0ee8bb09139957b. Reverting changes that caused estimation problems. --- inst/stan/affine_single_ind.stan | 3 +++ 1 file changed, 3 insertions(+) diff --git a/inst/stan/affine_single_ind.stan b/inst/stan/affine_single_ind.stan index e16a7a0..04425fc 100644 --- a/inst/stan/affine_single_ind.stan +++ b/inst/stan/affine_single_ind.stan @@ -76,6 +76,9 @@ parameters { real ind_y_0; real ind_const; real ind_beta_1; + + //Global level + real global_error_sigma; } // The model to be estimated. From ce7f5bb690d917c88f395bc2b2f8c37127bfd8c2 Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Tue, 7 Jan 2025 14:08:57 +1100 Subject: [PATCH 4/7] Changed parameter name in comment to reflect parameterisation. --- inst/stan/affine_single_ind.stan | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/inst/stan/affine_single_ind.stan b/inst/stan/affine_single_ind.stan index 04425fc..c79a94d 100644 --- a/inst/stan/affine_single_ind.stan +++ b/inst/stan/affine_single_ind.stan @@ -1,7 +1,7 @@ //Growth function functions{ //Growth function for use with Runge-Kutta method - //pars = (beta_0, beta_1, y_bar) + //pars = (ind_const, beta_1, y_bar) real DE_rk4(real y, array[] real pars){ //change number of pars return pars[1] - (pars[2] * (y-pars[3])); //growth function } From 2e7964b72ab16ca353cd79d10de6a9e25a53e009 Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Tue, 7 Jan 2025 14:12:13 +1100 Subject: [PATCH 5/7] Put y_bar back in. --- R/hmde_models.R | 1 + 1 file changed, 1 insertion(+) diff --git a/R/hmde_models.R b/R/hmde_models.R index 43ed631..dac41e9 100644 --- a/R/hmde_models.R +++ b/R/hmde_models.R @@ -117,6 +117,7 @@ hmde_affine_single_ind <- function(){ obs_index = NULL, time = NULL, int_method = NULL, + y_bar = NULL, prior_means = c(1,1), prior_sds = c(2,2), model = "affine_single_ind") From b1b23ff860773782aca5565f00439b3d76f5552f Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Tue, 7 Jan 2025 14:29:48 +1100 Subject: [PATCH 6/7] Returned y_bar to test. --- tests/testthat/test-hmde_models_affine.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/testthat/test-hmde_models_affine.R b/tests/testthat/test-hmde_models_affine.R index 72683b4..b6000b5 100644 --- a/tests/testthat/test-hmde_models_affine.R +++ b/tests/testthat/test-hmde_models_affine.R @@ -3,7 +3,7 @@ test_that("Model structures: affine", { # Single individual single_model <- hmde_model("affine_single_ind") expect_named(single_model, c("step_size", "n_obs", "y_obs", "obs_index", - "time", "int_method", "prior_means", + "time", "int_method", "y_bar", "prior_means", "prior_sds", "model")) expect_type(single_model, "list") expect_visible(single_model) From 3435682c09f1378b3d1b7a23563ec4ae2240cbc4 Mon Sep 17 00:00:00 2001 From: Tess-LaCoil Date: Tue, 7 Jan 2025 15:30:28 +1100 Subject: [PATCH 7/7] Removed constraint on ind_const to be non-negative in order to permit initialisation testing. --- inst/stan/affine_single_ind.stan | 9 +++------ 1 file changed, 3 insertions(+), 6 deletions(-) diff --git a/inst/stan/affine_single_ind.stan b/inst/stan/affine_single_ind.stan index c79a94d..21166ac 100644 --- a/inst/stan/affine_single_ind.stan +++ b/inst/stan/affine_single_ind.stan @@ -1,7 +1,7 @@ //Growth function functions{ //Growth function for use with Runge-Kutta method - //pars = (ind_const, beta_1, y_bar) + //pars = (const, beta_1, y_bar) real DE_rk4(real y, array[] real pars){ //change number of pars return pars[1] - (pars[2] * (y-pars[3])); //growth function } @@ -74,11 +74,8 @@ data { parameters { //Individual level real ind_y_0; - real ind_const; + real ind_const; real ind_beta_1; - - //Global level - real global_error_sigma; } // The model to be estimated. @@ -125,7 +122,7 @@ model { //Priors //Individual level - ind_const ~lognormal(log(prior_means[1]), prior_sds[1]); + ind_const ~normal(prior_means[1], prior_sds[1]); ind_beta_1 ~lognormal(log(prior_means[2]), prior_sds[2]); }