#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: ivanala...@yahoo.com.br/ivanala...@gmail.com
@
\end{signature}_______________________________________________
R-br mailing list
R-br@listas.c3sl.ufpr.br
https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça código
mínimo reproduzível.
