#
#
########## Multinomial Logistic Regression example ##########
#
#
# This script assumes you have worked through all the previous notes from
# the web page and you have downloaded, installed, and updated all available
# R packages.
# Load the following libraries if you have not already.
library(foreign)
library(Rcmdr)
# The MultiNomReg.sav data file contains one categorical outcome (y) with 3 categories,
# and three continuous predictors (x1, x2, x3).
# Load the "MultiNomReg.sav" data file naming is 'mdata1'.
mdata1 <- read.spss("http://researchsupport.unt.edu/class/Jon/R_SC/Module9/MultiNomReg.sav",
use.value.labels=TRUE, max.value.labels=Inf, to.data.frame=TRUE)
summary(mdata1)
# Next, change the outcome variable (which is currently numeric) to a factor.
mdata2 <- mdata1
mdata2$y <- as.factor(mdata2$y)
summary(mdata2)
## Conduct the multinomial analysis using the mlogit library; which allows you to specify the reference
## category of the outcome.
# Load the appropriate library.
library(mlogit)
# Then create a new 'form' of data which can be read by the mlogit function.
# To do this, we need to expand the outcome variable (y) much like we would
# for dummy coding a categorical variable for inclusion in standard multiple
# regression.
mdata3 <- mlogit.data(mdata2, varying=NULL, choice="y", shape="wide")
head(mdata3)
# Run the multinomial analysis with 'mlogit' function, specifying the reference category as "1".
model.1 <- mlogit(y ~ 1 | x1 + x2 + x3, data=mdata3, reflevel="1")
summary(model.1)
# To get the "exp(b)" values (ignore the intercepts which are the first 3 printed and are
# labeled "alt[T.2]", "alt[T.3]", "alt[T.4]"):
exp(coef(model.1))
# If desired, residuals can be produced.
residuals(model.1)
hist(residuals(model.1))
# Rcmdr can lead to different results because it uses the "nnet" package.
# If you specify the same model using Rcmdr, you get the following script; which
# produces slightly different results when compared to what is above.
MLM.1 <- multinom(y ~ x1 +x2 +x3, data=mdata2, trace=FALSE)
summary(MLM.1, cor=FALSE, Wald=TRUE)