Re: [R] Questing on fitting Baseline category Logit model

Wed, 14 Mar 2012 08:28:02 -0700

Not sure why VGAM::vglm doesn't work here, but most likely it is the small & zero counts cited on the page you quote below. This data set 
is very sparse.  You should communicate with the author of VGAM.

You can fit this model with nnet::multinom instead, something like

library(nnet)
# multinomial logit model
(mod1 <- multinom(food ~ lake+size+sex, data=Alligator, weights=count))

> multinom(food ~ lake+size+sex, data=Alligator, weights=count)
# weights:  35 (24 variable)
initial  value 352.466903
iter  10 value 270.397070
iter  20 value 268.958046
final  value 268.932740
converged
Call:
multinom(formula = food ~ lake + size + sex, data = Alligator,
    weights = count)

Coefficients:
(Intercept) lakeHancock lakeOklawaha lakeTrafford sizesmall sexmale fish 1.70178892 -0.5752403 0.5503569 -1.23679067 0.7303298 0.60639521 invert 0.53452560 -2.3557451 1.4635491 -0.08096449 2.0665999 0.14342792 other -0.01957203 0.1913919 0.5764707 0.32102428 1.0209285 0.35382356 reptile -1.15700455 0.5539263 3.0803954 1.82403973 0.1733300 -0.02116283 

Residual Deviance: 537.8655
AIC: 585.8655
>


On 3/14/2012 4:15 AM, Christofer Bogaso wrote:

Dear all,

I am facing some problem with how to fit a "Baseline category Logit
model" with R. Basically I am considering famous  "Alligator" data as
discussed by Agresti. This data can also be found here:

https://onlinecourses.science.psu.edu/stat504/node/174
(there is also an accompanying R file, however the underlying R code
could not load the data properly!!!)

Below are the stuffs what I have done so far:

My_Data<- structure(list(Number = c(7L, 4L, 1L, 0L, 0L, 0L, 0L, 1L, 5L,
2L, 16L, 3L, 3L, 0L, 2L, 1L, 2L, 2L, 3L, 3L, 2L, 13L, 2L, 7L,
0L, 6L, 0L, 0L, 1L, 0L, 3L, 0L, 9L, 1L, 1L, 0L, 0L, 1L, 2L, 0L,
3L, 8L, 7L, 6L, 1L, 6L, 0L, 3L, 1L, 5L, 2L, 0L, 4L, 1L, 1L, 0L,
1L, 0L, 4L, 0L, 13L, 9L, 10L, 0L, 0L, 0L, 2L, 1L, 2L, 2L, 3L,
8L, 9L, 1L, 1L, 0L, 0L, 0L, 1L, 1L), Food = structure(c(2L, 3L,
5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L,
1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L,
4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L,
2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L,
3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L, 2L, 3L, 5L, 1L, 4L), .Label = c("Bird",
"Fish", "Invertebrate", "Other", "Reptile"), class = "factor"),
     Size = structure(c(2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,
     2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L,
     1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,
     2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L,
     1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L,
     2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L), .Label = c("Large",
     "Small"), class = "factor"), Sex = structure(c(2L, 2L, 2L,
     2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
     1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L,
     1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
     2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L,
     2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
     1L, 1L), .Label = c("Female", "Male"), class = "factor"),
     Lake = structure(c(2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L,
     2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L,
     3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L,
     4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L, 4L,
     4L, 4L, 4L, 4L, 4L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L,
     1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L), .Label = c("George",
     "Hancock", "Oklawaha", "Trafford"), class = "factor")), .Names =
c("Number",
"Food", "Size", "Sex", "Lake"), row.names = c(NA, 80L), class = "data.frame")


library(VGAM)
vglm(Food~Size+Sex+Lake, data = My_Data, fam=multinomial, weights = Number)



However I am getting following error:

Error in if (max(abs(ycounts - round(ycounts)))>  smallno)
warning("converting 'ycounts' to integer in @loglikelihood") :
   missing value where TRUE/FALSE needed
In addition: Warning messages:
1: In checkwz(wz, M = M, trace = trace, wzepsilon = control$wzepsilon) :
   96 elements replaced by 1.819e-12
2: In checkwz(wz, M = M, trace = trace, wzepsilon = control$wzepsilon) :
   96 elements replaced by 1.819e-12
3: In checkwz(wz, M = M, trace = trace, wzepsilon = control$wzepsilon) :
   96 elements replaced by 1.819e-12
4: In checkwz(wz, M = M, trace = trace, wzepsilon = control$wzepsilon) :
   96 elements replaced by 1.819e-12
5: In checkwz(wz, M = M, trace = trace, wzepsilon = control$wzepsilon) :
   96 elements replaced by 1.819e-12
6: In checkwz(wz, M = M, trace = trace, wzepsilon = control$wzepsilon) :
   96 elements replaced by 1.819e-12

Can somebody points me why I am getting this error?

Thanks for you help




--
Michael Friendly     Email: friendly AT yorku DOT ca
Professor, Psychology Dept.
York University      Voice: 416 736-5115 x66249 Fax: 416 736-5814
4700 Keele Street    Web:   http://www.datavis.ca
Toronto, ONT  M3J 1P3 CANADA

______________________________________________
[email protected] 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.

Reply via email to