Appropriate GLMM distribution for ratings data that are bounded and discrete

I am using a linear mixed model to explain variation in an object's ratings. These ratings are bounded between 0 and 10 [EDIT: see note below], and take only discrete values (example histogram of the raw data below). Note that the distance between values can be thought of as consistent across the scale (unlike in a typical Likert scale) i.e. the difference between 1 and 2 is about the same as that between 8 and 9.

I know that a normal distribution can be a useful enough approximation for discrete data in some cases, but I suspect that I need to use a generalized linear model for this.

What sort of distribution/link function would be appropriate for this type of data?

Edit: I should note that I am using a mixed model with random effects to account for non-independence due to multiple measurements by some individuals. My apologies for not including this detail earlier.

Edit 2: the bounds are not 0 and 10 as I stated above, but 0 and 9.

• I suggest binomial family and logit link. In software I know 10 would be specified as the upper limit. Jul 18, 2019 at 14:25
• What about a model that deals explicitly with the ordinal nature of the ratings such as a proportional odds model? Jul 18, 2019 at 14:35
• @mkt Just because the ratings are equally spaced doesn't mean you can treat the rating as metric in my opinion (see here and here). Your comment gives me the impression that you want to avoid an ordinal model if possible and only use it when it is needed. Could you explain why you hope to use another model? Jul 18, 2019 at 15:09
• @COOLSerdash Thank you for those links, I will take a look. I'm open to switching, but you're right that I do have a mild bias towards avoiding ordinal models. This is simply because I'm less familiar with them and therefore less able to judge weaknesses/problems in model diagnostics. But if they are the best option here, I'll gladly use them.
– mkt
Jul 18, 2019 at 15:12
• It does seem odd that 9 is the mode while 10 doesn't occur at all. Jul 18, 2019 at 15:38

If you want a crude representation with a reasonable mean-variance relationship, a quasibinomial model would probably work OK: in R,

glm(cbind(score,10-score) ~ ..., family=quasibinomial, ...)


cbind(score,10-score) specifies to treat the response as a number of "successes" out of a maximum 10. Specifically, this means you'd be fitting a model with a predicted value of

$$\mu = \textrm{logistic}(\beta_0 + \beta_1 x_1 + ...)$$

where $$x_i$$ are your predictor variables and $$\mu$$ is the proportion out of 10. The assumed mean-variance relation is $$V \propto \mu(1-\mu)$$, which is reasonable since you expect the variance to decrease when the average score approaches 0 or the maximum value.

As pointed out in the comments, this does assume that the levels are in some sense "equally spaced" (although the logistic transformation does mean that it will take (for example) less difference in predictors to go from 4 to 5 or 5 to 6 than from 1 to 2 or from 9 to 10 ...) An ordinal model would add an extra $$n-1$$ parameters to quantify the distances between levels.

• Thanks very much, Ben! This makes sense to me. As an aside, would this be possible to implement in GLMMs? My understanding was that glmer & glmmTMB didn't allow quasi-families.
– mkt
Jul 19, 2019 at 7:54
• It appears to be implementable in glmmPQL, though I'm less familiar with that.
– mkt
Jul 19, 2019 at 8:45
• there's code in the GLMM FAQ that shows how to convert a summary table from (e.g.) binomial to quasibinomial by adjusting standard errors and p-values appropriately. Jul 19, 2019 at 12:43

You could also treat your outcome as an ordinal variable and fit a mixed-effects model for this ordinal response. Ordinal mixed models are available in the GLMMadaptive package; for an example, check the vignette Mixed Models for Ordinal Data.

Alternatively, and if you do not have any missing data in your outcome (or the missing data that you may have can be assumed to be missing completely at random), you could also consider using generalized estimating equations (GEEs), either with a binomial family or an ordinal one. These are available in the geepack in the functions geeglm(..., family = binomial()) and ordgee(), respectively. The GEEs work under quasi-likelihood and include an over-dispersion parameter automatically.

• (+1) Thank you for the suggestions! Much appreciated, I will look into those.
– mkt
Jul 19, 2019 at 9:25