# Using delta method to get confidence intervals for multinomial logit?

I am working with some choice modeling data and am interested in trying to potentially use the delta method with the multinomial logit model that I'm analyzing the data with. Here's an example:

First, I uploaded some data in R:

# First upload and format the data
colClasses = c(seat = "factor", price = "factor",
choice="integer"))
cbc.df$$eng <- factor(cbc.df$$eng, levels=c("gas", "hyb", "elec"))
cbc.df$$carpool <- factor(cbc.df$$carpool, levels=c("yes", "no"))
summary(cbc.df)


Format data for mlogit package depending on version:

if (packageVersion("mlogit") < "1.1") {
# for mlogit up through 1.0.3
cbc.mlogit <- mlogit.data(data=cbc.df, choice="choice", shape="long",
varying=3:6, alt.levels=paste("pos", 1:3),
id.var="resp.id")
} else {
# for mlogit starting with version 1.1
library(dfidx)      # install if needed

# add a column with unique question numbers, as needed in mlogit 1.1+
cbc.df\$chid <- rep(1:(nrow(cbc.df)/3), each=3)

# shape the data for mlogit
cbc.mlogit <- dfidx(cbc.df, choice="choice",
idx=list(c("chid", "resp.id"), "alt" ))
}


I think fit a basic multinomial logit model:

# fit the models
library(mlogit)
m1 <- mlogit(choice ~ 0 + seat + cargo + eng + price, data = cbc.mlogit)


Next, I want to predict the probability of making a choice given different combinations of the predictors. I found this function to do that:

# Predicting shares
predict.mnl <- function(model, data) {
# Function for predicting shares from a multinomial logit model
# model: mlogit object returned by mlogit()
# data: a data frame containing the set of designs for which you want to
#       predict shares.  Same format at the data used to estimate model.
data.model <- model.matrix(update(model$$formula, 0 ~ .), data = data)[ , -1] utility <- data.model%*%model$$coef
share <- exp(utility)/sum(exp(utility))
cbind(share, data)
}


So then I generated some combinations of predictors to generate predictions for:

> # and set up attributes needed later
> attrib <- list(seat  = c("6", "7", "8"),
+                cargo = c("2ft", "3ft"),
+                eng   = c("gas", "hyb", "elec"),
+                price = c("30", "35", "40"))
>
>
>
> (new.data <- expand.grid(attrib)[c(8, 1, 3, 41, 49, 26), ]) # find attrib at top
seat cargo  eng price
8     7   2ft  hyb    30
1     6   2ft  gas    30
3     8   2ft  gas    30
41    7   3ft  gas    40
49    6   2ft elec    40
26    7   2ft  hyb    35


And finally, get the predictions for these combinations:

> predict.mnl(m1, new.data)
share seat cargo  eng price
8  0.11273356    7   2ft  hyb    30
1  0.43336911    6   2ft  gas    30
3  0.31917819    8   2ft  gas    30
41 0.07281396    7   3ft  gas    40
49 0.01669280    6   2ft elec    40
26 0.04521237    7   2ft  hyb    35


This share variable is great, but what I'm really wanting is something like a 95% CI for it and I've heard that the delta method could potentially help. I just think that knowing that the first combination has an 11% share is meaningless without knowing the variability in the estimate.

Can anyone help?

• Can I ask did you solve this issue using mlogit and the delta method? Aug 16, 2021 at 12:30

Sort of not answering your question, but the brms R package covers the multinomial model. E.g. like this brm(bf(y | trials(size) ~ 1, muy2 ~ x), data = mydata, family = multinomial(), prior = myprior) (see e.g. here). The great thing - assuming you are happy with a Bayesian approach - is that you can then get predictions based on each MCMC sample.