[R-sig-ME] GLMM for underdispersed count data: Conway-Maxwell-Poisson and Ordinal

[Mollie Brooks]

I wanted to follow up on this discussion and let you know that glmmTMB recently changed the way that it does Conway-Maxwell-Poisson regression. So if you are using glmmTMB for this, you should probably install the new version from source. It now is parameterized to do regression on the mean rather than the approximate mean, which should improve accuracy and interpretability. It’s part of the master branch now, so the standard installation from source should work as described here https://github.com/glmmTMB/glmmTMB

cheers,
Mollie 

———————————
Mollie E. Brooks, Ph.D.
Postdoctoral Researcher
National Institute of Aquatic Resources
Technical University of Denmark

> On 15Dec 2016, at 1:52, Ben Bolker <bbolker at gmail.com> wrote:
> 
> 
> 
> 
> On 16-12-14 11:05 AM, Simone Santoro wrote:
>> Thank you all so much for it. Really a very useful discussion, hope it
>> may be so to others too.
>> 
>> El 14/12/2016 a las 16:06, Mollie Brooks escribió:
>>> Hi Simone,
>>> 
>>> For the glmmTMB model with the Conway-Maxwell Poisson distribution,
>>> the left side of the equation should technically by fledges rather
>>> than as.factor(fledges). However, it looks like glmmTMB doesn’t
>>> evaluate the as.factor() command and fits the model with fledges as
>>> the response anyway.
> 
>  I'd be careful with this conclusion.  I think these are *not* the same
> model.  For the fake data set (where all the true effects are zero) the
> results aren't that different, but what happens when you convert an
> integer value to a factor is that the unique values get converted to
> codes 1, 2, ...  This would potentially be disastrous.
> 
>  I got the clmm model to run with Rune's development version; it's a
> little hard to see whether the results are comparable or not since it's
> fitting a qualitatively different model ...
> 
> 
>>> 
>>> If you end up needing zero-inflation also, it can be specified using
>>> the ziformula command. See vignette("glmmTMB") or here
>>> https://github.com/glmmTMB/glmmTMB/blob/master/misc/salamanders.pdf
>>> for an example.
>>> 
>>> cheers,
>>> Mollie
>>> 
>>> ———————————
>>> Mollie E. Brooks, Ph.D.
>>> Postdoctoral Researcher
>>> National Institute of Aquatic Resources
>>> Technical University of Denmark
>>> 
>>>> On 9Dec 2016, at 19:40, Simone Santoro <santoro at ebd.csic.es
>>>> <mailto:santoro at ebd.csic.es>> wrote:
>>>> 
>>>> Hi,
>>>> 
>>>> Thank you all very much your hints. They have been really really
>>>> helpful for me. Below you may find a reproducible code to see how
>>>> three approaches fit a simulated data set (clmm::ordinal,
>>>> glmmTMB::glmmTMB, fitme:spaMM). Results seem to me qualitatively
>>>> similar but with clmm:ordinal I cannot use the three crossed random
>>>> effects because I get an error like this:
>>>> Error: no. random effects (=135) >= no. observations (=100)
>>>> 
>>>> set.seed(1234)
>>>> library(ordinal)
>>>> library(glmmTMB)
>>>> library(spaMM)
>>>> dati<- data.frame(fledges= rpois(100,10), habitatF=
>>>> as.factor(rbinom(100,1,0.5)), areaPatchFath= rnorm(100), poligF01=
>>>> as.factor(rbinom(100,1,0.5)),StdLayingDate= rnorm(100), ageFath1=
>>>> rpois(100,3), ageMoth1= rpois(100,3), year=
>>>> as.factor(rpois(100,200)), ringMoth= as.factor(rpois(100,200)),
>>>> ringFath= as.factor(rpois(100,200)))
>>>> str(dati)
>>>> 
>>>> system.time(Fitclm<- clmm(as.factor(fledges) ~
>>>> habitatF*(areaPatchFath+poligF01+StdLayingDate+ageFath1+ageMoth1)+(1|year)+(1|ringMoth)+(1|ringFath),data=dati,Hess=T))
>>>> # this way it works...
>>>> system.time(Fitclm1<- clmm(as.factor(fledges) ~
>>>> habitatF*(areaPatchFath+poligF01+StdLayingDate+ageFath1+ageMoth1)+(1|year)+(1|ringFath),data=dati,Hess=T))
>>>> summary(Fitclm1)
>>>> 
>>>> system.time(FitglmmTMB<- glmmTMB(as.factor(fledges) ~
>>>> habitatF*(areaPatchFath+poligF01+StdLayingDate+ageFath1+ageMoth1)+(1|year)+(1|ringMoth)+(1|ringFath),data=dati,family=
>>>> "compois"))
>>>> summary(FitglmmTMB)
>>>> 
>>>> system.time(FitglmmTMB<- glmmTMB(as.factor(fledges) ~
>>>> habitatF*(areaPatchFath+poligF01+StdLayingDate+ageFath1+ageMoth1)+(1|year)+(1|ringMoth)+(1|ringFath),data=dati,family=
>>>> "compois"))
>>>> summary(FitglmmTMB)
>>>> 
>>>> # This lasts much more (3-4')
>>>> system.time(Fitfitme<- fitme(fledges ~
>>>> habitatF*(areaPatchFath+poligF01+StdLayingDate+ageFath1+ageMoth1)+(1|year)+(1|ringFath)+(1|ringMoth),data=dati,COMPoisson(),method
>>>> = "ML"))
>>>> summary(Fitfitme)
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> 
>>>> El 08/12/2016 a las 4:32, Ben Bolker escribió:
>>>>>   One reference that uses ordinal regression in a similar situation
>>>>> (litter size of Florida panthers) is
>>>>> http://link.springer.com/article/10.1007/s00442-011-2083-0 ("Does
>>>>> genetic introgression improve female reproductive performance? A test on
>>>>> the endangered Florida panther")
>>>>> 
>>>>>  Not sure about the number-of-random-effects error: a reproducible
>>>>> example would probably be needed (smaller is better!)
>>>>> 
>>>>>  Ben Bolker
>>>>> 
>>>>> 
>>>>> On 16-12-06 08:41 AM, Simone Santoro wrote:
>>>>>> Dear all,
>>>>>> 
>>>>>> I am trying to find an appropriate GLMM (with temporal and individual 
>>>>>> crossed random effects) to model underdispersed count data (clutch 
>>>>>> size). I have found several possible ways of doing that. A good 
>>>>>> distribution for data like this would seem to be the 
>>>>>> Conway-Maxwell-Poisson but I have not found a way of using it within a 
>>>>>> GLMM in R (I have asked here 
>>>>>> <http://stats.stackexchange.com/questions/249738/how-to-define-the-nu-parameter-of-conway-maxwell-poisson-in-spamm-package> 
>>>>>> and here 
>>>>>> <http://stats.stackexchange.com/questions/249798/conway-maxwell-poisson-with-crossed-random-effects-in-r>).
>>>>>> I have seen that Ben Bolker suggested (here 
>>>>>> <https://stat.ethz.ch/pipermail/r-sig-mixed-models/2014q1/021945.html>and 
>>>>>> here 
>>>>>> <http://stats.stackexchange.com/questions/92156/how-to-handle-underdispersion-in-glmm-binomial-outcome-variable>) 
>>>>>> to use an ordinal model in cases like this(e.g. _ordinal:clmm_). I have 
>>>>>> tried this solution and the results I obtain makes (biological) sense to 
>>>>>> me. However, I wonder why but I cannot put all the three crossed random 
>>>>>> effects I have in the clmm model (_Error: no. random effects (=1254) >= 
>>>>>> no. observations (=854)_) whereas it is not a problem for the glmer 
>>>>>> model (the no. of levels of each single random effect does not exceed 854)*.
>>>>>> Beyond that, and that's what I would like to ask you, *I cannot find a 
>>>>>> reference to justify I used the ordinal model* to deal with 
>>>>>> underdispersed count data (referee will ask it for sure).
>>>>>> Best,
>>>>>> 
>>>>>> Simone
>>>>>> 
>>>>>> * FMglmer<- glmer(fledges ~ habitatF * (areaPatchFath + poligF01 + 
>>>>>> StdLayingDate + ageFath1 + ageMoth1) + (1|year) + (1|ringMoth) + 
>>>>>> (1|ringFath), data = datiDRS)
>>>>>>    FMclmm<- glmer(as.factor(fledges)~ habitatF * (areaPatchFath + 
>>>>>> poligF01 + StdLayingDate + ageFath1 + ageMoth1) + (1|year) + 
>>>>>> (1|ringMoth) + (1|ringFath), data = datiDRS)
>>>>>> 
>>>>>> 
>>>>>> 	[[alternative HTML version deleted]]
>>>>>> 
>>>>>> _______________________________________________
>>>>>> R-sig-mixed-models at r-project.org mailing list
>>>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>>>>> 
>>>> 
>>>> -- 
>>>> Simone Santoro
>>>> PhD
>>>> Department of Ethology and Biodiversity Conservation
>>>> Doñana Biological Station
>>>> Calle Américo Vespucio s/n
>>>> 41092 Seville - Spain
>>>> Phone no. +34 954 466 700 (ext. 1213)
>>>> http://www.researchgate.net/profile/Simone_Santoro
>>>> http://orcid.org/0000-0003-0986-3278
>>> 
>> 
>> -- 
>> Simone Santoro
>> PhD
>> Department of Ethology and Biodiversity Conservation
>> Doñana Biological Station
>> Calle Américo Vespucio s/n
>> 41092 Seville - Spain
>> Phone no. +34 954 466 700 (ext. 1213)
>> http://www.researchgate.net/profile/Simone_Santoro
>> http://orcid.org/0000-0003-0986-3278
>> 


	[[alternative HTML version deleted]]