How to implement a mixed-model with a beta distribution? I am interested in using a generalised linear mixed model with a response variable (values ranging from 0.001-0.999) that best fits a beta distribution when checked using the 'fitdistrplus' package and the 'descdist()' function in R.

I understand that with mixed-models it is the distribution of the residuals that is important, not the distribution of the raw data points. A QQ-plot of the model residuals also indicates non-normality.

I'm hesitant to transform the data. From what I have seen, the beta distribution is not listed as an option in 'glmer'. I also believe that the package 'betareg' does not allow inclusion of random effects. I've come across two other packages: 'GLMMadaptive' and 'glmmTMB', both of which permit specification of 'family=beta' in the model. The manuals for these packages are very high level and I'm struggling to understand whether either of these packages, or an alternative, would help me to run a mixed-model on this non-normally distributed data.
Can anyone advise about the use of these packages ('GLMMadaptive' or 'glmTMB') or recommend a package or code for a generalised linear mixed model that allows a beta distribution?
Thanks very much in advance!

Please note that there is no requirement, condition, or assumption regarding the distribution of the variables in any regression model. When the data are strictly positive and bounded then the beta distribution is often a very good choice.
GLMMadaptive and glmmTMB both allow for the beta distribution. Since you seem to be familiar with glmer then glmmTMB would be the easist choice for you since all you have to do is specify family = beta_family()
As for the residuals, since it's a beta model there is not expectation that the residuals would be normally distributed. The DHARMa package has some good functionality for assessing the residuals from a beta model.
