#=========================================
#     ANÁLISE DE DADOS CATEGÓRICOS
#     Modelo Politômico Nominal
#     Projeto Juliana-Dani
# 
#    Resposta: Tipo de Reabsorção - Tipo1: 0- ausente, 1- substituição e 2-inflamtória
#     Objetivo: avaliar o efeito das variáveis:
#         1- Idade (usei os pontos de corte de 11 e 16 anos), 
#         2- SAT - Antibiótico
#         3- TemposplintLog, 
#         4- TempoCont1Log, 
#         5- MeioG e 
#         6-8 Genéticas ( RANKL, OPG e IL10) 
#                      
#================================================================================
# 1. Descrição do Banco de Dados 
#================================================================================
dados<-read.table("dados_dani.txt",h=T)
head(dados)
dim(dados)
names(dados)
#
#   Análise Descritiva das Variáveis
#
table(dados$Tipo1)
hist(dados$Idade, xlab = "Idade", ylab = "Frequência", main = "")
table(dados$SAT)
table(dados$TempoSplintLog)
table(dados$TempoCont1)
hist(dados$TempoSplint, xlab = "tempo splint", ylab = "Frequência", main = "")
hist(dados$TempoCont1, , xlab = "Tempo cont", ylab = "Frequência", main = "")
table(dados$MEioArmaz)
table(dados$MeioG)
table(dados$RankGen)
table(dados$OPGGen)
table(dados$IL10Gen)
tapply(dados$Idade,dados$Tipo1,summary)
#
#
#================================================================================
# 2. Regressão Logística Politômica
#       2 logits com categoria de baseline 0 (ausência de Reab.) p. 268 of Agresti, 
#       Package VGAM
#====================================================================
library(VGAM)
#
#  vglm trata o maior nível como referência
#    tomamos o maior nível como sendo ausência de reabsorção 
dados$Tipo1 <- ifelse(dados$Tipo1==0,3,dados$Tipo1)
class(dados$Tipo1)

#modelo
fit01 <- vglm(Tipo1~1,family="multinomial",data=dados)
summary(fit01)
#
# Ajuste com MeioG
fit1 <- vglm(Tipo1~factor(MeioG),family=multinomial,
             data=dados) # 
summary(fit1)
# Ajuste com TempoSplintLog
#fit2 <- vglm(Tipo1~factor(TempoSplintLog) ,family="multinomial",
#             data=dados) # 
#summary(fit2)
# Ajuste com SAT - antibiótico 
fit3 <- vglm(Tipo1~SAT,family="multinomial",data=dados) # 
summary(fit3)
#
# Ajuste com TempoCont ERRO (considerar TempoCont1Log)???
fit4 <- vglm(Tipo1 ~ TempoCont1, family="multinomial", data=dados) # 
summary(fit4)
# Ajuste com RANK
fit5<-vglm(Tipo1~factor(RanklGen), family = "multinomial",
           data=dados) # 
summary(fit5)
# Ajuste com OPGGen sig
fit6<-vglm(Tipo1~factor(OPGGen), family = "multinomial",
           data=dados) # 
summary(fit6)
# Ajuste com IL10
fit7<-vglm(Tipo1~factor(IL10Gen), family = "multinomial",
           data=dados) # 
summary(fit7)
#Ajuste com Idade
fit8 <- vglm(Tipo1~Idade, family = "multinomial",
             data=dados) #
summary(fit8)
# Ajuste com IdadeG11
fit9 <- vglm(Tipo1~IdadeG11, family = "multinomial",
             data=dados) # 
summary(fit9)
# Ajuste com IdadeG16
fit10 <- vglm(Tipo1~IdadeG16, family = "multinomial",
             data=dados) # 
summary(fit10)
# Ajuste com OPG e TempoSplintLog + IdadeG11
fit11<-vglm(Tipo1~IdadeG11 + factor(OPGGen) + log(TempoCont1)+
              factor(MeioG) + factor(RanklGen), family = "multinomial",
             data=dados) # 
summary(fit11)
#
#
fit13<-vglm(Tipo1~log(TempoCont1), family = "multinomial",
            data=dados) # 
summary(fit13)
#
#