3

I'm trying to build an AIC table for my candidate set of models in R, which were run using mlogit. I've used glm and glmer in the past, and have always used the package AICcmodavg and aictab to extract values and create a model selection table. This package doesn't seem to work for mlogit, so I'm wondering if there are any other ways of creating an AIC table in R besides manual calculation using the log-likelihood value?

Example of mlogit model output:

Call:
mlogit(formula = Case ~ Dist_boulder + Mesohabitat + Depth + 
    Size + Size^2 | -1, data = reach.dc, method = "nr")

Frequencies of alternatives:
0 1 2 3 
1 0 0 0 

nr method
5 iterations, 0h:0m:0s 
g'(-H)^-1g = 1.19E-05 
successive function values within tolerance limits 

Coefficients :
                   Estimate Std. Error z-value Pr(>|z|)  
Dist_boulder      -0.052165   0.162047 -0.3219  0.74752  
Mesohabitatriffle -1.400752   0.612329 -2.2876  0.02216 *
Mesohabitatrun     0.302697   0.420181  0.7204  0.47128  
Depth              0.137524   0.162521  0.8462  0.39745  
Size               0.336949   0.145036  2.3232  0.02017 *
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Log-Likelihood: -86.627

example of models run (from my candidate set of 14)

 predation.reach<-mlogit(Case ~ Dist_boulder + Mesohabitat + Depth + Size + Size^2 | -1, data=reach.dc)
velocity.reach<-mlogit(Case ~ Mid_vel | -1, data=reach.dc)
spaces.reach<-mlogit(Case~ Embedded + Class | -1, data=reach.dc)
substrate.reach<-mlogit(Case ~ Class | -1, data=reach.dc)

defining candidate set list

cand.set.reach<-list(predation.reach, velocity.reach, spaces.reach, substrate.reach)
Community
  • 1
  • 1
  • `broom::augment()` will be the way to go, but not yet ... https://github.com/tidymodels/broom/issues/192 Does `bbmle::AICtab()` work ... ? – Ben Bolker Sep 15 '19 at 17:59

1 Answers1

2

bbmle::AICtab() appears to work.

library("mlogit")
m1 <- mlogit(formula = mode ~ price + catch | income,
  data = Fish,     
  alt.subset = c("charter", "pier", "beach"), method = "nr")
m2 <- update(m1, . ~ . - price)
bbmle::AICtab(m1,m2)
##    dAIC  df
## m1   0.0 6 
## m2 412.1 5 

By default bbmle::AICtab() gives only delta-AIC and the model degrees of freedom/number of parameters, but you can use optional arguments to get the absolute AIC, AIC weights, etc..

It also works with a list:

L <- list(m1,m2)
bbmle::AICtab(L)

In the tidyverse world,

library(broom)
L %>% purrr::map(augment) %>% bind_rows()

ought to work, but doesn't yet.

Ben Bolker
  • 211,554
  • 25
  • 370
  • 453