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@(protected)
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.