[R-sig-ME] Bivariate Poisson MCMCglmm

Jarrod Hadfield j.hadfield at ed.ac.uk
Tue Nov 27 15:23:06 CET 2012


Hi Katie,

Your model makes sense. The different malaria species have different  
abundances (trait-1) and their abundances change with age at different  
rates (trait:Age). The 2x2 covariance matrices specified by  
us(trait):random.term capture how much variance there is in the counts  
of each species (along the diagonal) between villages, between  
children within villages, and between observations within children.  
The latter being the overdispersion (rcov). The off-diagonal elements  
measure the degree to which the counts of the two species covary  
(negatively) at these different levels.  However, you have fixed the  
overdispersion to one in the prior for each species, and fixed the  
residual correlation between species to zero. You need to remove the  
fix argument. Also nu=2 in your G-stricture priors may have quite a  
strong effect particularly for those groupings where there is little  
information (e.g. village). I tend to use parameter expanded priors  
instead.

Predicting "mixed clinical episodes" is a bit more tricky.  Since,  
E[PFepisodes*PVepisodes]  
=COV(PFepisodes,PVepisodes)-E[PFepisodes]*E[PVepisodes], your model  
already predicts the product. However,  COV(PFepisodes,PVepisodes) is  
a constant in your model (i.e. the between village covariance is not  
variable) so it is perhaps not as flexible as you would like.  I think  
you would like to add predictors for COV(PFepisodes,PVepisodes)(?) but  
this is not straightforward in MCMCglmm. Discretising cases into  
single and multiple infection and using a binary glmm might be easier.

Cheers,

Jarrod
Quoting Katie Colborn <benton at stat.berkeley.edu> on Tue, 27 Nov 2012  
08:35:46 -0500:

> Hi:
>
> Thanks for the helpful comments from my last post. I am writing  
> today with a new question.
>
> I'd like to specify a bivariate Poisson model with two random  
> effects. The responses are counts of clinical episodes of two  
> different malaria species during an 8-week interval. There are a  
> handful of time-dependent covariates that I would like to include on  
> the RHS of the equation as well as a random effect for child  
> (repeated measures) and one for village (clustering). There's also  
> the complication of needing to add an offset for time-at-risk which  
> subtracts off the time after treatment when a child is not at risk  
> of becoming ill (I will worry about that later). There are 1867  
> observations from 264 children nested within 11 villages.
>
> My question is, can I fit this model in MCMCglmm or GLMER? I don't  
> believe I can fit it in GLMER, but I should be able to fit it in  
> MCMCglmm. I'd like to fit a saturated model, so I'd like an estimate  
> of the effects of the covariates on "mixed clinical episodes", i.e.  
> the interaction between counts of episodes between the two species.  
> The two counts are negatively correlated, which I hope can be  
> estimated using this package as well. Does anyone know of a good  
> manuscript that has done something similar? I think what I am truly  
> stuck on is writing down the model. How many parameters should there  
> be on the RHS? An intercept and a slope for each fixed effect for  
> each species? For example, if there is simply age on the RHS, do I  
> end up with 4 parameters: 2 intercepts and 2 slopes for the two  
> response categories? What about Y1*Y2?
>
> My first attempt at analyzing this is below (which I obviously  
> should not do until I figure out how to write down the model, but  
> perhaps this will help to understand what I am trying to do):
>
> prior.bivar <- list(R = list(V = diag(2), nu = 0, fix=1),
> 		 G = list(G1 = list(V = diag(2), n=2),
> 		 G2 = list(V = diag(2), n=2)))
>
> model.bivar <- MCMCglmm(cbind(PFepisodes, PVepisodes) ~ -1 + trait +  
> trait:Age,
> 		random = ~ us(trait):child + us(trait):village,
> 		rcov = ~ us(trait):units, prior = prior.bivar,
> 		family = c("poisson", "poisson"),
> 		nitt = 60000, burnin = 10000,
> 		thin=25, data = MixedDat)
>
> summary(model.bivar)
>
> Thanks in advance!
>
> Katie Colborn
> benton at stat.berkeley.edu
>
>
> 	[[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.



More information about the R-sig-mixed-models mailing list