[R-sig-ME] MCMCglmm with zibinomial

[Jarrod Hadfield]

Hi,


Below is a zero-inflated binomial example. I can only get your warning  
if the data are not positive integers. Are you sure this is not the  
case?

You can fit multiple correlation structures (e.g. animal and sire)  
using the ginverse argument.

Cheers,

Jarrod


library(VGAM)

n<-1000                                # number of observations
size<-rpois(n,1)+10            # number of trials per observation
prob<-rnorm(n, 1, sqrt(2))  # logit probability of success (with variance 2)
phi.logit<-(-1)  # logit probability of zero (ignoring the binomial  
distribution)

y<-rzibinom(n, size, plogis(prob), phi = plogis(phi.logit))

detach(package:VGAM) # detach: VGAM conflicts with everything

dat<-data.frame(success=y, failure=size-y)

prior<-list(R=list(V=diag(2), nu=0, fix=2))

# as in binary data the observation-level heterogeneity in  
zero-inflation cannot be estimated - the varaince of these effects  
I've arbitrarily set at 1.

m1<-MCMCglmm(c(success, failure)~trait-1, rcov=~idh(trait):units,  
data=dat, family="zibinomial", prior=prior)

# below posteriors with simulation value in red

hist(m1$Sol[,1])
abline(v=1, col="red")

c2 <- (16 * sqrt(3)/(15 * pi))^2  # correction fator because  
VGAM::rzbinom assumes the variance in logit zero-inflation  
probabilities is zero not one as in MCMCglmm

hist(m1$Sol[,2]/sqrt(1+c2))
abline(v=phi.logit, col="red")

hist(m1$VCV[,1])
abline(v=2, col="red")



Quoting "mjgeha at huskers.unl.edu" <mjgeha at huskers.unl.edu> on Tue, 7  
Feb 2012 18:33:49 +0000:

> I have been trying to fit a zibinomial distribution to a set of data I have.
>
> The response variable is in a single column and represents the  
> percentages of a certain incidence (ranging obviously between 0 and  
> 1).
>
> The data is assumed to be zero-inflated as > 40% of the data is made  
> up of zeroes.
>
> I have tried to fit the following:
>
> analysis.1<-MCMCglmm(response~1,random=~animal,rcov=~units,family="zibinomial",nitt=50000,burnin=5000,thin=25,data=source1,pedigree=ped.dat)
>
>
>
> I get the same error message as reported in a previous question  
> regarding the "zibinomial" :
>
> "Error in rowSums(data[, match(response.names[0:1 + nt], names(data))]) :
>   error in evaluating the argument 'x' in selecting a method for  
> function 'rowSums': Error in `[.data.frame`(data, ,  
> match(response.names[0:1 + nt], names(data))) :
>   undefined columns selected"
>
>
>
> I looked up the answer submitted to that question and tried creating  
>  "success" and "failure" variables and fit the following:
>
> MCMCglmm(cbind(success,failure)~  
> -1,random=~us(trait):animal,rcov=~us(trait):units,family=rep("zibinomial",2),nitt=50000,burnin=1000,thin=25,data=snap)
>
> Now I'm getting the following error message:
>
> "Error in MCMCglmm(cbind(success, failure) ~ - 1, random =  
> ~us(trait):animal,  :
>   binomial data must be non-negative integers"
>
> I'm sure I don't have any negative integers in my data set.
>
> I was wondering if it would be possible to get more  
> clarification/example/reference on fitting zibinomial in MCMCglmm.
>
>
>
> Also would it be possible to fit more than one pedigree relationship  
> matrix? The current example is for an animal model but what if we  
> had both an animal effect and a sire effect and we want to fit two  
> pedigree relationship matrices instead of just one?
>
>
>
> Any help on that issue is more than appreciated.
>
>
>
> Thanks,
>
>
>
> Mak
>
>
>
> Makram J. Geha,PhD
> Quantitative Geneticist
> Dow AgroSciences, LLC.
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-mixed-models at 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.