# How to interpret estimated posterior distribution of hyperparameters returned by brms

Hi, I am fitting a simple model with random effect as shown below:

`alpha_i` in the model represents random effect for different genes (`i` here represents genes). I assume `alpha_i` follows a normal distribution with a hierarchical prior `sigma_gene`. The ultimate goal is to estimate the posterior distribution of `sigma_gene`. Other variables are not related to this question so I skipped the introduction of them.

Here is how I specified prior for this model:

``````bpriors <- prior(normal(0, sigma_global_sd), class = "b", coef = "Intercept") +
prior(constant(1), class = "b", coef = "sc.prop") +
prior(normal(0, sigma_gene_sd), class = "sd") +
prior("target += cauchy_lpdf(sigma_global_sd | 0, 10)", check = FALSE) +
prior("target += cauchy_lpdf(sigma_gene_sd | 0, 10)", check = FALSE) +
prior(gamma(0.01, 0.01), class = "shape") +
prior(beta(1, 1), class = "zi")

stanvars <- stanvar(scode = "real<lower=0> sigma_global_sd; real<lower=0> sigma_gene_sd;",
block = "parameters")
``````

`sigma_gene_sd` in the code represents `sigma_gene` in the model above. Below is the brms result after fitting the model:

Based on the result and the data I have, I am pretty sure that the distribution of `sd_gene__Intercept` represents the estimate posterior distribution of the hyperparameter `sigma_gene` in the model, i.e., the posterior distribution of the variance of the normal distribution of the random effect. Then my first question is, if this is correct, then what is the distribution for `sigma_gene_sd`? What is the difference between `sd_gene__Intercept` and `sigma_gene_sd`?

Another question is that, if I run `prior_summary` on the model, here are the list of parameters and their corresponding priors:

What is the difference between the three parameters in the red square? Are they three different parameters or they represent the same thing?

Thank you!

1 Like

Hi @paul.buerkner, not sure if you are the person to ask but since no one answers this question, I was wondering if I could get some directions from you?

The priors all mean the same here as you simply have one parameter. Prior specifiction brms is â€śhierarchicalâ€ť in the way that the prior on the â€śmore generalâ€ť row will be vectorized to more specific cases if the specific cases donâ€™t have priors themselves. For more details, see ?set_prior

With regard to your question about the SDs, I donâ€™t know why you are manually adding these parameters. The one one alpha_i could be added via a random effect `(1 | gene)`. The `sigma_global` as a hyperparameter seems strange to me as you simply have one parameter, `alpha` to inform that one (Imagine trying to estimate a standard deviation one just a single observation).

1 Like

1. I manually adding parameters like `sigma_global_sd` and `sigma_gene_sd` because I saw this post about hierarchical prior specification. I think my situation is the same as the example you showed in that post so I manually specified these hyperparameters. I tried to run with default prior and I am still able to get posterior distribution of `b_Intercept` and `sd_gene__Intercept`. Then my question is, is the result the same between default priors and manually specified hierarchical priors? It seems they are similar at least in my case. If that is true, then is there any advantage of using manually specified hierarchical priors? And for the hierarchical prior case, what is the meaning of the posterior distribution of `sigma_gene_sd` and `sigma_global_sd` in the plot output? It seems that these two are the two things that are not estimated if I use default prior.

2. For your concern of `sigma_global_sd`: I think this is the same as the point above. `alpha` is simply a global intercept. We can specify prior for it and estimate its posterior distribution. Here I simply give it a normal prior with standard deviation as `sigma_global_sd`. I can go without manually specifying it and the posterior distribution will still be estimated. So again, it comes back to the questions above, i.e., what is the difference between default priors and manually specified hierarchical priors.

Can you share your `brms` model syntax/specification? That might make it easier to see how the priors relate to model parameters.

specifying sigma_global_sd doesnâ€™t make sense to me here. it will just follow its prior anyway so I would use a default prior in this case. I think you are making it unnecessarily complicated by estimating this hyperparameter.

you didnâ€™t provide your brm() model specification so I am unsure how that looked but for the gene part right now. I recommend adding (1 | gene) which will add the desired hyperparamter automatically so you donâ€™t need to fiddle with stanvars at all.

Thanks! For `sigma_global_sd`, I will just remove it and use default prior.

Here is how I specify the model:

``````      fit.m <- brm(
bf(sp.count ~ 0 + Intercept + sc.prop + (1 | gene) + offset(log(cs))),
data = data2fit.all,
prior = bpriors, stanvars = stanvars,                                      # this line passes the prior & modified parameter blocks into brm()
control = list(adapt_delta = 0.99, max_treedepth = 15),
save_pars = save_pars(all = TRUE),
family = zero_inflated_negbinomial(),
iter = 2000, chains = 4, cores = 4
)
``````

`bpriors` and `stanvars` are provided at the beginning of my question. I think it is nothing but what you suggested, together with the manually specified priors. For the sake of simplicity, I will probably remove `sigma_global_sd` and `sigma_gene_sd` and use the default priors. Right now I am just trying to understand what difference would the hierarchical prior provide and how to interpret the additional parameters.

I think the purpose of specifying those hierarchical priors explicitly is because I want the priors for both `alpha` and `alpha_i` to be normal distribution with mean 0 and unknown variance. In the default setting, the prior for `alpha` is flat and the prior for `alpha_i` is `student_t(3, 0, 2.5)`. I have seen people simply using `normal(0,10)` as priors for intercept and random effect but I kind of donâ€™t want to specify a fixed standard deviation (10 here). I think this is the main reason why I specify these hierarchical priors instead of using default priors.