MCMCglmm 'undefined columns' error

Hi Greg,

There are two ways of fitting multinomial models when the number of  
observations per unit is 1 (as in your case). You specified  
family="multinomial4" so MCMCglmm is expecting a 4 column response  
variable with the counts in each column. In your case each row should  
contain 3 zeros and 1 one.  The response should be specified something  
like:

cbind(category1,category2,category3,category4)~

Alternatively, you can pass a single column of factors as a response  
and use family="categorical" and keep your original  syntax:

count ~

Both models are parameterised as a log-odds ratio against a baseline  
category. In the first model  the final category is the baseline  
category. I've done this because for binomial traits ("multinomial2")  
the standard glm syntax  is to pass cbind(successes, failures) and the  
results are the log odds of succeeding vs failing (ie failing is the  
base-line category).   In the second model  the first factor level is  
the baseline category.  I've done this because for binary traits it is  
natural to have zero's as the base-line and interpret the parameters  
as predicting 1's (successes).

For your data there will be three response variables in both cases  
({cat1-cat4, cat2-cat4, cat3-cat4} or {cat2-cat1, cat3-cat1, cat4- 
cat1}).   You will want to use the reserved variable "trait" so that  
different intercepts and effects are fitted for the three responses.  
For example count~trait-1 or possibly count~trait+trait:cane-1 etc.

A multivariate error structure must also be defined. You could use:

~ idh(trait):units or ~ us(trait):units

However, for multinomial models with a single count you can't estimate  
these from the data so you may as well fix them.  I recommend fixing  
the residual variances to 1, and the covariances to either 0 or 0.5.   
This is specified in the prior:

prior=list(R=list(V=diag(3), n=1, fix=1), G=Gprior)

where V=diag(3) is  3X3 (number of categories -1) identity matrix and  
fix=1 fixes the matrix at these values. The random effects should also  
have priors which are also passed as a list (Gprior)

Its up to you what type of structure you use for vine. The vine  
variance structure is also 3X3 and you could use any of these options:

~ vine                             # 1 parameter    
v1=v2=v3=c12=c13=c23      Gprior=list(G1=list(V=diag(1), n=1))
~ vine:trait                     # 1 parameter   v1=v2=v3   
c12=c13=c23=0  Gprior=list(G1=list(V=diag(1), n=1))
~ vine+vine:trait           # 2 parameters v1=v2=v3  c12=c13=c23        
Gprior=list(G1=list(V=diag(1), n=1)), G2=list(V=diag(1), n=1))
~idh(vine)                     # 3 parameters  v1 v2 v3     
c12=c13=c23=0  Gprior=list(G1=list(V=diag(3), n=3))
~idh(vine)+vine:trait   # 4 parameters  v1 v2 v3    c12=c13=c23        
Gprior=list(G1=list(V=diag(3), n=3), G2=list(V=diag(1), n=1))
~ us(vine)                     # 6 parameters  v1 v2 v3    c12 c13  
c23          Gprior=list(G1=list(V=diag(3), n=3))

where v1 v2 and v3 are the variances for each logg-odds ratio due to  
vine, and  c12 c13 c23  are the covarainces between them.  I've made  
up the prior (co)variances but I've put in the minimum degree of  
belief parameter (n) that is required to make the priors proper.

One obvious problem is that the multinomial model in this form looses  
the information about the ordering 0<1<2<3.

Cheers,

Jarrod

MCMCglmm 'undefined columns' error

Thread (2 messages)