
#Boa noite senhores! Após vários dias estudando sobre dados categóricos e pesquisando (vários fóruns com respostas do dono da library - todas incipientes) como interpretar o 'output' da função MCMCglmm venho recorrer a ajuda dos senhores, caso alguém já tenha trabalhado com função MCMCglmm no caso de variáveis dependentes ordinais. Já li e reli todos os pdf's do pacote, coursenotes do Jarrod, enfim, estou esgotado. Para esclarecer a base de dados, os tratamentos (denominado FASES) consistem de 3 níveis (1-pre, 2-própolis, 3-vincris), cujo as variável resposta ordinal tem 3 categorias (1-normal, 2-agudo, 3-cronico). Vide tabela! du <- transform(read.table('http://dl.dropbox.com/u/33619290/Dados/Dtest.txt',h=T),FASES=factor(FASES),ALT.RENAIS=ordered(ALT.RENAIS)) summary(du) library('MCMCglmm') du <- subset(du, ALT.RENAIS != 'NA') tabela <- table(du[,c(2,4)]) tabela colnames(tabela) <- c('Normal','Aguda','Crônica') rownames(tabela) <- c('Pre','Propolis','Vincr') tabela #Um modelo misto seria: set.seed(1) mod1 <- MCMCglmm(ALT.RENAIS ~-1+FASES, random= ~ ANIMAIS, family='ordinal',pl=TRUE,data=du) summary(mod1) #Aí começa o sofrimento, pois a documentação é insuficiente neste caso. Segundo o próprio Jarrod (fóruns), as médias a posteriori dos coeficientes das covariáveis estão na escala probit. Segundo o meu estudo, estes coeficientes são os scores da distribuição normal padrão. Mais os escores não deveriam corresponder aos cutpoints? Neste caso, teríamos j(categorias da variável resposta)-1 cutpoints, ou seja, 2 cutpoints. O output só me apresenta um cutpoint. Como é possível então calcular as probabilidades com apenas um cutpoint? Segundo a documentação(Vignettes, pag 22), se P(y=k)=F(yk|l(vlatente),1) - F(yk-1|l,1), esse 'um' provavelmente seria a categoria 'um' da variável dependente? Enfim senhores, como posso extrair as probabilidades para as fases referentes a cada categoria da variável dependente? Desde já agradeço a atenção de todos. #embora um pouco lógico, os resultados são absurdos! latentv <- mean(mod1$Liab) cutpoint <- mean(mod1$CP) pnorm(-(latentv), 0, sqrt(2)) pnorm(cutpoint - (latentv),0, sqrt(2)) - pnorm((latentv),0, sqrt(2)) 1- pnorm(cutpoint - (latentv),0, sqrt(2)) #isso teria algum resultado lógico até certo ponto, mais seria apenas referente a categoria 1 da variável dependente. E as outras? bre <- c(mean(mod1$Liab),mean(mod1$Sol[,1]),mean(mod1$Sol[,2]),mean(mod1$Sol[,3])) pnorm(bre[2])-pnorm(bre[1]) pnorm(bre[3])-pnorm(bre[2]) pnorm(bre[4])-pnorm(bre[3])#probabilidade negativa nãooo! \begin{signature} <<>>= Prof. Dr. Ivan Bezerra Allaman Universidade Estadual de Santa Cruz Departamento de Ciências Exatas e Tecnológicas Ilhéus/BA - Brasil Fone: +55 73 3680-5076 E-mail: ivanalaman@yahoo.com.br/ivanalaman@gmail.com @ \end{signature}