thr3ads.net - R help - [R] Questing on fitting Baseline category Logit model [Mar 2012]

If this information is useful, please help other people find it:
Share via:

Christofer Bogaso

2012-Mar-14 08:15 UTC

[R] Questing on fitting Baseline category Logit model

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

2012-Mar-14 15:19 UTC

head link

[R] Questing on fitting Baseline category Logit model

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

Reasonably Related Threads

Search for more possibly parallel threads

R help - Mar 2012 - Questing on fitting Baseline category Logit model

[R] Questing on fitting Baseline category Logit model

[R] Questing on fitting Baseline category Logit model

Reasonably Related Threads