# How to translate (relate) **mgcv** `gam()` and **brms** `brm()` results

I am keen to get a detailed grasp of the relationship of some non-linear models analysed with mgcv vs. brms. The syntax is virtually the same, apart from the functions themselves `gam()` vs. `brm()`. But the summary tables differ. To illustrate this, I am using the first toy-example from mgcv:

``````require(mgcv)
require(brms)

# example from mgcv
set.seed(2)
dat <- gamSim(1, n=400, dist='normal', scale=2)

# gam from mgcv
summary(gam1 <- gam(y ~ s(x0) + s(x1) + s(x2) + s(x3), data=dat))
# Parametric coefficients:
#             Estimate Std. Error t value Pr(>|t|)
# (Intercept)  7.83328    0.09878    79.3   <2e-16

# Approximate significance of smooth terms:
#         edf Ref.df      F  p-value
# s(x0) 2.500  3.115  6.921 0.000131
# s(x1) 2.401  2.984 81.858  < 2e-16
# s(x2) 7.698  8.564 88.158  < 2e-16
# s(x3) 1.000  1.000  4.343 0.037818

# brm from brms
summary(brm1 <- brm(y ~ s(x0) + s(x1) + s(x2) + s(x3), data=dat))
# Smooth Terms:
#            Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
# sds(sx0_1)     2.18      1.19     0.74     5.11 1.00     1856     2533
# sds(sx1_1)     2.05      1.11     0.70     4.88 1.00     1928     2416
# sds(sx2_1)    19.74      5.03    12.16    32.31 1.00     1142     1741
# sds(sx3_1)     0.90      0.91     0.03     3.23 1.00     1468     1731

# Population-Level Effects:
#           Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
# Intercept     7.83      0.10     7.64     8.03 1.00     5138     3185
# sx0_1        -1.06      3.77    -8.79     6.21 1.00     2113     2201
# sx1_1         9.31      3.74     1.90    17.20 1.00     2019     2222
# sx2_1        35.47     11.85    11.58    58.22 1.00     2581     2655
# sx3_1        -0.61      2.44    -4.89     5.45 1.00     1881     1648

# Family Specific Parameters:
#       Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
# sigma     1.98      0.07     1.84     2.13 1.00     4643     3062
``````

So, all four smooths are significant according to the `gam`. Equally, by `brm`, CIs do not contain zeros. Nicely, `x3` is flagged weak by both, so to say.

However, how should I interpret `Population-Level Effects` from `brm`? Are these simple comparisons with the Intercept? Is there more to it? Interestingly enough, besides x3 x0 also does not differ from the Intercept.

• Operating System: macOS Monterey (12.5.1)
• brms Version: 2.17.0

Thanks!

@striatum hello, I think the answer is to stop trying to interpret the coefficientsā¦ If you wanted to compare the predictions from each of these models then Iād examine each of the smooth functions directly, thereās a clean shortcut for this from `brms` with `conditional_smooths()` and of course you can show any predictions you want with `fitted()` or `predict()`.

I believe that both of these models are equivalently thin-plate regression splines using `s()`, by default. The default priors in `brms` are not all flat though, so you wouldnāt expect the estimates to be exactly the same, even for the same model.

1 Like

Hey @AWoodward, thanks so much for this! It makes perfect sense. What I like to know/have is a sort of initial ātriageā of the smooths worth further inspection with `conditional_smooths()`. For that, I suppose, CIs should suffice. (By the way, very often, one needs to make reviewers happy rather than present/summarise what is actually meaningful. š š)

See this post and the article by Tristan Mahr linked there. As @AWoodward said, the population-level effects are not very interpretable.
Also, brms treats the smooths as varying effects (aka random effects) and the hyperparameter that controls the wiggliness is like the sd in a varying effects model (that is what those `sds()` are). So this is more like gamm4.