Skip to content

MCMCglmm troubles

3 messages · Jarrod Hadfield, Michael Lawrence

#
Hi folks,

I'm trying to do a mixed multinomial model with 3 response variables,
but I'm getting the error "please use idh() or us() error structure".
I believe that idh() and us() are used to specify the nature of the
interaction between multiple random variables, but I don't have
multiple random variables, so I'm not sure why I'm getting that error.
Some example data are below, but note that in my real data, the dv's
are indeed non-independent (they are in fact the number of trials
classified into each of three categories).

a = expand.grid(
	sid = 1:20
	, condition = factor(1:2)
)
a$group = factor(a$sid%%2)
a$sid = factor(a$sid)
a$dv1 = rnorm(nrow(a))
a$dv2 = rnorm(nrow(a))
a$dv3 = rnorm(nrow(a))

MCMCglmm(
	cbind(dv1,dv2,dv3)~condition*group
	, random = ~ sid
	, family = 'multinomial3'
	, data = a
)
6 days later
#
Hi Mike,

Sorry for the delay on this.  I'm not sure you want "multinomial3" as  
your family given that the three response variables are normally  
distributed. Specifying family=rep("gaussian", 3) would be more  
appropriate.  You will still get an error message because the default   
residual term is ~units which specifies a residual for each row of the  
response. For univariate models this is equivalent to IID residuals.  
For multi-response models it is more usual to allow different residual  
variances for each response, and often residual covariances between  
responses. rcov=~idh(trait):units fits different variances across  
responses (trait)  but sets the covariances to zero. rcov=  
~us(trait):units  estimates the covariances too. There are other less  
general possibilities too (see CourseNotes) which may be usefull. For  
example  rcov=~trait:units fits a common variance across responses and  
fixes the covariances to zero.

Having random=~sid also assumes that each sid effect is identical  
across responses so you may want to relax this assumption as above.

Cheers,

Jarrod
On 17 Jun 2010, at 13:43, Mike Lawrence wrote:

            

  
    
#
Ah, that makes sense now. Oh, and I think it is indeed the
"multinomial3" family that I need to model; I used rnorm() as a
shortcut to generate data to show the error I encountered, but the
real data are counts.
On Thu, Jun 24, 2010 at 5:46 AM, Jarrod Hadfield <j.hadfield at ed.ac.uk> wrote: