Assumptions of linearity and homoskedasticity when the DV is a Likert scale - how to interpret the (optimized) residual plots?

Question

I am currently building mixed models for several analyses where each the DV is ordinal, specifically a 7-point Likert scale. Now, I am fully aware of the controversies of using Likert scale as a DV for a linear model, however after some deliberation, at least for the sake of my current project I have decided to proceed with this.

Of course, this entails problems with interpreting the residual plots, needed for the assumptions of homoskedasticity and linearity, since all the residuals are fitted onto descrete lines, as seen on the plot below:

I have found an interesting way to optimize the residual plot in a way that allows for easier interpretation of Likert scale-derived residuals, that someone posted online. Specifically, they recommend introducing some disturbance to the DV before plotting. I have done so, using this code (in R):

Y1 <- df$Y + runif(N,-.5,.5)

where N is the number of observations.

Having plotted the residuals of the model containing the new Y1, I indeed found the new plots to be much more "typical". Nonetheless, I struggle somewhat with understanding if the assumptions of linearity and homoskedasticity are met, perhaps because the plots are still not as good as the typical plot with a continuous DV, or perhaps simply because I've reached the point where I completely lost any confidence regarding my knowledge in stats. Either way, an input of a another person would be greatly appreciated...

Below I attach the residual plots for two of the 10 models I run. The first one looked the best to my ignorant student I out of all of them, while the second one look the worst, therefore I figured it might be more informative to include both and see what you think.

Plot 1 is below. For this DV, the data was normally distributed, there were no outliers, the mean was pretty much in the middle of the range of values. It seems pretty random to me, however the tilt is still present, and I am simply unsure if this suggests lack of linearity and/or homoskedasticity or not...

And here's Plot 2. The data here does have many issues. The data was not normally distributed (kurtosis value over 9), there were multiple outliers, and the mean was close to 0 -> floor effect. How would you interpret this in terms of the two assumptions? And, what exactly here would indicate heterostedasticity and lack of linearity? Would the very strong floor effect itself be the thing that shows that the assumption of linearity is not met? And, while there's no fanning or curve that I know would indicate heteroskedasticity, can I really assume homoskedasticity?

Here's also some more info about my dataset/analysis, perhaps it'll be useful (perhaps not):

• I am not doing anything about the outliers in my dataset, since they did not look erroneous upon inspecton and exist only due to the floor effect observed, unfortunately, for a majority of my DVs
• The sample is very large, there's 5056 observations per each DV (repeated measures, 79 subjects x 64 stimuli)

Thank you in advance for the input!

Answer 1

First, as far as your "optimization" goes, it is not really "optimization," and further, "there is an app for that": Use "jitter". It will be better than your approach because it does not add as much noise. The purpose of adding noise is not to validate any assumptions; instead, it just allows to to see the data density more clearly when the data are discrete, as in your case.

Second, as far as the assumptions go, every problem you see in the residual plot is a function of the discreteness/boundedness of your DV. And for that reason alone, your linearity and homoscedasticity assumptions are obviously violated, for essentially the same reasons that they are violated for binary DVs. Also obviously, normality is violated as well for the same reasons. There is no reasons to look at kurtosis to tell you that, although the somewhat large kurtosis does alert you to the fact that there are rare, extreme response values in your data.

You can improve the diagnostic value of your plots by adding the LOESS smooth (of the original, not jittered) data to the jittered (predict, resid) plot to better diagnose deviation from linearity. Deviation from a flat line indicates nonlinearity. To better diagnose constant variance, you can add the LOESS smooth (of the original, not jittered) data to the jittered (predict, abs(resid)) plot. Here again, deviation from the flat line indicates a problem, although the problem is heteroscedasticity in this plot.

You might be able to get by with just using ordinary linear regression if these plots look ok, but it will definitely be sub-optimal. In particular, the predicted conditional distributions are normal when you do that, but your distributions are in reality, discrete and highly skewed. So your predictions will be scientifically unrealistic.

Another, probably better option, especially considering that there are plenty of observations available for estimating the additional parameters, is ordinal regression. It also makes assumptions, but you can assess their validity by estimating an even more highly parameterized model such as multinomial logistic regression, then by comparing fit measures such as penalized likelihood, and also by comparing estimated conditional distributions. If the estimated conditional distributions are similar for both models, then the simpler ordinal regression model wins, according to Occam's razor.