I am working on a dataset of vegetation characteristics recorded over three seasons in two different areas used by buffalo.

I have run generalised linear mixed models to determine whether vegetation biomass changes seasonally in different areas used by individual buffalo. I have run every possible combination of model and used AIC to identify the one with the best fit, which includes a fixed effect of season.

Within the best-fitting model, I would like to know which seasons differ from each other, and for this I used a post-hoc Tukey test from the multcomp package to compare pairs of seasons, and identified season pairs that differed significantly based on p-values.

I have been told that using the AIC model selection and the post-hoc tests based on p-values is mixing two types of analysis methods and is therefore not viable. Is this correct?

I have searched the literature and cannot find any statistical papers that refer to this specific question, but I have found some research papers that have used this method. Can anyone direct me to relevant publications that can help me with this?

If indeed I cannot mix the two methods, does anyone have a suggestion for another way to look at differences between factor levels of a fixed effect?

Edit Thank you for the response, but I just want to make sure that the question I am asking is clear. I got this response from a different forum, which may help to clarify the question that I am asking:

"This might depend on how you ran the model comparisons. If you had Season as e.g. a dummy-coded variable or something like that, and put all levels of it in as one step of the model (or e.g. just did it the R way, putting 'Season' in as its own predictor), then indeed model comparison doesn't tell you which levels of Season are different from each other, and therefore you still need some kind of post-hoc tests (or at least just inspecting the model coefficients, but indeed these should likely be corrected for multiple comparisons) to see which levels are different from which."

I am trying to make sure that I am interpreting answers correctly, but for that to happen I have to be sure that my question is clear as well.

Many thanks for all the responses



1 Answer 1


If you choose your model based on the actual observations (using , etc.), your inferential statistics will be invalid for all the reasons gung has described here. Don't do it.

If indeed I cannot mix the two methods, does anyone have a suggestion for another way to look at differences between factor levels of a fixed effect?

Specify your model in advance, before observing any data.

Of course, you can conduct a pilot study and try different models on pilot data. This is completely valid. However, you cannot do valid inferential statistics on this dataset. Instead, specify your model based on your pilot data, and then collect independent new data on which you can perform inferential statistics with the prespecified model.

(Theoretically, you could possibly also account for your model selection in running inferential statistics on the same dataset you also selected your model for. However, I don't know how I would go about this, since this would in turn depend heavily on the specific dataset in question.)

  • $\begingroup$ Reasonable ways to account for the model selection that could be combined with AIC based criteria include to use model averaging or bootstrapping the model selection process (and averaging across the outcomes across the different bootstrap samples) and a number of other approaches. Of course there's other approaches that are not AIC based e.g. LASSO etc. $\endgroup$
    – Björn
    Jun 14, 2016 at 11:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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