Good point. Two modes are bad news and probably mean that you need the tighter priors or respecify the model to avoid the modes (do the two modes result in identical inferences?)

As for modelling the correlation - given the non-identifiability, this is unlikely to help much, but I wrote it before seeing your previous post, so I am posting it anyways :-)

You would need to figure out how do the sigmas bind together - which requires understanding the data generating process (which I do not). Two guesses below, might not actually work…

A) if you can assume (as you suggested) that sigmaA^2+sigmaW^2=unknown const you can model this directly with something like:

```
parameters {
simplex[2] sigmas_raw;
real<lower=0> sigmas_sum;
}
transformed data {
real sigmaA = sigmas_raw[1] * sigmas_sum;
real sigmaW = sigmas_raw[2] * sigmas_sum;
}
```

+ you would need a prior on sigmas_sum

B) Model the correlation explicitly (more of a pseudocode - I’ve never done this and don’t understand it deeply, it’s just something that could work):

```
real sigmaA_mean ~ [original prior for sigmaA];
real sigmaW_mean ~ [original prior for sigmaW];
cov_matrix sigma_cov ~ inv_wishart( [expected correlation] );
vector[2] sigmas ~ multi_normal({sigmaA, sigmaW}, sigma_cov);
```