#
#   Experimento em Blocos Completamente Aleatorizados
#     Três soluções Enxaguantes aleatorizadas por dia (em 4 dias)
#      Resposta: tempo até aparecimento de bactérias (horas)
#
 rm(list=ls(all=TRUE))
#
#############################################################
# instalando alguns pacotes necess?rios para a an?lise:
#
#install.packages("multcomp", dependencies=TRUE)
#install.packages("multcompView", dependencies=TRUE)
#install.packages("agricolae", dependencies=TRUE)
#install.packages("nortest", dependencies=TRUE)
#install.packages("car", dependencies=TRUE)
#
# carrecando os pacotes instalados:
#
# require(multcomp)
# require(multcompView)
 require(agricolae)
# require(nortest)
 require(car)
#
#############################################################
#
# lendo o banco de dados:
 dados=read.table("dados1_bloco.txt", header=TRUE)
 dados
#
#verificando a classe de cada vari?vel do nosso objeto dados:
 sapply(dados,class)
#
# obtendo as m?dias amostrais associadas a cada n?vel do fator de interesse:
 medias.trat <- tapply(dados$resp, dados$trat, mean)
 medias.bloco<- tapply(dados$resp, dados$bloco, mean)
 medias.trat
 medias.bloco
#
# explorando graficamente:
 plot(resp~trat,data=dados)
 dev.new()
 plot(resp~bloco, data=dados)
#
# Obtendo a tabela ANOVA:
 modelo=aov(resp~trat+bloco,data=dados)
 summary(modelo)
 model.tables(modelo, ty="means")
#
# valor tabelado da distribui??o F:
 qf(0.05, 3, 6, lower.tail = FALSE)
#
 class(modelo)
 names(modelo)
#
##########################################################################
##########################################################################
# Análise de resíduos
#
 trat.hat=modelo$coeff
 trat.hat   # coeficientes estimados
#
 y.hat=modelo$fitted.values
 y.hat     # valores preditos
#
 residuo=modelo$residuals
 residuo   # resíduos
#
############################################
#
# Testes de homogeneidade da variância
 bartlett.test(resp~trat, data=dados)
 leveneTest(resp~trat, data=dados)
#
 require(lawstat)  # pedir o pacote lawstat
 levene.test(dados$resp, dados$trat)
# Testes de normalidade:
 shapiro.test(residuo)
# ad.test(residuo)
#
  par(mfrow=c(2,2))
 plot(modelo)
 qqnorm(residuo)
 qqline(residuo, col = 2)
#
 win.graph()
 hist(residuo)
#
 win.graph()
 plot(residuo)
#
#########################################################################
# Comparações múltiplas:
#
# Como as suposições foram atendidads, necessitamos
# de comparações múltiplas para encontrar quais grupos são diferentes
# e quais são iguais. Vamos utilizar os métodos de Bonferroni e Tukey.
#
#  pairwise.t.test(dados$resp, dados$trat, p.adj = "bonf")
#
 ex01.tu <- TukeyHSD(modelo,dados$trat,conf.level=0.95)
 plot(ex01.tu)
 ex01.tu
#
# Concluimos que existe diferença média de cerca de 15 horas entre o 1 e 2 com relação
# ao terceiro
#
# Teste da m?nima diferen?a significativa de Fisher
# print(LSD.test(modelo, "trat", group=FALSE))
# print(LSD.test(modelo, "trat", group=TRUE))
#
# Teste de Bonferroni
 print(LSD.test(modelo, "trat", group=FALSE, p.adj="bonferroni"))
 print(LSD.test(modelo, "trat", group=TRUE, p.adj="bonferroni"))
#
# teste da m?xima diferen?a significativa de Tukey
 print(HSD.test(modelo, "trat", group=FALSE))
 print(HSD.test(modelo, "trat", group=TRUE))
#
#
# fun??o alternativa para o teste de Tukey:
 dif=TukeyHSD(modelo, dados$trat,ordered = TRUE)
 dif
 win.graph()
 plot(dif)
#
##############
#
# mcHSD=glht(modelo, linfct = mcp(trat = "Tukey"))
#summary(mcHSD)
#comp=LSD.test(modelo, "trat")
#
############################################
#
#   Quadrado Latino
#
#    Tempo de Montagem de um componente elétrico
#      Y: tempo (minutos)
#      Fator: 4 métodos de montagem
#      Bloco1: operador
#      Bloco2: ordem de montagem
################################################################
#############################################################
#
# lendo o banco de dados:
 dados1=read.table("quadrado_latino.txt", header=TRUE)
 dados1
 names(dados1)
#
#verificando a classe de cada vari?vel do nosso objeto dados:
 sapply(dados1,class)
 dados1$bloco1<-as.factor(dados1$bloco1)
 dados1$bloco2<-as.factor(dados1$bloco2)
#
# obtendo as m?dias amostrais associadas a cada n?vel do fator de interesse:
 medias.trat <- tapply(dados1$tempo, dados1$fator, mean)
 medias.bloco1<- tapply(dados1$tempo, dados1$bloco1, mean)
 medias.bloco2<- tapply(dados1$tempo, dados1$bloco2, mean)
 medias.trat
 medias.bloco1
 medias.bloco2
#
# explorando graficamente:
 plot(tempo~fator,data=dados1)
 dev.new()
 plot(tempo~bloco1, data=dados1, xlab="operador")
 plot(tempo~bloco2, data=dados1, xlab="ordem")
#
#
# Obtendo a tabela ANOVA:
 modelo1=aov(tempo~fator+bloco1+bloco2,data=dados1)
 summary(modelo1)
 model.tables(modelo1, ty="means")
#
 par(mfrow=c(2,2))
 plot(modelo1)
# Testes de homogeneidade da variância
 bartlett.test(tempo~fator, data=dados1)
 leveneTest(tempo~fator, data=dados1)
#
# Testes de normalidade:
 shapiro.test(residuals(modelo1))
#
#  Comparações Múltiplas
#
 ex02.tu <- TukeyHSD(modelo1,dados1$fator,conf.level=0.95)
 plot(ex02.tu)
# 
