2
$\begingroup$

I am trying to determine which evolutionary model is best for my discrete data using the function fitDiscrete() from the geiger package.

These are the values that I get for the number of parameters (k), maximum log likelihood (lnL), AIC and AICc for each model:

     k   lnL        AIC       AICc
ER   1   -115.8006  233.6012  233.6637
ARD  90  -85.98459  351.9692  -303.2308
SYM  45  -97.23202  284.4640  491.464

The same dataset (n = 66), tree and single trait with 10 levels were used to create each model. The only difference is the evolutionary model fitted (equal rates (ER), all rates different (ARD) and symmetrical rates (SYM))

I am having trouble interpreting these results, however.

To start, for AIC, I'm fairly sure that I should select the model with the smallest AIC score, i.e the ER model.

For lnL, however, I have seen that the model with the "largest value" should be selected with this being interpreted as the value closest to 0 (https://www.r-phylo.org/wiki/HowTo/Ancestral_State_Reconstruction), i.e. the ARD model. I realise though that lnL values tend to be biased towards models with higher k values. To address this, I did do a likelihood test as suggested by the website above (chi-squared test), which came to p < 0.001. This would suggest that the ARD model should be preferred over the ER model, which contradicts what the AIC scores are telling me.

As for AICc, again, the "smallest value" should be selected but the negative sign mixed in with the positive ones has thrown me. Is this the smallest absolute value or the value value closest to negative infinity?

So, all in all, how can I tell which model should be preferred?

Improve this question
$\endgroup$
6
  • $\begingroup$ How many degrees of freedom did you use in your chi-square test? With 89 degrees of freedom, the difference between the ER and ARD log likelihoods doesn't look significant to me... And is the ER model nested in the ARD model? $\endgroup$
    – jbowman
    Commented Nov 12, 2020 at 19:18
  • $\begingroup$ I've had to double check the website above but the df should be the difference between the values for k so 89 here. The likelihood test is calculated through the following command (from the webiste above): 1-pchisq(2*abs(ER(lnL) - ARD(lnL)), 1). ER should be nested in ARD but the code isn't clear: ER <- fitDiscrete(phy = lattree3, dat = f, model = "ER") and ARD <- fitDiscrete(phy = lattree3, dat = f, model = "ARD"). They contain the same variables but ARD is more complex so I'm going to guess ER is nested within. Sorry I couldn't be more certain! I'm very new to this! $\endgroup$
    – Carolina Karoullas
    Commented Nov 12, 2020 at 20:09
  • $\begingroup$ The degrees of freedom of the chi-square test is the difference in the dimensionality of the parameters, in this case, 89: en.wikipedia.org/wiki/Wilks%27_theorem , the example should make it clear. This should help reconcile the difference you are seeing! ... so pchisq(2*abs(ER(lnL) - ARD(lnL)), 89). $\endgroup$
    – jbowman
    Commented Nov 12, 2020 at 20:26
  • $\begingroup$ Ahh that would be it! I've just tried it and it all works now! I knew the df should be 89 but the example didn't specify where to insert that and I didn't realise that their model only had 1 df! Thank you :) $\endgroup$
    – Carolina Karoullas
    Commented Nov 12, 2020 at 22:23
  • $\begingroup$ The AICc value for ARD cannot be correct since AICc is always larger than the corresponding AIC. The AICc adds an additional penalty for parameter number. See here: en.wikipedia.org/wiki/… $\endgroup$
    – abstrusiosity
    Commented Nov 12, 2020 at 23:07

0

Reset to default

Your Answer

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