// generated with brms 2.13.9 data { int N; // number of observations vector[N] Y; // response variable int K_logK1; // number of population-level effects matrix[N, K_logK1] X_logK1; // population-level design matrix int K_logBPnd; // number of population-level effects matrix[N, K_logBPnd] X_logBPnd; // population-level design matrix int K_logVnd; // number of population-level effects matrix[N, K_logVnd] X_logVnd; // population-level design matrix int K_logk4; // number of population-level effects matrix[N, K_logk4] X_logk4; // population-level design matrix int K_logvB; // number of population-level effects matrix[N, K_logvB] X_logvB; // population-level design matrix // covariate vectors for non-linear functions vector[N] C_1; vector[N] C_2; vector[N] C_3; vector[N] C_4; vector[N] C_5; vector[N] C_6; vector[N] C_7; vector[N] C_8; vector[N] C_9; vector[N] C_10; vector[N] C_11; int K_sigma; // number of population-level effects matrix[N, K_sigma] X_sigma; // population-level design matrix // data for splines int Ks_sigma; // number of linear effects matrix[N, Ks_sigma] Xs_sigma; // design matrix for the linear effects // data for spline s(t_tac) int nb_sigma_1; // number of bases int knots_sigma_1[nb_sigma_1]; // number of knots // basis function matrices matrix[N, knots_sigma_1[1]] Zs_sigma_1_1; // data for group-level effects of ID 1 int N_1; // number of grouping levels int M_1; // number of coefficients per level int J_1[N]; // grouping indicator per observation // group-level predictor values vector[N] Z_1_sigma_1; // data for group-level effects of ID 2 int N_2; // number of grouping levels int M_2; // number of coefficients per level int J_2[N]; // grouping indicator per observation // group-level predictor values vector[N] Z_2_logK1_1; vector[N] Z_2_logBPnd_2; vector[N] Z_2_logVnd_3; vector[N] Z_2_logk4_4; int NC_2; // number of group-level correlations // data for group-level effects of ID 3 int N_3; // number of grouping levels int M_3; // number of coefficients per level int J_3[N]; // grouping indicator per observation // group-level predictor values vector[N] Z_3_logK1_1; vector[N] Z_3_logBPnd_2; vector[N] Z_3_logVnd_3; vector[N] Z_3_logk4_4; int NC_3; // number of group-level correlations // data for group-level effects of ID 4 int N_4; // number of grouping levels int M_4; // number of coefficients per level int J_4[N]; // grouping indicator per observation // group-level predictor values vector[N] Z_4_logVnd_1; vector[N] Z_4_logk4_2; vector[N] Z_4_logvB_3; int NC_4; // number of group-level correlations // data for group-level effects of ID 5 int N_5; // number of grouping levels int M_5; // number of coefficients per level int J_5[N]; // grouping indicator per observation // group-level predictor values vector[N] Z_5_logvB_1; int prior_only; // should the likelihood be ignored? } transformed data { int Kc_sigma = K_sigma - 1; matrix[N, Kc_sigma] Xc_sigma; // centered version of X_sigma without an intercept vector[Kc_sigma] means_X_sigma; // column means of X_sigma before centering for (i in 2:K_sigma) { means_X_sigma[i - 1] = mean(X_sigma[, i]); Xc_sigma[, i - 1] = X_sigma[, i] - means_X_sigma[i - 1]; } } parameters { vector[K_logK1] b_logK1; // population-level effects vector[K_logBPnd] b_logBPnd; // population-level effects vector[K_logVnd] b_logVnd; // population-level effects vector[K_logk4] b_logk4; // population-level effects vector[K_logvB] b_logvB; // population-level effects vector[Kc_sigma] b_sigma; // population-level effects real Intercept_sigma; // temporary intercept for centered predictors vector[Ks_sigma] bs_sigma; // spline coefficients // parameters for spline s(t_tac) // standarized spline coefficients vector[knots_sigma_1[1]] zs_sigma_1_1; real sds_sigma_1_1; // standard deviations of spline coefficients vector[M_1] sd_1; // group-level standard deviations vector[N_1] z_1[M_1]; // standardized group-level effects vector[M_2] sd_2; // group-level standard deviations matrix[M_2, N_2] z_2; // standardized group-level effects cholesky_factor_corr[M_2] L_2; // cholesky factor of correlation matrix vector[M_3] sd_3; // group-level standard deviations matrix[M_3, N_3] z_3; // standardized group-level effects cholesky_factor_corr[M_3] L_3; // cholesky factor of correlation matrix vector[M_4] sd_4; // group-level standard deviations matrix[M_4, N_4] z_4; // standardized group-level effects cholesky_factor_corr[M_4] L_4; // cholesky factor of correlation matrix vector[M_5] sd_5; // group-level standard deviations vector[N_5] z_5[M_5]; // standardized group-level effects } transformed parameters { // actual spline coefficients vector[knots_sigma_1[1]] s_sigma_1_1; vector[N_1] r_1_sigma_1; // actual group-level effects matrix[N_2, M_2] r_2; // actual group-level effects // using vectors speeds up indexing in loops vector[N_2] r_2_logK1_1; vector[N_2] r_2_logBPnd_2; vector[N_2] r_2_logVnd_3; vector[N_2] r_2_logk4_4; matrix[N_3, M_3] r_3; // actual group-level effects // using vectors speeds up indexing in loops vector[N_3] r_3_logK1_1; vector[N_3] r_3_logBPnd_2; vector[N_3] r_3_logVnd_3; vector[N_3] r_3_logk4_4; matrix[N_4, M_4] r_4; // actual group-level effects // using vectors speeds up indexing in loops vector[N_4] r_4_logVnd_1; vector[N_4] r_4_logk4_2; vector[N_4] r_4_logvB_3; vector[N_5] r_5_logvB_1; // actual group-level effects // compute actual spline coefficients s_sigma_1_1 = sds_sigma_1_1 * zs_sigma_1_1; r_1_sigma_1 = (sd_1[1] * (z_1[1])); // compute actual group-level effects r_2 = (diag_pre_multiply(sd_2, L_2) * z_2)'; r_2_logK1_1 = r_2[, 1]; r_2_logBPnd_2 = r_2[, 2]; r_2_logVnd_3 = r_2[, 3]; r_2_logk4_4 = r_2[, 4]; // compute actual group-level effects r_3 = (diag_pre_multiply(sd_3, L_3) * z_3)'; r_3_logK1_1 = r_3[, 1]; r_3_logBPnd_2 = r_3[, 2]; r_3_logVnd_3 = r_3[, 3]; r_3_logk4_4 = r_3[, 4]; // compute actual group-level effects r_4 = (diag_pre_multiply(sd_4, L_4) * z_4)'; r_4_logVnd_1 = r_4[, 1]; r_4_logk4_2 = r_4[, 2]; r_4_logvB_3 = r_4[, 3]; r_5_logvB_1 = (sd_5[1] * (z_5[1])); } model { // initialize linear predictor term vector[N] nlp_logK1 = X_logK1 * b_logK1; // initialize linear predictor term vector[N] nlp_logBPnd = X_logBPnd * b_logBPnd; // initialize linear predictor term vector[N] nlp_logVnd = X_logVnd * b_logVnd; // initialize linear predictor term vector[N] nlp_logk4 = X_logk4 * b_logk4; // initialize linear predictor term vector[N] nlp_logvB = X_logvB * b_logvB; // initialize non-linear predictor term vector[N] mu; // initialize linear predictor term vector[N] sigma = Intercept_sigma + Xc_sigma * b_sigma + Xs_sigma * bs_sigma + Zs_sigma_1_1 * s_sigma_1_1; for (n in 1:N) { // add more terms to the linear predictor nlp_logK1[n] += r_2_logK1_1[J_2[n]] * Z_2_logK1_1[n] + r_3_logK1_1[J_3[n]] * Z_3_logK1_1[n]; } for (n in 1:N) { // add more terms to the linear predictor nlp_logBPnd[n] += r_2_logBPnd_2[J_2[n]] * Z_2_logBPnd_2[n] + r_3_logBPnd_2[J_3[n]] * Z_3_logBPnd_2[n]; } for (n in 1:N) { // add more terms to the linear predictor nlp_logVnd[n] += r_2_logVnd_3[J_2[n]] * Z_2_logVnd_3[n] + r_3_logVnd_3[J_3[n]] * Z_3_logVnd_3[n] + r_4_logVnd_1[J_4[n]] * Z_4_logVnd_1[n]; } for (n in 1:N) { // add more terms to the linear predictor nlp_logk4[n] += r_2_logk4_4[J_2[n]] * Z_2_logk4_4[n] + r_3_logk4_4[J_3[n]] * Z_3_logk4_4[n] + r_4_logk4_2[J_4[n]] * Z_4_logk4_2[n]; } for (n in 1:N) { // add more terms to the linear predictor nlp_logvB[n] += r_4_logvB_3[J_4[n]] * Z_4_logvB_3[n] + r_5_logvB_1[J_5[n]] * Z_5_logvB_1[n]; } for (n in 1:N) { // add more terms to the linear predictor sigma[n] += r_1_sigma_1[J_1[n]] * Z_1_sigma_1[n]; } for (n in 1:N) { // apply the inverse link function sigma[n] = exp(sigma[n]); } for (n in 1:N) { // compute non-linear predictor values mu[n] = twotcm_log_stan(nlp_logK1[n] , nlp_logVnd[n] , nlp_logBPnd[n] , nlp_logk4[n] , nlp_logvB[n] , C_1[n] , C_2[n] , C_3[n] , C_4[n] , C_5[n] , C_6[n] , C_7[n] , C_8[n] , C_9[n] , C_10[n] , C_11[n]); } // priors including all constants target += normal_lpdf(b_logK1[1] | -2.5, 0.5); target += normal_lpdf(b_logK1[2] | 0, 0.3); target += normal_lpdf(b_logK1[3] | 0, 0.3); target += normal_lpdf(b_logK1[4] | 0, 0.3); target += normal_lpdf(b_logK1[5] | 0, 0.3); target += normal_lpdf(b_logK1[6] | 0, 0.3); target += normal_lpdf(b_logK1[7] | 0, 0.3); target += normal_lpdf(b_logK1[8] | 0, 0.3); target += normal_lpdf(b_logK1[9] | 0, 0.3); target += normal_lpdf(b_logK1[10] | 0, 0.05); target += normal_lpdf(b_logK1[11] | 0, 0.05); target += normal_lpdf(b_logBPnd[1] | 2, 0.5); target += student_t_lpdf(b_logBPnd[2] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[3] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[4] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[5] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[6] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[7] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[8] | 3, 0, 0.3); target += student_t_lpdf(b_logBPnd[9] | 3, 0, 0.3); target += normal_lpdf(b_logBPnd[10] | 0, 0.1 ); target += normal_lpdf(b_logBPnd[11] | 0, 0.1 ); target += normal_lpdf(b_logBPnd[12] | 0, 0.1); target += normal_lpdf(b_logBPnd[13] | 0, 0.1); target += normal_lpdf(b_logBPnd[14] | 0, 0.05); target += normal_lpdf(b_logBPnd[15] | 0, 0.05); target += normal_lpdf(b_logBPnd[16] | 0, 0.05); target += normal_lpdf(b_logBPnd[17] | 0, 0.05); target += normal_lpdf(b_logBPnd[18] | 0, 0.05); target += normal_lpdf(b_logBPnd[19] | 0, 0.05); target += normal_lpdf(b_logBPnd[20] | 0, 0.05); target += normal_lpdf(b_logBPnd[21] | 0, 0.05); target += normal_lpdf(b_logBPnd[22] | 0, 0.05); target += normal_lpdf(b_logBPnd[23] | 0, 0.05); target += normal_lpdf(b_logBPnd[24] | 0, 0.05); target += normal_lpdf(b_logBPnd[25] | 0, 0.05); target += normal_lpdf(b_logBPnd[26] | 0, 0.05); target += normal_lpdf(b_logBPnd[27] | 0, 0.05); target += normal_lpdf(b_logBPnd[28] | 0, 0.05); target += normal_lpdf(b_logBPnd[29] | 0, 0.05); target += normal_lpdf(b_logVnd | -1, 0.5); target += normal_lpdf(b_logk4 | -4, 0.5); target += normal_lpdf(b_logvB[1] | -4, 1); target += normal_lpdf(b_logvB[2] | 0, 0.05); target += normal_lpdf(b_logvB[3] | 0, 0.05); target += normal_lpdf(b_sigma[1] | 0, 0.2); target += normal_lpdf(b_sigma[2] | 0, 0.2); target += normal_lpdf(b_sigma[3] | 0, 0.2); target += normal_lpdf(b_sigma[4] | 0, 0.2); target += normal_lpdf(b_sigma[5] | 0, 0.2); target += normal_lpdf(b_sigma[6] | 0, 0.2); target += normal_lpdf(b_sigma[7] | 0, 0.2); target += normal_lpdf(b_sigma[8] | 0, 0.2); target += normal_lpdf(Intercept_sigma | -5, 1); target += student_t_lpdf(sds_sigma_1_1 | 3, 0, 2.5) - 1 * student_t_lccdf(0 | 3, 0, 2.5); target += std_normal_lpdf(zs_sigma_1_1); target += normal_lpdf(sd_1 | 0, 0.3) - 1 * normal_lccdf(0 | 0, 0.3); target += std_normal_lpdf(z_1[1]); target += student_t_lpdf(sd_2 | 3, 0, 0.2) - 4 * student_t_lccdf(0 | 3, 0, 0.2); target += std_normal_lpdf(to_vector(z_2)); target += lkj_corr_cholesky_lpdf(L_2 | 1); target += normal_lpdf(sd_3[1] | 0, 0.01) - 1 * normal_lccdf(0 | 0, 0.01); target += normal_lpdf(sd_3[2] | 0, 0.02) - 1 * normal_lccdf(0 | 0, 0.02); target += normal_lpdf(sd_3[3] | 0, 0.01) - 1 * normal_lccdf(0 | 0, 0.01); target += normal_lpdf(sd_3[4] | 0, 0.01) - 1 * normal_lccdf(0 | 0, 0.01); target += std_normal_lpdf(to_vector(z_3)); target += lkj_corr_cholesky_lpdf(L_3 | 1); target += normal_lpdf(sd_4 | 0, 0.1) - 3 * normal_lccdf(0 | 0, 0.1); target += std_normal_lpdf(to_vector(z_4)); target += lkj_corr_cholesky_lpdf(L_4 | 1); target += student_t_lpdf(sd_5 | 3, 0, 0.2) - 1 * student_t_lccdf(0 | 3, 0, 0.2); target += std_normal_lpdf(z_5[1]); // likelihood including all constants if (!prior_only) { target += normal_lpdf(Y | mu, sigma); } } generated quantities { // actual population-level intercept real b_sigma_Intercept = Intercept_sigma - dot_product(means_X_sigma, b_sigma); // compute group-level correlations corr_matrix[M_2] Cor_2 = multiply_lower_tri_self_transpose(L_2); vector[NC_2] cor_2; // compute group-level correlations corr_matrix[M_3] Cor_3 = multiply_lower_tri_self_transpose(L_3); vector[NC_3] cor_3; // compute group-level correlations corr_matrix[M_4] Cor_4 = multiply_lower_tri_self_transpose(L_4); vector[NC_4] cor_4; // extract upper diagonal of correlation matrix for (k in 1:M_2) { for (j in 1:(k - 1)) { cor_2[choose(k - 1, 2) + j] = Cor_2[j, k]; } } // extract upper diagonal of correlation matrix for (k in 1:M_3) { for (j in 1:(k - 1)) { cor_3[choose(k - 1, 2) + j] = Cor_3[j, k]; } } // extract upper diagonal of correlation matrix for (k in 1:M_4) { for (j in 1:(k - 1)) { cor_4[choose(k - 1, 2) + j] = Cor_4[j, k]; } } }