[R-sig-ME] Calculating repeatability with ordinal data
Samantha Patrick
spatrick at cebc.cnrs.fr
Tue May 22 16:22:26 CEST 2012
Hi
Thanks for the comments - sorry for the delayed response, I had to learn
how to use MCMCglmm!
Following Shinichi's advice I have been using an adapted version of the
code in Dean et al 2012 (Am Nat). To start with I have been trying to
fit the data as categorical, as the original script is based on this
distribution.
Simplified code:
Bird = factor with 1183 levels with 1638 observations
Year = factor with 4 levels
scoremax = ordinal factor with 5 levels (0-4)
year<-as.factor(year)
# setting prior probabilities for random effects
IJ<-(1/5)*(diag(4) + matrix(1,4,4))
priortest1 <- list(R=list(V = IJ, nu = 0.002), G = list(G1 = list(V = 1,
nu = 0.002), G2 = list(V = 1, nu = 0.002)))
model1b<-MCMCglmm(scoremax~trait-1, random =~Bird +year, rcov =
~us(trait):units, data = Data, prior = priortest1, family = "categorical")
# Getting "adjusted" repeatability - Equation A10
summary(model1$VCV[,1]/(model1$VCV[,1]+2/3+pi^2/3))
posterior.mode(model1$VCV[,1]/(model1$VCV[,1]+2/3+pi^2/3))
This does run but the model fit is pretty bad. I haven't tried to
optimise it as it is important that my data is ordered so I am trying to
work on the script below:
priortest2 <- list(R=list(V = 1 , fix=1), G = list(G1=list(V=1, nu=1,
alpha.mu=0, alpha.V=100),G2=list(V=1, nu=1, alpha.mu=0, alpha.V=100)))
## this seems to be a commonly used prior for ordinal data so I was
using it as a starting point but have not got as far as to test whether
it is more suitable than priortest1 due to the problem below
model1b<-MCMCglmm(scoremax~trait-1, random =~Bird + year, rcov =
~us(trait):units, data = Data, prior = priortest2, family = "ordinal")
At this stage I get an error message which says:
Error in `contrasts<-`(`*tmp*`, value = contr.funs[1 + isOF[nn]]) :
contrasts can be applied only to factors with 2 or more levels
I have tried to create dummy data to seek out the problem. I thought
the problem may stem from the model trying to estimate contrasts between
each level of the ordinal response ie. 0-1, 1-0, 0-2, 2-0, 0-3, 3-0, 0-4
etc. Within the data individuals move between all scores but not both
ways (i.e. There are individuals who change between 1 and 2, but not
necessarily individuals that move from 1 to 2 AND 2 to 1). However, I
have not fitted observation number in the model (the temporal order of
observations is not specified) and altering the order of observation
does not seem to fix this. I have simplified the model and run without
year, but I get the same error.
I can happily post more data if it is helpful.
Any advice would be greatly appreciated!
Thanks
Sam
Le 02/05/2012 15:43, Jarrod Hadfield a écrit :
> Hi,
>
> The residual variance can never be estimated from ordinal data. Most
> programs will set it to zero, some programs allow you to set it at
> anything (MCMCglmm). I have not seen repeatabilities for ordinal data
> but I presume you can add the variance of the relevant distribution
> (pi^2/3 or 1 for logit/probit) to the denominator in order to get an
> intra-class correlation. I'm not confident how this would be
> interpreted for an ordinal response though. Note that if you do use
> MCMCglmm you need to include the constrained residual variance in the
> denominator (i.e. to get the denominator you need to add two to the
> "Bird" variance if you have constrained the residual variance to one).
>
> Cheers,
>
> Jarrod
>
>
> Quoting Samantha Patrick <spatrick at cebc.cnrs.fr> on Wed, 02 May 2012
> 13:11:09 +0200:
>
>> Hi
>>
>> Thank you for the tip - I had editted my code! However, sadly it
>> doesn't solve the problem of how to estimate the residual variance.
>>
>> Many Thanks
>>
>> Sam
>>
>>
>> Le 02/05/2012 12:38, Federico Calboli a écrit :
>>> On 2 May 2012, at 11:08, Samantha Patrick wrote:
>>>
>>>> Hi
>>>>
>>>> Firstly I sent this message last week but can find no evidence that
>>>> it actually sent.
>>>> However, if this is a double posting I am very sorry.
>>>>
>>>> I am currently working with repeated measures for individuals and I am
>>>> trying to quantify individual repeatability. Normally, for continuous
>>>> distributions, I would use a mixed model and calculate the variance
>>>> explained by individual divided by the total variance. However my
>>>> individual scores
>>>> are ordinal and I have been using the clmm function in the Ordinal
>>>> package:
>>>>
>>>> example of data:
>>>> ID Bird Sex scoremax year
>>>> 622 BS8831 M 2 2008
>>>> 623 BS8831 M 1 2010
>>>> 624 BS8831 M 1 2011
>>>> 625 BS9065 M 1 2010
>>>> 626 BS9065 M 3 2011
>>>> 627 BS19724 F 4 2010
>>>> 628 BS19724 F 5 2010
>>>> 629 BS21302 F 1 2010
>>>> 630 BS25376 F 1 2011
>>>> 631 BS9184 F 2 2009
>>>> 632 BS19989 M 3 2011
>>>> 633 BS21617 M 4 2008
>>>> 634 BS21617 M 2 2009
>>>> 635 BS21617 M 1 2010
>>>>
>>>> where scoremax ranges from 1-5, and there are 1188 birds and 1638
>>>> observations.
>>>>
>>>> scoremax<-as.factor(scoremax)
>>> does
>>>
>>> scoremax = as.ordered(scoremax)
>>>
>>> make any difference in your results?
>>>
>>> I ask this because as.factor() does not, strictly speaking, create
>>> an ordered factor.
>>>
>>> BW
>>>
>>> F
>>>
>>>
>>>
>>>
>>>> bird<-as.factor(bird)
>>>> year<-as.factor(year)
>>>> fmm1<- clmm(scoremax~year+ (1|bird), link = c("probit"), Hess =TRUE)
>>>>
>>>> summary(fmm1)
>>>>
>>>> but this only gives the variance estimate for bird, with no residual
>>>> estimate. Some investigations reveal that using an ordinal regression
>>>> in MCMCglmm will also not estimate the residual variance, and it seems
>>>> you need to constrain this value. I have been unable to find any
>>>> posts
>>>> about repeatability in ordinal data.
>>>>
>>>> My questions I guess are:
>>>> Is using a mixed model appropriate for calculating the
>>>> repeatability of
>>>> ordinal data (and if not does anyone know any other methods)?
>>>>
>>>> If it is, does anyone have any hints on how to calculate the
>>>> residual variance,
>>>> to enable repeatability estimates to be calculated.
>>>>
>>>> Many Thanks
>>>>
>>>> Sam
>>>>
>>>> --
>>>>
>>>> Dr Samantha Patrick
>>>> Post Doctoral Fellow
>>>> Centre d'Etudes Biologiques de Chizé - CNRS
>>>> 79360 Villiers-en-Bois
>>>> France
>>>> T:+33 549 097 846
>>>> M:+33 675 603 451
>>>> Skype: sammy_patrick
>>>> http://www.cebc.cnrs.fr/Fidentite/patrick/patrick.htm
>>>> http://www.researchgate.net/profile/Samantha_Patrick/
>>>>
>>>> _______________________________________________
>>>> R-sig-mixed-models at r-project.org mailing list
>>>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>> --
>>> Federico C. F. Calboli
>>> Neuroepidemiology and Ageing Research
>>> Imperial College, St. Mary's Campus
>>> Norfolk Place, London W2 1PG
>>>
>>> Tel +44 (0)20 75941602 Fax +44 (0)20 75943193
>>>
>>> f.calboli [.a.t] imperial.ac.uk
>>> f.calboli [.a.t] gmail.com
>>>
>>>
>>>
>>
>> --
>>
>> Dr Samantha Patrick
>> Post Doctoral Fellow
>> Centre d'Etudes Biologiques de Chizé - CNRS
>> 79360 Villiers-en-Bois
>> France
>> T:+33 549 097 846
>> M:+33 675 603 451
>> Skype: sammy_patrick
>> http://www.cebc.cnrs.fr/Fidentite/patrick/patrick.htm
>> http://www.researchgate.net/profile/Samantha_Patrick/
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>>
>
>
>
--
Dr Samantha Patrick
Post Doctoral Fellow
Centre d'Etudes Biologiques de Chizé - CNRS
79360 Villiers-en-Bois
France
T:+33 549 097 846
M:+33 675 603 451
Skype: sammy_patrick
http://www.cebc.cnrs.fr/Fidentite/patrick/patrick.htm
http://www.researchgate.net/profile/Samantha_Patrick/
More information about the R-sig-mixed-models
mailing list