[R-sig-ME] Comparison between the ouputs of the lmer and the brm (from the brms package) functions

[Vincent Bremhorst]

Dear all,

I fitted the same mixed model using two different functions : lmer and brm.

The estimation of the standard deviation of the random effect and the estimation of the standard errors of the intercept differ. Both estimates are higher with the Bayesian procedure.
Since  I use non-informative prior in the brm specification, I would expect similar results.
The other estimates are similar for both procedures.

Do you have any idea what's happen here?
Thanks for your help,
Vincent Bremhorst.

lmer (model assumptions are met):

res <- lmer(dTmeanoff ~ habitat + (1|week), data=trh)


Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']

Formula: dTmeanoff ~ habitat + (1 | week)

   Data: trh



REML criterion at convergence: 312.5



Scaled residuals:

    Min      1Q  Median      3Q     Max

-3.5949 -0.5149 -0.0181  0.4792  2.2995



Random effects:

 Groups   Name        Variance Std.Dev.

 week     (Intercept) 0.06142  0.2478

 Residual             1.93221  1.3900

Number of obs: 89, groups:  week, 4



Fixed effects:

            Estimate Std. Error      df t value Pr(>|t|)

(Intercept)   0.6200     0.2788 12.5997   2.224   0.0451 *

habitatu     -0.8101     0.3503 83.0542  -2.312   0.0232 *

habitatw     -1.6366     0.3698 83.2151  -4.425  2.9e-05 ***

---

Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1



Correlation of Fixed Effects:

         (Intr) habitt

habitatu -0.638

habitatw -0.605  0.481

brm (convergence of the posterior chains ok)


> priors<- c(prior(normal(0,100), class="b"),

+            prior(normal(0, 100), class="Intercept"),

+            prior(exponential(0.1), class="sigma"),

+            prior(exponential(0.1), class="sd", group= "week")

+

+ )



> fit1<- brm(dTmeanoff~habitat+(1|week), data=trh,

+            prior=priors,

+            iter=4000, warmup=2000,chains=2,

+            family = gaussian(),#no logit function applied

+            file = "output.rds",

+            sample_prior = "yes",

+            control = list(adapt_delta = .9))



Family: gaussian

  Links: mu = identity; sigma = identity

Formula: dTmeanoff ~ habitat + (1 | week)

   Data: trh (Number of observations: 89)

Samples: 2 chains, each with iter = 4000; warmup = 2000; thin = 1;

         total post-warmup samples = 4000



Group-Level Effects:

~week (Number of levels: 4)

              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS

sd(Intercept)     0.50      0.48     0.02     1.82 1.00      900     1048



Population-Level Effects:

          Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS

Intercept     0.63      0.36    -0.11     1.40 1.00     1159      901

habitatu     -0.82      0.35    -1.51    -0.13 1.00     2235     2581

habitatw     -1.65      0.37    -2.37    -0.89 1.00     2091     2500



Family Specific Parameters:

      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS

sigma     1.41      0.11     1.21     1.65 1.00     2863     2199



	[[alternative HTML version deleted]]