I'd like to know if a method I'm trying to use for analysis is valid (statistically speaking).

Here's the deal :

My dataset has a few quantitative variables and I'm trying to see if a qualitative variable has a significant effect on these, as well as determine which category has significant effect and quantify it.

To that end, I'm using the lm() function in R, with quant-mean(quant)~0+qual as the given equation. The 0+ gives me a model with no intercept to avoid having first category as reference, which allows me to have the summary() give all coefficients compared to 0 instead of a comparison to the reference class. The quant-mean(quant) is to have the model coefficients correspond to the difference between the mean of their class and the global mean.

To start the analysis, I use the anova() function on the model to see if a significant portion of the Sum Sq. is explained by the qualitative variable. If I'm not wrong, this function should tell me of a significant effect if at least one of the categories significantly differs from the mean. Then, if need be I'll use the summary to see which category(ies) are causing this and in what way.

The trouble is that with the "no intercept" model, I lose a degree of freedom, so a significant effect as a whole may not be detected through an anova of the proposed lm, versus one of the basic quant~qual lm.

It doesn't seem right to build the "basic" lm to check for effect significance and then build the other model to identify which categories cause this, but I'm not sure if I should be taking the risk of directly building the "no intercept" model... Any suggestions ?

  • $\begingroup$ Why don't you do fit <- aov(quant ~ qual); summary(fit); TukeyHSD(fit)? I don't understand why you add the mean as an offset (usually an offset should be on the left-hand side). Also, "which category has significant effect" needs a reference. Right now you use the mean as a reference which seems strange. $\endgroup$ – Roland May 24 '18 at 7:41
  • $\begingroup$ Yes, I'm trying to identify if some categories cause the quantitative variable to be higher or lower than the rest (as a whole). That's why I'm using the mean as reference. The whole point is to be able to say that "category A has an average <quant. variable> which is <coef.> lower than the rest". $\endgroup$ – RoB May 24 '18 at 8:14
  • $\begingroup$ I think you should simply use "sum to zero" contrasts: atyre2.github.io/2016/09/03/sum-to-zero-contrasts.html $\endgroup$ – Roland May 24 '18 at 8:40
  • $\begingroup$ Thanks for the info, although I fail to see how this is possible with a model using only a qualitative variable ? $\endgroup$ – RoB May 24 '18 at 8:53
  • $\begingroup$ Contrast only concern qualitative variables. So, I don't understand why you think there would be an issue? $\endgroup$ – Roland May 24 '18 at 9:10

Well, I now realize my approach was wrong (thanks @Roland for pointing me in the right direction). To obtain the desired results ("Categ. A is lower than mean..."), the "sum to zero" contrasts were indeed the right way to go. What I failed to realize at the time was that the intercept of such a model was the general mean, and that we lose a category in the model's summary.

To get the last coefficient, one would then simply modify the contrast matrix appropriately.

Since I wrote a bit of code to do this, I'll post it here (hopefully, it helps someone)

k = length(levels(qual))

cont1 = contr.sum(k)
colnames(cont1) = (levels(qual))[1:k-1]   #so summary stays readable

cont2 = rbind(cont1[4,], cont1[-nrow(cont1), ]) # shift matrix down a line
colnames(cont2) = (levels(catrefs2))[2:k]

contrasts(qual) <- cont1 
g1 = glm(quant~qual)
summary(g1) # gives comparisons to mean for coefs 1 to k-1

contrasts(qual) <- cont2 
g2 = glm(quant~qual)
summary(g2) # gives comparisons for coefs 2 to k
| cite | improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.