title: Annual count index from GAM looking at long-term trends by site
body: <p>I'm interested in estimating a shared, global trend over time for counts monitored at several different sites using generalized additive models (gams). I've read this great introduction to hierarchical gams (hgams) by <a href="https://peerj.com/articles/6876/" rel="nofollow noreferrer">Pederson et al. (2019)</a>, and I believe I can setup the model as follows (the Pederson et al. (2019) GS model),</p>
  <pre><code>fit_model = gam(count ~ s(year, m = 2) + s(year, site, bs = 'fs', m = 2),
                  data = count_df,
                  family = nb(link = 'log'),
                  method = 'REML')
  </code></pre>
  <p>I can plot the partial effect smooths, look at the fit diagnostics, and everything looks reasonable. My question is how to extract a non-centered annual relative count index? My first thought would be to add the estimated intercept (the average count across sites at the beginning of the time series) to the <code>s(year)</code> smooth (the shared global smooth). But I'm not sure if the uncertainty around that smooth already incorporates uncertainty in the estimated intercept? Or if I need to add that in? All of this was possible thanks to the amazing R libraries <code>mgcv</code>, <code>gratia</code>, and <code>dplyr</code>.</p>