How to check normality assumption in Mixed-effects models? I'm dealing with longitudinal data, and to take into account the dependence of observations within the cluster, I shall rely on a generalized linear mixed model. I have a continuous response variable, and I'd like to fit a Gaussian mixed model. However, plotting the density of the response (even after a log transformation) does not seem to be normal. It has two local maxima (where the second local is also a global maximum).
Is it appropriate to work with a Gaussian model?

There are some issues when incorporating mixed effects, as you have two sources of residual variation, stemming from your level 1 and level 2 effects. It's been a while since I looked into this, but if I recall correctly, there is a debate going on regarding the appropriateness of different types of residuals. I think, Santos Nobre and da Motta Singer (2007) give a good overview over the challenges in modelling as well as show the most commonly used methods.
If you're working in R, I suggest looking into the HLMdiag package. I remember finding it particularly helpful when diagnosing mixed models. For a Bayesian approach, I looked into DHARMa, which might also be worth checking out if it applies to your use case.
