ex <- read.table("dadosexemplo.txt", header=T)
ex

#fornecedores efeito fixo
#lotes efeito aleatório

dim(ex)
names(ex)

is.factor(ex$forn)
is.factor(ex$lot)
#forn e lot não foram lidos como fatores
#está lendo como quantidades numérica e
#não como indicadores dos níveis dos fatores

ex$forn <- as.factor(ex$forn)
ex$lot <- as.factor(ex$lot)

is.factor(ex$resp) 
is.numeric(ex$resp)# a variável é resposta é numérica

summary(ex) #resumo dos dados

#Deveríamos fazer alguns gráficos descritivos

#Ajuste do modelo
ex.av <- aov(resp ~ forn/lot, data=ex) # / indica efeitos aninhados
summary(ex.av)
#Embora os elementos do quadro de análise de variância estejam
#corretos o teste F para efeito dos fornecedores está ERRADO. 
#A análise acima considerou todos os efeitos como fixos e portanto 
#dividiu os quadrados médios dos efeitos pelo quadrado médio do resíduo.
#Como lotes é um efeito aleatório deveríamos dividir o quadrado médio 
#do termo lot pelo quadrado médio de forn:lot


ex.av1 <- aov(resp ~ forn/lot + Error(forn) , data=ex)
summary(ex.av1)
#Agora o teste F errado não é mais mostrado, mas o teste correto 
#também não foi feito!


ex.anova <- anova(ex.av)
is.list(ex.anova)

names(ex.anova)
ex.anova$Df
ex.anova$Mean
Fcalc <- ex.anova$Mean[1]/ex.anova$Mean[2]
Fcalc
pvalor <- 1 - pf(Fcalc, ex.anova$Df[1], ex.anova$Df[2])
pvalor

#Outra possível e elegante solução no R para este problema
#é utilizar a função lme do pacote nlme.
require(nlme)

#A seguir criamos uma variável para indicar o efeito aleatório
#que neste exemplo chamamos de ex$fa utilizando a função getGroups. 
ex$fa <- getGroups(ex, ~ 1|forn/lot, level=2)

#Ajuste do modelo
ex.av <- lme(resp ~ forn, random=~1|fa, data=ex)
ex.av
anova(ex.av)
intervals(ex.av)
#mostra os intervalos de confiança para os componentes de variância

summary(ex.av)