Please also provide the following information in addition to your question:
-
Operating System:
Ubuntu 19.10 -
brms Version:
2.12.0
I am relatively new to Bayesian inference and brms has made the transition very approachable!
I am currently trying to analyze some multinomial data (8 possible categorical in each row) with a hierarchical model.
The data are generated from 16 groups of people (8 groups in “control”, 8 groups in "intervention) across 5 sessions, led by 2 different facilitators.
The model is taking a very long time to run. Most recently, I ran the model with 2 chains, for 2000 iterations with 1000 warmup iterations. The first chain took 19 hours, the second took 50 hours and still had a couple of divergent transitions and low ESS.
I’m wondering if Stan’s difficulty is from how I’ve specified the model, or if it’s just the complexity of the data leading to identifiability issues.
Is reparameterization something I could/should consider? Is this possible using the brms interface, or does it require directly coding in Stan?
My brms model is as follows:
library(tidyverse)
library(brms)
## Model formula
What_Multi_fmla <- bf(
Codes | trials(Count) ~ Condition +
Session_Type +
Facilitator +
Condition:Session_Type +
Facilitator:Session_Type +
(1 | Condition / Group)
)
# Parse data
Data <- read.csv("What_Wider.csv") %>%
mutate(Tran_Type = ifelse(Session == "Instruction", NA, as.character(Tran_Type))) %>%
unite("Session_Type", c(Session, Tran_Type)) %>%
mutate(
Group = factor(Group),
Session_Type = str_remove(Session_Type, "_NA"),
Session_Type = factor(
Session_Type,
levels = c(
"Instruction",
"Initial_Homecare",
"Initial_Debrief",
"Follow-up_Homecare",
"Follow-up_Debrief"
)
)
) %>%
droplevels() %>%
data.frame()
## Bind outcome columns to use as dependent variable in multinomial analysis
What_Wider$Codes <-
with(
What_Wider,
cbind(
NoCode,
Building,
Caregiver,
Patient,
Patientpsych,
SDOH,
Canmeds,
Expertise,
Recommendations
)
)
The generated Stan code is as follows:
functions {
/* multinomial-logit log-PMF
* Args:
* y: array of integer response values
* mu: vector of category logit probabilities
* Returns:
* a scalar to be added to the log posterior
*/
real multinomial_logit_lpmf(int[] y, vector mu) {
return multinomial_lpmf(y | softmax(mu));
}
}
data {
int<lower=1> N; // number of observations
int<lower=2> ncat; // number of categories
int Y[N, ncat]; // response array
int trials[N]; // number of trials
int<lower=1> K_muBuilding; // number of population-level effects
matrix[N, K_muBuilding] X_muBuilding; // population-level design matrix
int<lower=1> K_muCaregiver; // number of population-level effects
matrix[N, K_muCaregiver] X_muCaregiver; // population-level design matrix
int<lower=1> K_muPatient; // number of population-level effects
matrix[N, K_muPatient] X_muPatient; // population-level design matrix
int<lower=1> K_muPatientpsych; // number of population-level effects
matrix[N, K_muPatientpsych] X_muPatientpsych; // population-level design matrix
int<lower=1> K_muSDOH; // number of population-level effects
matrix[N, K_muSDOH] X_muSDOH; // population-level design matrix
int<lower=1> K_muCanmeds; // number of population-level effects
matrix[N, K_muCanmeds] X_muCanmeds; // population-level design matrix
int<lower=1> K_muExpertise; // number of population-level effects
matrix[N, K_muExpertise] X_muExpertise; // population-level design matrix
int<lower=1> K_muRecommendations; // number of population-level effects
matrix[N, K_muRecommendations] X_muRecommendations; // population-level design matrix
// data for group-level effects of ID 1
int<lower=1> N_1; // number of grouping levels
int<lower=1> M_1; // number of coefficients per level
int<lower=1> J_1[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_1_muBuilding_1;
// data for group-level effects of ID 2
int<lower=1> N_2; // number of grouping levels
int<lower=1> M_2; // number of coefficients per level
int<lower=1> J_2[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_2_muBuilding_1;
// data for group-level effects of ID 3
int<lower=1> N_3; // number of grouping levels
int<lower=1> M_3; // number of coefficients per level
int<lower=1> J_3[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_3_muCaregiver_1;
// data for group-level effects of ID 4
int<lower=1> N_4; // number of grouping levels
int<lower=1> M_4; // number of coefficients per level
int<lower=1> J_4[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_4_muCaregiver_1;
// data for group-level effects of ID 5
int<lower=1> N_5; // number of grouping levels
int<lower=1> M_5; // number of coefficients per level
int<lower=1> J_5[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_5_muPatient_1;
// data for group-level effects of ID 6
int<lower=1> N_6; // number of grouping levels
int<lower=1> M_6; // number of coefficients per level
int<lower=1> J_6[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_6_muPatient_1;
// data for group-level effects of ID 7
int<lower=1> N_7; // number of grouping levels
int<lower=1> M_7; // number of coefficients per level
int<lower=1> J_7[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_7_muPatientpsych_1;
// data for group-level effects of ID 8
int<lower=1> N_8; // number of grouping levels
int<lower=1> M_8; // number of coefficients per level
int<lower=1> J_8[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_8_muPatientpsych_1;
// data for group-level effects of ID 9
int<lower=1> N_9; // number of grouping levels
int<lower=1> M_9; // number of coefficients per level
int<lower=1> J_9[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_9_muSDOH_1;
// data for group-level effects of ID 10
int<lower=1> N_10; // number of grouping levels
int<lower=1> M_10; // number of coefficients per level
int<lower=1> J_10[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_10_muSDOH_1;
// data for group-level effects of ID 11
int<lower=1> N_11; // number of grouping levels
int<lower=1> M_11; // number of coefficients per level
int<lower=1> J_11[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_11_muCanmeds_1;
// data for group-level effects of ID 12
int<lower=1> N_12; // number of grouping levels
int<lower=1> M_12; // number of coefficients per level
int<lower=1> J_12[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_12_muCanmeds_1;
// data for group-level effects of ID 13
int<lower=1> N_13; // number of grouping levels
int<lower=1> M_13; // number of coefficients per level
int<lower=1> J_13[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_13_muExpertise_1;
// data for group-level effects of ID 14
int<lower=1> N_14; // number of grouping levels
int<lower=1> M_14; // number of coefficients per level
int<lower=1> J_14[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_14_muExpertise_1;
// data for group-level effects of ID 15
int<lower=1> N_15; // number of grouping levels
int<lower=1> M_15; // number of coefficients per level
int<lower=1> J_15[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_15_muRecommendations_1;
// data for group-level effects of ID 16
int<lower=1> N_16; // number of grouping levels
int<lower=1> M_16; // number of coefficients per level
int<lower=1> J_16[N]; // grouping indicator per observation
// group-level predictor values
vector[N] Z_16_muRecommendations_1;
int prior_only; // should the likelihood be ignored?
}
transformed data {
int Kc_muBuilding = K_muBuilding - 1;
matrix[N, Kc_muBuilding] Xc_muBuilding; // centered version of X_muBuilding without an intercept
vector[Kc_muBuilding] means_X_muBuilding; // column means of X_muBuilding before centering
int Kc_muCaregiver = K_muCaregiver - 1;
matrix[N, Kc_muCaregiver] Xc_muCaregiver; // centered version of X_muCaregiver without an intercept
vector[Kc_muCaregiver] means_X_muCaregiver; // column means of X_muCaregiver before centering
int Kc_muPatient = K_muPatient - 1;
matrix[N, Kc_muPatient] Xc_muPatient; // centered version of X_muPatient without an intercept
vector[Kc_muPatient] means_X_muPatient; // column means of X_muPatient before centering
int Kc_muPatientpsych = K_muPatientpsych - 1;
matrix[N, Kc_muPatientpsych] Xc_muPatientpsych; // centered version of X_muPatientpsych without an intercept
vector[Kc_muPatientpsych] means_X_muPatientpsych; // column means of X_muPatientpsych before centering
int Kc_muSDOH = K_muSDOH - 1;
matrix[N, Kc_muSDOH] Xc_muSDOH; // centered version of X_muSDOH without an intercept
vector[Kc_muSDOH] means_X_muSDOH; // column means of X_muSDOH before centering
int Kc_muCanmeds = K_muCanmeds - 1;
matrix[N, Kc_muCanmeds] Xc_muCanmeds; // centered version of X_muCanmeds without an intercept
vector[Kc_muCanmeds] means_X_muCanmeds; // column means of X_muCanmeds before centering
int Kc_muExpertise = K_muExpertise - 1;
matrix[N, Kc_muExpertise] Xc_muExpertise; // centered version of X_muExpertise without an intercept
vector[Kc_muExpertise] means_X_muExpertise; // column means of X_muExpertise before centering
int Kc_muRecommendations = K_muRecommendations - 1;
matrix[N, Kc_muRecommendations] Xc_muRecommendations; // centered version of X_muRecommendations without an intercept
vector[Kc_muRecommendations] means_X_muRecommendations; // column means of X_muRecommendations before centering
for (i in 2:K_muBuilding) {
means_X_muBuilding[i - 1] = mean(X_muBuilding[, i]);
Xc_muBuilding[, i - 1] = X_muBuilding[, i] - means_X_muBuilding[i - 1];
}
for (i in 2:K_muCaregiver) {
means_X_muCaregiver[i - 1] = mean(X_muCaregiver[, i]);
Xc_muCaregiver[, i - 1] = X_muCaregiver[, i] - means_X_muCaregiver[i - 1];
}
for (i in 2:K_muPatient) {
means_X_muPatient[i - 1] = mean(X_muPatient[, i]);
Xc_muPatient[, i - 1] = X_muPatient[, i] - means_X_muPatient[i - 1];
}
for (i in 2:K_muPatientpsych) {
means_X_muPatientpsych[i - 1] = mean(X_muPatientpsych[, i]);
Xc_muPatientpsych[, i - 1] = X_muPatientpsych[, i] - means_X_muPatientpsych[i - 1];
}
for (i in 2:K_muSDOH) {
means_X_muSDOH[i - 1] = mean(X_muSDOH[, i]);
Xc_muSDOH[, i - 1] = X_muSDOH[, i] - means_X_muSDOH[i - 1];
}
for (i in 2:K_muCanmeds) {
means_X_muCanmeds[i - 1] = mean(X_muCanmeds[, i]);
Xc_muCanmeds[, i - 1] = X_muCanmeds[, i] - means_X_muCanmeds[i - 1];
}
for (i in 2:K_muExpertise) {
means_X_muExpertise[i - 1] = mean(X_muExpertise[, i]);
Xc_muExpertise[, i - 1] = X_muExpertise[, i] - means_X_muExpertise[i - 1];
}
for (i in 2:K_muRecommendations) {
means_X_muRecommendations[i - 1] = mean(X_muRecommendations[, i]);
Xc_muRecommendations[, i - 1] = X_muRecommendations[, i] - means_X_muRecommendations[i - 1];
}
}
parameters {
vector[Kc_muBuilding] b_muBuilding; // population-level effects
real Intercept_muBuilding; // temporary intercept for centered predictors
vector[Kc_muCaregiver] b_muCaregiver; // population-level effects
real Intercept_muCaregiver; // temporary intercept for centered predictors
vector[Kc_muPatient] b_muPatient; // population-level effects
real Intercept_muPatient; // temporary intercept for centered predictors
vector[Kc_muPatientpsych] b_muPatientpsych; // population-level effects
real Intercept_muPatientpsych; // temporary intercept for centered predictors
vector[Kc_muSDOH] b_muSDOH; // population-level effects
real Intercept_muSDOH; // temporary intercept for centered predictors
vector[Kc_muCanmeds] b_muCanmeds; // population-level effects
real Intercept_muCanmeds; // temporary intercept for centered predictors
vector[Kc_muExpertise] b_muExpertise; // population-level effects
real Intercept_muExpertise; // temporary intercept for centered predictors
vector[Kc_muRecommendations] b_muRecommendations; // population-level effects
real Intercept_muRecommendations; // temporary intercept for centered predictors
vector<lower=0>[M_1] sd_1; // group-level standard deviations
vector[N_1] z_1[M_1]; // standardized group-level effects
vector<lower=0>[M_2] sd_2; // group-level standard deviations
vector[N_2] z_2[M_2]; // standardized group-level effects
vector<lower=0>[M_3] sd_3; // group-level standard deviations
vector[N_3] z_3[M_3]; // standardized group-level effects
vector<lower=0>[M_4] sd_4; // group-level standard deviations
vector[N_4] z_4[M_4]; // standardized group-level effects
vector<lower=0>[M_5] sd_5; // group-level standard deviations
vector[N_5] z_5[M_5]; // standardized group-level effects
vector<lower=0>[M_6] sd_6; // group-level standard deviations
vector[N_6] z_6[M_6]; // standardized group-level effects
vector<lower=0>[M_7] sd_7; // group-level standard deviations
vector[N_7] z_7[M_7]; // standardized group-level effects
vector<lower=0>[M_8] sd_8; // group-level standard deviations
vector[N_8] z_8[M_8]; // standardized group-level effects
vector<lower=0>[M_9] sd_9; // group-level standard deviations
vector[N_9] z_9[M_9]; // standardized group-level effects
vector<lower=0>[M_10] sd_10; // group-level standard deviations
vector[N_10] z_10[M_10]; // standardized group-level effects
vector<lower=0>[M_11] sd_11; // group-level standard deviations
vector[N_11] z_11[M_11]; // standardized group-level effects
vector<lower=0>[M_12] sd_12; // group-level standard deviations
vector[N_12] z_12[M_12]; // standardized group-level effects
vector<lower=0>[M_13] sd_13; // group-level standard deviations
vector[N_13] z_13[M_13]; // standardized group-level effects
vector<lower=0>[M_14] sd_14; // group-level standard deviations
vector[N_14] z_14[M_14]; // standardized group-level effects
vector<lower=0>[M_15] sd_15; // group-level standard deviations
vector[N_15] z_15[M_15]; // standardized group-level effects
vector<lower=0>[M_16] sd_16; // group-level standard deviations
vector[N_16] z_16[M_16]; // standardized group-level effects
}
transformed parameters {
vector[N_1] r_1_muBuilding_1; // actual group-level effects
vector[N_2] r_2_muBuilding_1; // actual group-level effects
vector[N_3] r_3_muCaregiver_1; // actual group-level effects
vector[N_4] r_4_muCaregiver_1; // actual group-level effects
vector[N_5] r_5_muPatient_1; // actual group-level effects
vector[N_6] r_6_muPatient_1; // actual group-level effects
vector[N_7] r_7_muPatientpsych_1; // actual group-level effects
vector[N_8] r_8_muPatientpsych_1; // actual group-level effects
vector[N_9] r_9_muSDOH_1; // actual group-level effects
vector[N_10] r_10_muSDOH_1; // actual group-level effects
vector[N_11] r_11_muCanmeds_1; // actual group-level effects
vector[N_12] r_12_muCanmeds_1; // actual group-level effects
vector[N_13] r_13_muExpertise_1; // actual group-level effects
vector[N_14] r_14_muExpertise_1; // actual group-level effects
vector[N_15] r_15_muRecommendations_1; // actual group-level effects
vector[N_16] r_16_muRecommendations_1; // actual group-level effects
r_1_muBuilding_1 = (sd_1[1] * (z_1[1]));
r_2_muBuilding_1 = (sd_2[1] * (z_2[1]));
r_3_muCaregiver_1 = (sd_3[1] * (z_3[1]));
r_4_muCaregiver_1 = (sd_4[1] * (z_4[1]));
r_5_muPatient_1 = (sd_5[1] * (z_5[1]));
r_6_muPatient_1 = (sd_6[1] * (z_6[1]));
r_7_muPatientpsych_1 = (sd_7[1] * (z_7[1]));
r_8_muPatientpsych_1 = (sd_8[1] * (z_8[1]));
r_9_muSDOH_1 = (sd_9[1] * (z_9[1]));
r_10_muSDOH_1 = (sd_10[1] * (z_10[1]));
r_11_muCanmeds_1 = (sd_11[1] * (z_11[1]));
r_12_muCanmeds_1 = (sd_12[1] * (z_12[1]));
r_13_muExpertise_1 = (sd_13[1] * (z_13[1]));
r_14_muExpertise_1 = (sd_14[1] * (z_14[1]));
r_15_muRecommendations_1 = (sd_15[1] * (z_15[1]));
r_16_muRecommendations_1 = (sd_16[1] * (z_16[1]));
}
model {
// initialize linear predictor term
vector[N] muBuilding = Intercept_muBuilding + Xc_muBuilding * b_muBuilding;
// initialize linear predictor term
vector[N] muCaregiver = Intercept_muCaregiver + Xc_muCaregiver * b_muCaregiver;
// initialize linear predictor term
vector[N] muPatient = Intercept_muPatient + Xc_muPatient * b_muPatient;
// initialize linear predictor term
vector[N] muPatientpsych = Intercept_muPatientpsych + Xc_muPatientpsych * b_muPatientpsych;
// initialize linear predictor term
vector[N] muSDOH = Intercept_muSDOH + Xc_muSDOH * b_muSDOH;
// initialize linear predictor term
vector[N] muCanmeds = Intercept_muCanmeds + Xc_muCanmeds * b_muCanmeds;
// initialize linear predictor term
vector[N] muExpertise = Intercept_muExpertise + Xc_muExpertise * b_muExpertise;
// initialize linear predictor term
vector[N] muRecommendations = Intercept_muRecommendations + Xc_muRecommendations * b_muRecommendations;
// linear predictor matrix
vector[ncat] mu[N];
for (n in 1:N) {
// add more terms to the linear predictor
muBuilding[n] += r_1_muBuilding_1[J_1[n]] * Z_1_muBuilding_1[n] + r_2_muBuilding_1[J_2[n]] * Z_2_muBuilding_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muCaregiver[n] += r_3_muCaregiver_1[J_3[n]] * Z_3_muCaregiver_1[n] + r_4_muCaregiver_1[J_4[n]] * Z_4_muCaregiver_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muPatient[n] += r_5_muPatient_1[J_5[n]] * Z_5_muPatient_1[n] + r_6_muPatient_1[J_6[n]] * Z_6_muPatient_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muPatientpsych[n] += r_7_muPatientpsych_1[J_7[n]] * Z_7_muPatientpsych_1[n] + r_8_muPatientpsych_1[J_8[n]] * Z_8_muPatientpsych_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muSDOH[n] += r_9_muSDOH_1[J_9[n]] * Z_9_muSDOH_1[n] + r_10_muSDOH_1[J_10[n]] * Z_10_muSDOH_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muCanmeds[n] += r_11_muCanmeds_1[J_11[n]] * Z_11_muCanmeds_1[n] + r_12_muCanmeds_1[J_12[n]] * Z_12_muCanmeds_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muExpertise[n] += r_13_muExpertise_1[J_13[n]] * Z_13_muExpertise_1[n] + r_14_muExpertise_1[J_14[n]] * Z_14_muExpertise_1[n];
}
for (n in 1:N) {
// add more terms to the linear predictor
muRecommendations[n] += r_15_muRecommendations_1[J_15[n]] * Z_15_muRecommendations_1[n] + r_16_muRecommendations_1[J_16[n]] * Z_16_muRecommendations_1[n];
}
for (n in 1:N) {
mu[n] = [0, muBuilding[n], muCaregiver[n], muPatient[n], muPatientpsych[n], muSDOH[n], muCanmeds[n], muExpertise[n], muRecommendations[n]]';
}
// priors including all constants
target += normal_lpdf(b_muBuilding[1] | 0,5);
target += normal_lpdf(b_muBuilding[2] | 0,5);
target += normal_lpdf(b_muBuilding[3] | 0,5);
target += normal_lpdf(b_muBuilding[4] | 0,5);
target += normal_lpdf(b_muBuilding[5] | 0,5);
target += normal_lpdf(b_muBuilding[6] | 0,5);
target += normal_lpdf(b_muBuilding[7] | 0,5);
target += normal_lpdf(b_muBuilding[8] | 0,5);
target += normal_lpdf(b_muBuilding[9] | 0,5);
target += normal_lpdf(b_muBuilding[10] | 0,5);
target += normal_lpdf(b_muBuilding[11] | 0,5);
target += normal_lpdf(b_muBuilding[12] | 0,5);
target += normal_lpdf(b_muBuilding[13] | 0,5);
target += normal_lpdf(b_muBuilding[14] | 0,5);
target += normal_lpdf(b_muCaregiver[1] | 0,5);
target += normal_lpdf(b_muCaregiver[2] | 0,5);
target += normal_lpdf(b_muCaregiver[3] | 0,5);
target += normal_lpdf(b_muCaregiver[4] | 0,5);
target += normal_lpdf(b_muCaregiver[5] | 0,5);
target += normal_lpdf(b_muCaregiver[6] | 0,5);
target += normal_lpdf(b_muCaregiver[7] | 0,5);
target += normal_lpdf(b_muCaregiver[8] | 0,5);
target += normal_lpdf(b_muCaregiver[9] | 0,5);
target += normal_lpdf(b_muCaregiver[10] | 0,5);
target += normal_lpdf(b_muCaregiver[11] | 0,5);
target += normal_lpdf(b_muCaregiver[12] | 0,5);
target += normal_lpdf(b_muCaregiver[13] | 0,5);
target += normal_lpdf(b_muCaregiver[14] | 0,5);
target += normal_lpdf(b_muPatient[1] | 0,5);
target += normal_lpdf(b_muPatient[2] | 0,5);
target += normal_lpdf(b_muPatient[3] | 0,5);
target += normal_lpdf(b_muPatient[4] | 0,5);
target += normal_lpdf(b_muPatient[5] | 0,5);
target += normal_lpdf(b_muPatient[6] | 0,5);
target += normal_lpdf(b_muPatient[7] | 0,5);
target += normal_lpdf(b_muPatient[8] | 0,5);
target += normal_lpdf(b_muPatient[9] | 0,5);
target += normal_lpdf(b_muPatient[10] | 0,5);
target += normal_lpdf(b_muPatient[11] | 0,5);
target += normal_lpdf(b_muPatient[12] | 0,5);
target += normal_lpdf(b_muPatient[13] | 0,5);
target += normal_lpdf(b_muPatient[14] | 0,5);
target += normal_lpdf(b_muPatientpsych[1] | 0,5);
target += normal_lpdf(b_muPatientpsych[2] | 0,5);
target += normal_lpdf(b_muPatientpsych[3] | 0,5);
target += normal_lpdf(b_muPatientpsych[4] | 0,5);
target += normal_lpdf(b_muPatientpsych[5] | 0,5);
target += normal_lpdf(b_muPatientpsych[6] | 0,5);
target += normal_lpdf(b_muPatientpsych[7] | 0,5);
target += normal_lpdf(b_muPatientpsych[8] | 0,5);
target += normal_lpdf(b_muPatientpsych[9] | 0,5);
target += normal_lpdf(b_muPatientpsych[10] | 0,5);
target += normal_lpdf(b_muPatientpsych[11] | 0,5);
target += normal_lpdf(b_muPatientpsych[12] | 0,5);
target += normal_lpdf(b_muPatientpsych[13] | 0,5);
target += normal_lpdf(b_muPatientpsych[14] | 0,5);
target += normal_lpdf(b_muSDOH[1] | 0,5);
target += normal_lpdf(b_muSDOH[2] | 0,5);
target += normal_lpdf(b_muSDOH[3] | 0,5);
target += normal_lpdf(b_muSDOH[4] | 0,5);
target += normal_lpdf(b_muSDOH[5] | 0,5);
target += normal_lpdf(b_muSDOH[6] | 0,5);
target += normal_lpdf(b_muSDOH[7] | 0,5);
target += normal_lpdf(b_muSDOH[8] | 0,5);
target += normal_lpdf(b_muSDOH[9] | 0,5);
target += normal_lpdf(b_muSDOH[10] | 0,5);
target += normal_lpdf(b_muSDOH[11] | 0,5);
target += normal_lpdf(b_muSDOH[12] | 0,5);
target += normal_lpdf(b_muSDOH[13] | 0,5);
target += normal_lpdf(b_muSDOH[14] | 0,5);
target += normal_lpdf(b_muCanmeds[1] | 0,5);
target += normal_lpdf(b_muCanmeds[2] | 0,5);
target += normal_lpdf(b_muCanmeds[3] | 0,5);
target += normal_lpdf(b_muCanmeds[4] | 0,5);
target += normal_lpdf(b_muCanmeds[5] | 0,5);
target += normal_lpdf(b_muCanmeds[6] | 0,5);
target += normal_lpdf(b_muCanmeds[7] | 0,5);
target += normal_lpdf(b_muCanmeds[8] | 0,5);
target += normal_lpdf(b_muCanmeds[9] | 0,5);
target += normal_lpdf(b_muCanmeds[10] | 0,5);
target += normal_lpdf(b_muCanmeds[11] | 0,5);
target += normal_lpdf(b_muCanmeds[12] | 0,5);
target += normal_lpdf(b_muCanmeds[13] | 0,5);
target += normal_lpdf(b_muCanmeds[14] | 0,5);
target += normal_lpdf(b_muExpertise[1] | 0,5);
target += normal_lpdf(b_muExpertise[2] | 0,5);
target += normal_lpdf(b_muExpertise[3] | 0,5);
target += normal_lpdf(b_muExpertise[4] | 0,5);
target += normal_lpdf(b_muExpertise[5] | 0,5);
target += normal_lpdf(b_muExpertise[6] | 0,5);
target += normal_lpdf(b_muExpertise[7] | 0,5);
target += normal_lpdf(b_muExpertise[8] | 0,5);
target += normal_lpdf(b_muExpertise[9] | 0,5);
target += normal_lpdf(b_muExpertise[10] | 0,5);
target += normal_lpdf(b_muExpertise[11] | 0,5);
target += normal_lpdf(b_muExpertise[12] | 0,5);
target += normal_lpdf(b_muExpertise[13] | 0,5);
target += normal_lpdf(b_muExpertise[14] | 0,5);
target += normal_lpdf(b_muRecommendations[1] | 0,5);
target += normal_lpdf(b_muRecommendations[2] | 0,5);
target += normal_lpdf(b_muRecommendations[3] | 0,5);
target += normal_lpdf(b_muRecommendations[4] | 0,5);
target += normal_lpdf(b_muRecommendations[5] | 0,5);
target += normal_lpdf(b_muRecommendations[6] | 0,5);
target += normal_lpdf(b_muRecommendations[7] | 0,5);
target += normal_lpdf(b_muRecommendations[8] | 0,5);
target += normal_lpdf(b_muRecommendations[9] | 0,5);
target += normal_lpdf(b_muRecommendations[10] | 0,5);
target += normal_lpdf(b_muRecommendations[11] | 0,5);
target += normal_lpdf(b_muRecommendations[12] | 0,5);
target += normal_lpdf(b_muRecommendations[13] | 0,5);
target += normal_lpdf(b_muRecommendations[14] | 0,5);
target += cauchy_lpdf(sd_1[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_1[1] | 0, 1);
target += cauchy_lpdf(sd_2[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_2[1] | 0, 1);
target += cauchy_lpdf(sd_3[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_3[1] | 0, 1);
target += cauchy_lpdf(sd_4[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_4[1] | 0, 1);
target += cauchy_lpdf(sd_5[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_5[1] | 0, 1);
target += cauchy_lpdf(sd_6[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_6[1] | 0, 1);
target += cauchy_lpdf(sd_7[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_7[1] | 0, 1);
target += cauchy_lpdf(sd_8[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_8[1] | 0, 1);
target += cauchy_lpdf(sd_9[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_9[1] | 0, 1);
target += cauchy_lpdf(sd_10[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_10[1] | 0, 1);
target += cauchy_lpdf(sd_11[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_11[1] | 0, 1);
target += cauchy_lpdf(sd_12[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_12[1] | 0, 1);
target += cauchy_lpdf(sd_13[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_13[1] | 0, 1);
target += cauchy_lpdf(sd_14[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_14[1] | 0, 1);
target += cauchy_lpdf(sd_15[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_15[1] | 0, 1);
target += cauchy_lpdf(sd_16[1] | 0,1)
- 1 * cauchy_lccdf(0 | 0,1);
target += normal_lpdf(z_16[1] | 0, 1);
// likelihood including all constants
if (!prior_only) {
for (n in 1:N) {
target += multinomial_logit_lpmf(Y[n] | mu[n]);
}
}
}
generated quantities {
// actual population-level intercept
real b_muBuilding_Intercept = Intercept_muBuilding - dot_product(means_X_muBuilding, b_muBuilding);
// actual population-level intercept
real b_muCaregiver_Intercept = Intercept_muCaregiver - dot_product(means_X_muCaregiver, b_muCaregiver);
// actual population-level intercept
real b_muPatient_Intercept = Intercept_muPatient - dot_product(means_X_muPatient, b_muPatient);
// actual population-level intercept
real b_muPatientpsych_Intercept = Intercept_muPatientpsych - dot_product(means_X_muPatientpsych, b_muPatientpsych);
// actual population-level intercept
real b_muSDOH_Intercept = Intercept_muSDOH - dot_product(means_X_muSDOH, b_muSDOH);
// actual population-level intercept
real b_muCanmeds_Intercept = Intercept_muCanmeds - dot_product(means_X_muCanmeds, b_muCanmeds);
// actual population-level intercept
real b_muExpertise_Intercept = Intercept_muExpertise - dot_product(means_X_muExpertise, b_muExpertise);
// actual population-level intercept
real b_muRecommendations_Intercept = Intercept_muRecommendations - dot_product(means_X_muRecommendations, b_muRecommendations);
}
I’m slowly attempting to work backwards to learn Stan from brms, but this relatively complex model has me at a loss.
Are there any corrections or optimizations that might help the sample run more efficiently in this case?
I have also attached the data:
What_Wider.csv (331.3 KB)