Dear all,

Thank you for the clarifications on the correct use and coding of nested random effects. I decided to try the MCMCglmm package for my analysis, but encountered some other problems.

First: dataset of 57 points (each point has a unique ID (a,b,c,...) and is nested in a certain site (26 sites, Site1,Site2,...). Points have been visited 5 times (v1,v2,...,v5). For simplicity, I want to relate the number of birds counted to one environmental variable, called pca1. Using most default settings, the model should look like:

MCMCglmm(abundance~data$pca1,family="zipoisson",random=~site+point+visit,data=data,verbose=TRUE)

However, this results in the following error:"Error in MCMCglmm(abundance ~ data$pca1, family = "zipoisson", random = ~site +  :   R-structure does not define unique residual for each data point". [[ Note that if I specify another distribution, e.g. an inappropriate Poisson, this error does not occur]]. As far as I understand, this error relates to the rcov part of the model and I updated the syntax to:

MCMCglmm(abundance~data$pca1,family="zipoisson",rcov=~us(trait):units,random=~site+point+visit,data=data,verbose=TRUE)

However, this results in another error:
MCMC iteration = 0
E[MH acceptance ratio] = 0.000431
MCMC iteration = 1000
E[MH acceptance ratio] = 0.439460
MCMC iteration = 2000
E[MH acceptance ratio] = 0.447993
Error in MCMCglmm(abundance ~ data$pca1, family = "poisson", random = ~site +  : 
  ill-conditioned G/R structure: use proper priors if you haven't or rescale data if you have

>From this I understand that I should define appropriate priors. Based on the Plodia tutorial, I tried 3 priors:
halfV = var(abundance)/2
prior=list(R=list(n=1,V=halfV),G=list(G1=list(n=1,V=halfV))) 
##halfV = 0.08639657##
prior = list(R=list(V=1,n=1e+06),G=list(G1=list(V=1,n=0)))
prior = list(R=list(V=1,n=0,fix=1),G=list(G1=list(V=1, +   n=0)))

However, all of them result in the following error:
Error in MCMCglmm(abundance ~ data$pca1, family = "zipoisson", rcov = ~us(trait):units,  : 
  ill-conditioned G/R structure: use proper priors if you haven't or rescale data if you have

Could anyone offer some advice on what my mistakes are? Or suggested papers to read, as this is my first encounter with Bayesian based analysis.

Many thanks,

Diederik


Diederik Strubbe
Evolutionary Ecology Group
Department of Biology, University of Antwerp
Universiteitsplein 1
B-2610 Antwerp, Belgium
http://webhost.ua.ac.be/deco
tel : 32 3 820 23 85



-----Original Message-----
From: Jarrod Hadfield [mailto:jhadfiel@staffmail.ed.ac.uk]
Sent: Thu 23-4-2009 9:33
To: Strubbe Diederik
Cc: r-sig-mixed-models@r-project.org; boristmr@yahoo.com
Subject: Re: [R-sig-ME] About glmer with a quasi-Poisson
 


Diederik Strubbe
Evolutionary Ecology Group
Department of Biology, University of Antwerp
Universiteitsplein 1
B-2610 Antwerp, Belgium
http://webhost.ua.ac.be/deco
tel : 32 3 820 23 85



-----Original Message-----
From: Jarrod Hadfield [mailto:jhadfiel@staffmail.ed.ac.uk]
Sent: Thu 23-4-2009 9:33
To: Strubbe Diederik
Cc: r-sig-mixed-models@r-project.org; boristmr@yahoo.com
Subject: Re: [R-sig-ME] About glmer with a quasi-Poisson
 
Hi,

This model can be fitted using MCMCglmm (and others I expect) if it is  
coded correctly. As Doug mentioned in an earlier email, calling  
different points in different sites the same name is not really  
necessary any more. If they are called different names then:

random=~site+point

is the same as the glmer model you use. If they are called the same thing then

random=~site+site:point

is the same as the glmer model below, as can be seen in the glmer summary.

Coded this way many packages can probably fit the type of model you envisage.

Cheers,

Jarrod





Quoting Strubbe Diederik <diederik.strubbe@ua.ac.be>:

> Hey Mike,
>
> I am working on a similar dataset as you. I have point counts of   
> birds in 26 sites, in total 57 points [number of point counts per   
> site varying from 1 to 5]. These points were visited 5 times ( in   
> the same year). Most of the counts are 0 (> 90%), and the dataset is  
>  thus zero-inflated - which, as far as I know is a special case of   
> overdispersion. There  are 3 continuous explanatory variables under   
> consideration. I have been looking at different packages to fit a   
> model with
> 1)a  zero-inflated poison distribution
> 2)a nested random effect to account for the fact that I have points   
> located within sites
> 3)a random effect to account for the 5 visits (repeated measures)
>
> However, none of the packages that I tried is able to do this ( I   
> tried zeroinfl {pscl}, fmr {gnlm}, MCMCglmm {MCMCglmm} and   
> glmm.ADMB{glmmADMB}). The problem seems to be the nested random   
> effect (1|site/point in my case).
>
> The only function that I found able to handle this nested random   
> design and a distribution that comes close to be 'good' for my data   
> is glmer with a quasi-poisson distribution. My  model syntax is:
>
> test <-   
> glmer(data$abundance~data$pca1+data$pca2+data$Gynoxis+(1|site/point)  
>  + (1|visit),family=quasipoisson(link="log"))
>
> and yields the folowing output
>
> Generalized linear mixed model fit by the Laplace approximation
> Formula: data$abundance ~ data$pca1 + data$pca2 + data$Gynoxis + (1   
> |      site/point) + (1 | visit)
>    AIC   BIC logLik deviance
>  158.1 187.1 -71.05    142.1
> Random effects:
>  Groups     Name        Variance  Std.Dev.
>  point:site (Intercept) 0.0468101 0.216356
>  site       (Intercept) 0.0000000 0.000000
>  visit      (Intercept) 0.0012396 0.035208
>  Residual               0.0187386 0.136889
> Number of obs: 276, groups: point:site, 57; site, 26; visit, 5
>
> Fixed effects:
>              Estimate Std. Error t value
> (Intercept)  -3.20340    0.10701 -29.936
> data$pca1    -0.02629    0.03100  -0.848
> data$pca2    -0.15055    0.02907  -5.180
> data$Gynoxis -0.02246    0.01320  -1.702
>
> Correlation of Fixed Effects:
>             (Intr) dt$pc1 dt$pc2
> data$pca1   -0.239
> data$pca2   -0.098 -0.071
> data$Gynoxs -0.878  0.330  0.205
>
> I was not aware of the comments of Ben Bolker on glmer and   
> quasi-Poisson, and am thus also interested in any responses on this.  
>  Am I right to conclude that, ecept glmer, none of the packages   
> mentioned above is capable to handle nested random effects and   
> zero-inflated Poisson data?
>
> Cheers and many thanks,
>
> Diederik
>
>
>
> Diederik Strubbe
> Evolutionary Ecology Group
> Department of Biology, University of Antwerp
> Universiteitsplein 1
> B-2610 Antwerp, Belgium
> http://webhost.ua.ac.be/deco
> tel : 32 3 820 23 85
>
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>



-- 
The University of Edinburgh is a charitable body, registered in
Scotland, with registration number SC005336.





	[[alternative HTML version deleted]]


