As we noted earlier and as is clearly stated in the docs, multinomRob is estimating an OVERDISPERSED multinomial model. And in your models here the overdispersion parameter is not identified; you need more observations. Walter pointed out using the print.level trick to get the coefs for the standard MNL model, but when the model the function is actually trying to estimate is not identified, that trick will not work.
As I also previously noted, it is a simple matter to change the multinomMLE function to estimate the standard MNL model. Since you don't want to change that file and since nnet's multinom function doesn't have some functionality people need, we'll add a "MLEonly" function to multinomRob which will allow you to do what you want. We'll post a new version on my webpage later tonight: http://sekhon.berkeley.edu/robust. And after some testing, we'll forward the new version to CRAN. Jas. ======================================= Jasjeet S. Sekhon Associate Professor Travers Department of Political Science Survey Research Center UC Berkeley http://sekhon.berkeley.edu/ V: 510-642-9974 F: 617-507-5524 ======================================= Roger Levy writes: > Walter Mebane wrote: > > Roger, > > > > > Error in if (logliklambda < loglik) bvec <- blambda : > > > missing value where TRUE/FALSE needed > > > In addition: Warning message: > > > NaNs produced in: sqrt(sigma2GN) > > > > That message comes from the Newton algorithm (defined in source file > > multinomMLE.R). It would be better if we bullet-proofed it a bit > > more. The first thing is to check the data. I don't have the > > multinomLogis() function, so I can't run your code. > > Whoops, sorry about that -- I'm putting revised code at the end of the > message. > > > But do you really > > mean > > > > > for(i in 1:length(choice)) { > > and > > > dim(counts) <- c(length(choice),length(choice)) > > > > Should that be > > > > for(i in 1:n) { > > and > > dim(counts) <- c(n, length(choice)) > > > > or instead of n, some number m > length(choice). As it is it seems to > > me you have three observations for three categories, which isn't going > > to work (there are five coefficient parameters, plus sigma for the > > dispersion). > > I really did mean for(i in 1:length(choice)) -- once again, the proper > code is at the end of this message. > > Also, I notice that I get the same error with another kind of data, > which works for multinom from nnet: > > > library(nnet) > library(multinomRob) > dtf <- data.frame(y1=c(1,1),y2=c(2,1),y3=c(1,2),x=c(0,1)) > summary(multinom(as.matrix(dtf[,1:3]) ~ x, data=dtf)) > summary(multinomRob(list(y1 ~ 0, y2 ~ x, y3 ~ x), data=dtf,print.level=128)) > > > The call to multinom fits the following coefficients: > > Coefficients: > (Intercept) x > y2 0.6933809622 -0.6936052 > y3 0.0001928603 0.6928327 > > but the call to multinomRob gives me the following error: > > multinomRob(): Grouped MNL Estimation > [1] "multinomMLE: -loglik initial: 9.48247391895106" > Error in if (logliklambda < loglik) bvec <- blambda : > missing value where TRUE/FALSE needed > In addition: Warning message: > NaNs produced in: sqrt(sigma2GN) > > > Does this shed any light on things? > > > Thanks again, > > Roger > > > > > > *** > > set.seed(10) > library(multinomRob) > multinomLogis <- function(vector) { > x <- exp(vector) > z <- sum(x) > x/z > } > > n <- 20 > choice <- c("A","B","C") > intercepts <- c(0.5,0.3,0.2) > prime.strength <- rep(0.4,length(intercepts)) > counts <- c() > for(i in 1:length(choice)) { > u <- intercepts[1:length(choice)] > u[i] <- u[i] + prime.strength[i] > counts <- c(counts,rmultinomial(n = n, pr = multinomLogis(u))) > } > dim(counts) <- c(length(choice),length(choice)) > counts <- t(counts) > row.names(counts) <- choice > colnames(counts) <- choice > data <- data.frame(Prev.Choice=choice,counts) > > for(i in 1:length(choice)) { > data[[paste("last",choice[i],sep=".")]] <- > ifelse(data$Prev.Choice==choice[i],1,0) > } > > multinomRob(list(A ~ last.A , > B ~ last.B , > C ~ last.C - 1 , > ), > data=data, > print.level=128) > > > > I obtained this output: > > > Your Model (xvec): > A B C > (Intercept)/(Intercept)/last.C 1 1 1 > last.A/last.B/NA 1 1 0 > > [1] "multinomRob: WARNING. Limited regressor variation..." > [1] "WARNING. ... A regressor has a distinct value for only one > observation." > [1] "WARNING. ... I'm using a modified estimation algorithm (i.e., > preventing LQD" > [1] "WARNING. ... from modifying starting values for the affected > parameters)." > [1] "WARNING. ... Affected parameters are TRUE in the following table." > > A B C > (Intercept)/(Intercept)/last.C FALSE FALSE TRUE > last.A/last.B/NA TRUE TRUE FALSE > > > > multinomRob(): Grouped MNL Estimation > [1] "multinomMLE: -loglik initial: 70.2764843511374" > Error in if (logliklambda < loglik) bvec <- blambda : > missing value where TRUE/FALSE needed > In addition: Warning message: > NaNs produced in: sqrt(sigma2GN) ______________________________________________ R-help@stat.math.ethz.ch mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
