Residuals not normal due to random effects structure I am new to mixed-effect models and statistics, and I need to do a project for my thesis using mixed-effect modeling techniques.
My research aims to model how cardiovascular risk scores change over time. In this case, my mixed-effect model is
 lmer(log(Score) ~ Time + I(Time ^2) + (1+ Time + I(Time ^2) | Patients))

I have log transformed the score to improve the model fit. Quadratic function of time shows a better model fit. The random effect is allowing intercept and time to vary among patients.
I got an issue with the residuals of my best-fitted mixed-effect model. Q-Q plot shows the residuals are not normal on the individual level. However, the residuals look normal at the population level. If I fit a linear model (without the random effect component), the residuals also look normal. Hence, I suspected the issue was due to the random effect structure. I've been playing around with different random effect structures and correlations in the lme4 package in R, but neither could fix the residual issue.
I really need help to get through.

 A: it's a good call to see that you can look at residuals marginal (population-level) and conditional (on the RE).
The bad news here though is that what is relevant for your case are the conditional residuals, and they don't look normal. You shouldn't adjust the RE structure to correct this - it's not the fault of the REs, it can just happen that your residual error is non-normal.
In this case, you have at least 3 options to deal with distributional problem

*

*try other transformations

*beta regression

*quantile regressions

I would probably try if beta solves the problem (because this is the natural choice), and else try a quantile regression.
Additional comment: my personal opinion is that the (1+ Time + I(Time ^2) | Patients)) structure creates a weird setup for modelling between-patient variability in temporal trends, because of the dependency between Time and Time^2. I would prefer using factor smooths or factor smooth interactions, see e.g. Multiply a seasonal spline by a region-specific coefficient in Generalized Additive Model (GAM) in R. You can do so keeping your population-level estimate as a polynomial or move to a spline as well.
