Effect of set_rescor(TRUE) on multivariate models with horseshoe prior

Interfaces brms
1

How does estimating residual correlations affect model fitting?

Does setting set_rescor(TRUE) change how effects are fit? In a simple multivariate model, it has no effect (left-side models in output plots).

However, with a horseshoe prior on b, set_rescor(TRUE) has a large impact on the modelled effects (right-side models in output plots).

  • Operating System: Ubuntu 20.10
  • brms Version: version 2.15.0

require(brms)
require(ggpubr)

n <- 18

set.seed(21)

dat <- tibble( Group = rep(rep(0:1,each=n/2),2),
               res1 = rnorm( 2*n, 0, 1 ),
               res2 = rnorm( 2*n, 0, 1 ),
               res3 = rnorm( 2*n, 0, 1 ),
               res4 = rnorm( 2*n, 0, 1 ) ) %>%
  mutate( y1 = Group + res1,
          y2 = Group + y1 + res2,
          y3 = y2 + res3,
          y4 = y1 + y2 + y3 + res4,
          
          y1 = (y1 - mean(y1))/sd(y1),
          y2 = (y2 - mean(y2))/sd(y2),
          y3 = (y3 - mean(y3))/sd(y3),
          y4 = (y4 - mean(y4))/sd(y4),
          
          Group = factor(Group) ) %>%
  as.data.frame

#************************************************************************************
#* USING BRMS
#************************************************************************************

bf0 <- bf( mvbind(y1,y2,y3,y4) ~ Group) + set_rescor(FALSE)
bf1 <- bf( mvbind(y1,y2,y3,y4) ~ Group) + set_rescor(TRUE)

b_priors <- c( set_prior('horseshoe(1)', class='b', resp='y1'),
               set_prior('horseshoe(1)', class='b', resp='y2'),
               set_prior('horseshoe(1)', class='b', resp='y3'),
               set_prior('horseshoe(1)', class='b', resp='y4') )

brms_0_0 <- brm( bf0,
                 control = list(adapt_delta=0.99),
                 cores = 4,
                 seed = 1,
                 data = dat )

brms_1_0 <- brm( bf1,
                 control = list(adapt_delta=0.99),
                 cores = 4,
                 seed = 1,
                 data = dat )

# regularization
brms_0_1 <- brm( bf0,
                 prior = b_priors,
                 control = list(adapt_delta=0.999),
                 cores = 4,
                 seed = 1,
                 data = dat )

brms_1_1 <- brm( bf1,
                 prior = b_priors,
                 control = list(adapt_delta=0.999),
                 cores = 4,
                 seed = 1,
                 data = dat )

ggarrange( mcmc_plot(brms_0_0) + coord_cartesian(xlim=c(-1.5,1.5)),
           mcmc_plot(brms_0_1) + coord_cartesian(xlim=c(-1.5,1.5)),
           mcmc_plot(brms_1_0) + coord_cartesian(xlim=c(-1.5,1.5)), 
           mcmc_plot(brms_1_1) + coord_cartesian(xlim=c(-1.5,1.5)), ncol=2, nrow=2)
1 Like
2

I did not initially check, but the interaction of set_rescor(TRUE) with priors does not require the horseshoe. Simply using e.g. b ~ normal(0,1) does the same. So the effect of set_rescor(TRUE) depends on whether or not priors are set on b.

1 Like
3

Hi, sorry for not getting to you earlier.

Stan sampler (and hence brms) uses always the same method to fit the effects (Hamiltonian Monte-Carlo). However, adding residual correlation changes the model. In full generality even small changes to a model can sometimes alter the estimates for all parameters. On the other hand, sometimes a very different model can yield exactly the same estimates.

Here’s the image produced on my computer with your code (with horseshoe prior):

obrazek
obrazek741×856 11.4 KB

There indeed is a slight difference when adding both residual correlations and priors that is bigger then when adding only one of those. However, when I change the seed to 5555 (for data generation and fits) I get almost no difference:

obrazek
obrazek741×856 11.1 KB

Similarly with a different seed (2164655)

obrazek
obrazek741×856 11.3 KB

To me this looks like fitting residual correlation has very little impact on your model (at least for this particular way of simulating data). Additionally, I can see how with such a small dataset, priors can influence the fit noticeablye. In particular, when you fit the residual correlation, you usually slightly decrease the amount of information influencing the other parameters (as you now have more parameters to fit). I find it plausible that this can let the zero-centered prior have a somewhat bigger influence on the posterior than in the fit without priors.

I don’t use residual correlations too often, but when I did use them, the difference from an uncorrelated also was not very big, so I am not completely surprised. But I admit I don’t understand the underlying mechanisms very well.

Hope that helps at least a little bit.

1 Like