Dealing with a large number of random effects in gam/bam I have large datasets and want to model a binomial outcome using gam or bam, but because of the large number of random effects in my datasets, r has memory errors. I'm looking for advice on how I might run my models using less memory.
Example: dataset with n = 60451 observations; varOUTCOME is binomial, var1 and var2 are categorical, var3 is linear numeric and var4RE is the integer variable for random effects. In this set n = 51660 participants, indicating about 10k observations are from repeat participants (participants are present in data 1 to 8 times). We know the ICC is high enough that we need to account for this repeating.
We have tried gam and bam with varOUTCOME ~ var1 + var2 + s(var3, by var1) + s(var4RE, bs='re').
We are now trying this but it also requires hours to run:
gamm4(varOUTCOME ~ var1 + var2  + s(var3, by var1) + s(var4RE, bs='re') + random = ~(1|var4RE),
data = thedata, verbose = TRUE, family = binomial(link='logit'))
As an added note, using the less optimal s(var3) in lieu of the interaction with var1 does not help.
I have searched for/read posts on 'bam, gam and random effects' but not found a solution. Apologies if I missed something.
Thank you in advance for any suggestions on how to handle the large number of random effects.

The default s() smooth in mgcv, the thin-plate regression spline, can be computationally expensive with very large numbers of observations. The default basis dimension might just be too big in your case. Quoting from the manual for gam():

The choice of the basis dimension (k in the s, te, ti and t2 terms) is something that should be considered carefully (the exact value is not critical, but it is important not to make it restrictively small, nor very large and computationally costly).

Specifying a smaller basis dimension might solve your problem with the first model. Your second model is double-modeling the random effect, as you include the same random effect both in a smooth and in an lme4-type random effect. It thus doesn't seem surprising that it never converges. Choose one way to handle the random effect, not both.
You should be able to accomplish what you need without using this type of GAM. You have only 1 continuous predictor (interacting with a categorical predictor) and 1 random effect. Adequately flexible modeling of the continuous predictor might be done with restricted cubic splines (e.g., the ns() function in the R splines package, part of the standard R distribution, or the rcs() function in the rms package). You specify the flexibility of the model via the number of "knots" between which you fit a set of cubic polynomials, constrained to fit smoothly at each knot location and to extend linearly beyond the outermost knots. Unless there's a lot of wiggliness that you need to capture in the model, something on the order of 5 knots often suffices; see Section 2.4 of Frank Harrell's course notes or book. Then you use glmer() from lme4 to handle the random effect.
If you really want to use gam(), you can specify cubic splines in the smooth specification.
