I have reason to believe that different treatments have influence on the variance of replicate measurements within a treatment. However I can’t figure out how my brm model formula should be in order to take this difference in variances into account.

To make things maybe a bit clearer I have simulated data which can be seen in below figure. This data is normally distributed with muA = 1, muB = 1.3 and muC = 0.8; the sd for all groups was 0.1. The sd for duplicate samples varies by group and were respectively sd_dupA = 0.05, sd_dupB = 0.01 and sd_dupC = 0.09.

My data has three treatment groups (A,B and C), each treatment group has 6 samples and each sample is measured twice (Dupl). I have tried the following formula which do not seem to give me estimates on the different sd’s for the duplicates.

Edit: Correction to previous post - you’re right, if you want within Sample variation (i.e. between duplicates) then you don’t want (1|Sample/Dupl), but I think you do want bf(Value ~ Treatment + (1|Sample), sigma ~ Treatment). You should be able to obtain different sigma for each level of Treatment, where sigma is the within Sample, i.e. between Dupl, variability.

Hey thanks for your help, I think this is working! However I have a hard time interpreting the results. It seems that I have to exponentiate the sigma_TreatmentN parameters in order to get them on the right scale but I can’t get my head around of why that is. Could you clarify this?

Family: gaussian
Links: mu = identity; sigma = log
Formula: Value ~ Treatment + (1 | Sample)
sigma ~ Treatment
Data: d (Number of observations: 3000)
Draws: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
total post-warmup draws = 4000
Group-Level Effects:
~Sample (Number of levels: 1500)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.99 0.03 0.94 1.05 1.01 927 1681
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept -0.04 0.06 -0.15 0.07 1.00 2351 3139
sigma_Intercept 0.05 0.03 -0.01 0.11 1.00 2564 2718
Treatment2 -1.03 0.07 -1.17 -0.88 1.00 831 1810
Treatment3 1.62 0.11 1.40 1.84 1.00 4006 3139
sigma_Treatment2 -1.05 0.04 -1.13 -0.96 1.00 2356 3045
sigma_Treatment3 0.98 0.04 0.90 1.05 1.00 3380 3307
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Maybe I can jump in to address this quickly:
In order to keep \sigma positive all the computation is done on the unconstrained \log-scale. This means that also the coefficients will be on the \log-scale (\log(\sigma) = \sigma_{intercept}+\sigma_{Treatment_{k}}) . For your interpretation it means that you have to think about sigma_Intercept, sigma_Treatment2 and sigma_Treatment3 as multiplicative effects, because summing them up would be the equivalence of \sigma = e^{\sigma_{intercept}} \cdot e^{\sigma_{Treatment_{k}}}. Let’s say for Treatment 2 it would mean you have a \sigma of e^{0.05}\cdot e^{-1.05} = 0.368.