Thank you all for your answers! I’m still a bit puzzled about why some \beta parameters have high Rhat, since this repeatedly occurs when moving from yearly to monthly data. I initially thought it was due to the increase in the number of points at which the GP is evaluated, but another model contradicted this (see this). Even with the same model structure, i.e., not accounting for within-year/month variation, high Rhat appears as soon as I use monthly data.
@Bob_Carpenter : I don’t think scalable spectral methods apply here, since, if I understand correctly, they don’t provide uncertainty intervals.
@ruarai I tried different priors, but this didn’t resolve the treedepth issue. Adapting the prior for the lengthscale makes sense and is something I will keep in mind for other variants. Since this model is already simplified, I want to fully understand its behaviour before adding more complexity.
@avehtari and @Niko : Thanks for pointing out partial non-centered parameterizations, I wasn’t aware of this. The blog post looks very interesting, but I’m not sure it would solve my issue, as I don’t see funnel plots or strong correlations between the \beta parameters and either the lengthscale or the marginal standard deviation.