Validity of using a combination of Wilcoxon tests and Spearman rho as an alternative to GLM?

I recently read a paper where the authors used either (a) a generalised linear model or (b) a combination of Wilcoxon rank sum tests and Spearman’s rank correlation, depending on whether the data followed a normal distribution or not. They also used a Bonferroni correction when applicable, to account for multiple comparisons.

The responses were two continuous variables, each analysed separately. They had 3 categorical variables and 1 continuous variable as independent variables. (a) For the first response variable, data followed a normal distribution, so a GLM was fitted assuming a normal distribution with identity link. Inference was (apparently) made by significance of model terms, i.e. whether or not the independent variable was significant in the model or not. (b) For the second response variable, the data did not follow a normal distribution, so a Wilcoxon rank sum tests and a Spearman’s rank correlation test were used instead.

My questions are:

1) Are these two approaches appropriate and equally valid alternatives, depending on (non)normality?

2) Even in case of non-normality, wouldn't a different GLM type be better than a combination of Wilcoxon rank sum tests and a Spearman’s rank correlation? (Not sure which type of GLM would be appropriate though, as the response is still continuous...gamma perhaps?)

3) In a GLM framework, aren't we worried about non-normality of residuals, rather than raw data (which is not always the same thing)?

4) Doesn't a GLM (generalized linear model) with "normal distribution and identity link" simply mean a general linear model (i.e. an ordinary linear model)?

5) Is this approach of using p-values from a full model (including all variables) appropriate, or would a variable selection (using AIC or similar) be more appropriate?

• I have a hard time seeing how these approaches even address the same problem. A GLM is a regression model (of a conditional response) whereas the Wilcoxon test is a particular hypothesis test of a difference in location and Spearman rank correlation is a descriptive statistic of co-variation! It's a little like asking whether riding a bike or going to the beach are valid alternatives to eating an apple. – whuber Dec 7 '18 at 20:27
• Can you provide a link to the paper? – Chris Dec 7 '18 at 21:01
• In addition to the good answers: doing different tests based on outcome of normality tests or naively doing a test on a model buildt using some model selection procedure (e.g. AIC based) invalidates p-values. Using a single pre-specified full model is by far the easiest way to get valid p-values assuming that your model assumptions are correct. If some model selection (whether normality testing, deciding what variables to include etc.) was done, the procedure for getting p-values needs to take that into account. It's possible to do, but requires extra contortions. – Björn Dec 10 '18 at 10:15
• @whuber, yes, and that is where my confusion comes from. I was surprised by this choice of "If normal, let's do a GLM, if non-normal, let's do the Wilcoxon instead", which is why I'm asking this question. Chris, I considered that, but I did not wish to publicly call out the authors (also because I know some of them), but rather inquired about the approach itself. – Tilen Dec 11 '18 at 14:53
• @Björn, are you saying that to determine which predictors are likely drivers of the response, one would either A) fit a pre-specified full model and judge "importance" (yes or no) of predictors based on p-values, or B) use something like AIC to determine which of these predictors of the full model should stay in, but not C) both? Did I get that right? If that is the case, would you opt for A or B to assess which variables seem "important" for the response? – Tilen Dec 15 '18 at 11:51

To clarify terminology, a GLM with Gaussian distribution and identity link is identical to the general linear model or linear regression. So I will simply refer to it as linear regression going forward.

Additionally, the normality assumption in regression is on the error term. Since the residuals are estimated errors, we verify this assumption using the residuals or variants of the residuals, like studentized deleted residuals. Claiming the distributional assumption applies to the error term works in linear regression. However, for other GLMs like logistic regression, we do not assume the error term is binomial. A more general assumption for linear regression is that the data are normally distributed with a mean that depends on the predictors and a variance that is constant. Since the mean of the distribution depends on the predictors, to verify normality, we extract these means, hence giving all of the data the same center/location. This is what residuals are and then we can plot these residuals to check normality. Beyond linear regression, such empirical verification becomes more complicated.

It happens that the normality assumption in linear regression is one of the assumptions we do not need to care much about. The least squares coefficients do not depend on normality. However, classical inference using t-tests relies on the normality assumption. But from central limit theorem, we know this assumption does not matter too much once the sample sizes are large enough.

Rather than moving away from linear regression because of non-normality, I might move to another GLM if I could hypothesize a theoretical reason for why another distribution would be more appropriate. For example, if my data were count so non-negative and they were largely low counts with a few high values, I might start to consider a Poisson distribution. If the data were binary, then Binomial distribution is reasonable. The normal distribution is a reasonable choice if all we are willing to assume is that the data have finite variance and their range is $$(-\infty, \infty)$$. Usually, we can make more assumptions about the data than these two assumptions.

Additionally, the linear regression can handle multiple predictors predicting a single outcome. Wilcoxon/Spearman are for single predictor situations, so they are not comparable to linear regression in this way. They are directly comparable to specialized linear regressions like the Pearson correlation and the t test. However, they ask different questions of the data. Wilcoxon tests are stochastic dominance tests comparing two groups, while Spearman test is a measure of monotonic relation (like an ordinal linear correlation). So even when they are similar to linear regression, they ask different questions of the data.

Finally, with regard to your final question, modeling is hard. Model selection is a difficult problem. A reasonable approach is to attempt different models then observe their varying implications. If their implications are similar, that makes life easy. If they differ, that can be interesting for research purposes. With Bayesian modeling, we have a few more options in that we do not simply have to select a model, we can combine different models to account for our uncertainty in our final model.

Edit

An error on my part from the first paragraph: general linear model also includes MANOVA; this differentiates general linear model from the standard linear regression. All other members of the general linear model family are linear regression models or generalized linear models with gaussian distribution and identity link function.

• With respect to model selection, what I meant was: to determined which of the predictors appears to have an effect on the response, would it be better to a) fit a full model and judge which predictors are "important" based on the resulting p-values, or would I be better of b) judge the "importance" of predictors via something like AIC (assuming my choice of variables is already reduced to a sensible subset, not including everything I can think of)? – Tilen Dec 15 '18 at 11:35
• @Tilen there is no simple answer to this particular question. You would be better off reading some of the good model selection questions on this website as a starting point. Search for model selection, AIC, model averaging, ... – Heteroskedastic Jim Dec 15 '18 at 12:33
• I did in fact, and have asked related questions on this website as well (which can be seen through my profile). Still rather confused though. I'm sure I may have missed something, but I did spend a substantial amount of time on these issues. Many thanks to everyone for the answers provided though. – Tilen Dec 15 '18 at 19:54
• @Tilen I don't think you've missed anything. It is a difficult question with correct answer. Just different perspectives, some more thought out than others. There are wrong answers though. – Heteroskedastic Jim Dec 15 '18 at 19:59

The model that is a generalization of Wilcoxon, Kruskal-Wallis, and, to some extent, Spearman's $$\rho$$, is the proportional odds ordinal logistic semiparametric model. For short we call it the proportional odds (PO) model. The PO model handles the full generality of the linear model plus more, because it is invariant to how $$Y$$ is transformed. This circumvents the need to guess about whether residuals have a normal distribution. The PO model is robust with regard to $$Y$$ and is 0.95 as efficient as the linear model if normality holds. Otherwise it can be more efficient than the linear model.

Software needs to make it fast to handle the case where there is a large number $$k$$ of distinct $$Y$$ values, as the number of intercepts in the PO model is $$k-1$$. The R rms package orm function is fast for $$k$$ up to about 6000.

See my detailed case study in my RMS course notes for modeling continuous $$Y$$ using ordinal regression.

Once you are comfortable with ordinal response models, you no longer even need the special cases such as Wilcoxon. And ironically the PO model handles extreme ties better than the Wilcoxon test.

1) Are these two approaches appropriate and equally valid alternatives, depending on (non)normality?

Definitely not. Regression, in general, seeks to describe the relationship between typical values of a random variable given known information. In their comment, whuber describes the (very different) purpose of the other two tests.

2) Even in case of non-normality, wouldn't a different GLM type be better than a combination of Wilcoxon rank sum tests and a Spearman’s rank correlation? (Not sure which type of GLM would be appropriate though, as the response is still continuous...gamma perhaps?)

Definitely so. We would need more information about the variable to provide specific guidance. Residual plots are a good way to convey the presumed lack of normality.

3) In a GLM framework, aren't we worried about non-normality of residuals, rather than raw data (which is not always the same thing)?

Correct. This is almost never the same thing. Residual plots are a good tool to explore this in linear regression. As Jim mentions, deviations from normality are not in and of themselves a great concern when you have a moderately sized sample and care only about the regression coefficients or mean predicted values.

On the other hand, a more important condition that residual plots also allow you to assess is a mean value of zero for the residuals. If the residuals stray too far at times, it suggests the model may be misspecified.

4) Doesn't a GLM (generalized linear model) with "normal distribution and identity link" simply mean a general linear model (i.e. an ordinary linear model)?

Not quite. A "general linear model" refers to a, well, general version of linear models. This generality also allows for multiple outcomes to be analyzed simultaneously. Generalized linear models only handle univariate outcomes. The normal-identity GLM could be described as "multiple linear regression". That makes it both a particular case of a GLM and of the general linear model.

5) Is this approach of using p-values from a full model (including all variables) appropriate, or would a variable selection (using AIC or similar) be more appropriate?

A p-value only allows you to say whether a statistically significant difference was observed, not that it wasn't, but it fails at communicating anything beyond that. I suggest using the confidence intervals to describe the findings provided by the model, as they show the range of plausible values for the estimated parameters. Even if a variable is significant, a wide confidence interval can suggest that little useful information about that particular relationship was gleaned. Conversely, a non-significant variable with a narrow confidence interval suggests that the presumed relationship is either absent or small in magnitude.

Using AIC for variable selection addresses a very different question. Selecting the AIC-best model will, on average, give you the model that makes the best predictions while also being as small as possible.

• Nicely written reply. I forgot GLM includes MANOVA. A pedantic point: AIC-based model selection is supposed to optimize good out of sample predictions under certain conditions. – Heteroskedastic Jim Dec 8 '18 at 0:52
• Let me see if I get this right. What I thought was the GENERALIZED linear model (which may be Poisson one, binomial, gamma, etc.) is not the same as a GENERAL linear model (which is simply a linear model, ordinary least squares, which may be simple or multiple predictors). I thought that the GENERALIZED linear model (GLM) with normal distribution and identity link function is essentially just a GENERAL linear model (LM, ordinary least squares). Both GLM and LM can have a single or multiple predictors. From Jim's answer, I got a feeling that's the case. But you suggest it's not? – Tilen Dec 15 '18 at 11:40
• By "multiple" outcomes, you mean multiple response variables? So, would a multinomial model be a sub-set of the general linear model then...? – Tilen Dec 15 '18 at 11:42
• Yes, I mean multiple outcomes. The definition put forth in McCullagh's monograph introducing GLMs specifies that the response variable must be a vector. Under the assumptions GLM works under, you can decompose a multinomial response so that each "draw" behaves a single distribution, and if you choose to represent the outcome as a numeric one (eg first option is 1, second is 2 and so on), you can produce a link function that take a single linear predictor. More explanation here. – zipzapboing Dec 15 '18 at 12:13
• As for the General Linear Model, it doesn't allow for a "link" function and directly relates the expectation of the potentially multivariate responses with the linear predictor. If this isn't clear you're not alone in your confusion. As to why these distinctions are necessary, I guess it's a matter of both historical usage and the different techniques that are used to fit each type of model. – zipzapboing Dec 15 '18 at 12:19