Exemplo: Introdução a análise de experimentos (cont.)
Neste outro exemplo ainda mostramos exemplos simplificados de recursos para análise de experimentos.
A ênfase aqui continua sendo ilustrar o uso da linguagem. Mais detalhes sobre recursos para análise serão fornecidos em outros tutoriais.
Experimento em blocos casualizados
Vamos agora analisar o experimento em blocos ao acaso. Os dados são fornecidos na tabela abaixo.
Conteúdo de óleo de S. linicola, em percentagem, em vários estágios de crescimento (Steel \& Torrie, 1980, p.199).
| Estágios | Blocos | |||
|---|---|---|---|---|
| I | II | III | IV | |
| Estágio 1 | 4,4 | 5,9 | 6,0 | 4,1 |
| Estágio 2 | 3,3 | 1,9 | 4,9 | 7,1 |
| Estágio 3 | 4,4 | 4,0 | 4,5 | 3,1 |
| Estágio 4 | 6,8 | 6,6 | 7,0 | 6,4 |
| Estágio 5 | 6,3 | 4,9 | 5,9 | 7,1 |
| Estágio 6 | 6,4 | 7,3 | 7,7 | 6,7 |
Inicialmente vamos entrar com os dados no R.
Há várias possíveis maneiras de fazer isto.
Vamos aqui usar a função scan() e entrar com os dados por linha da tabela.
Digitamos o comando abaixo e e função scan() recebe os dados. Depois de digitar o último dado digitamos
ENTER em um campo em branco e a função encerra a entrada de dados retornando para o prompt do programa.
OBS: Note que, sendo um programa escrito na língua inglesa, os decimais devem ser indicados por . e não por vírgulas.
y <- scan() #!1: 4.4 #!2: 5.9 #!3: 6.0 #!... #!24: 6.7 #!25: #!Read 24 items
Agora vamos montar um data.frame com os dados e os indicadores de blocos e tratamentos.
ex02 <- data.frame(estag = factor(rep(1:6, each=4)), bloco=factor(rep(1:4, 6)), resp=y)
Note que usamos a função factor para indicar que as variáveis blocos e estag são níveis de fatores e não
valores numéricos.
Vamos agora explorar um pouco os dados.
names(ex02) summary(ex02) attach(ex02) plot(resp ~ estag + bloco) interaction.plot(estag, bloco, resp) interaction.plot(bloco, estag, resp) ex02.mt <- tapply(resp, estag, mean) ex02.mt ex02.mb <- tapply(resp, bloco, mean) ex02.mb plot.default(estag, resp) points(ex02.mt, pch="x", col=2, cex=1.5) plot.default(bloco, resp) points(ex02.mb, pch="x", col=2, cex=1.5)
Nos gráficos e resultados acima procuramos captar os principais aspectos dos dados bem como verificar se não há interação entre blocos e tratamentos, o que não deve acontecer neste tipo de experimento.
A seguir vamos ajustar o modelo e obter outros resultados, incluindo a análise de resíduos e testes para verificar a validades dos pressupostos do modelo.
ex02.av <- aov(resp ~ bloco + estag)
anova(ex02.av)
names(ex02.av)
par(mfrow=c(2,2))
plot(ex02.av)
par(mfrow=c(2,1))
residuos <- (ex02.av$residuals)
plot(ex02$bloco,residuos)
title("Resíduos vs Blocos")
plot(ex02$estag,residuos)
title("Resíduos vs Estágios")
par(mfrow=c(2,2))
preditos <- (ex02.av$fitted.values)
plot(residuos,preditos)
title("Resíduos vs Preditos")
respad <- (residuos/sqrt(anova(ex02.av)$"Mean Sq"[2]))
boxplot(respad)
title("Resíduos Padronizados")
qqnorm(residuos,ylab="Residuos", main=NULL)
qqline(residuos)
title("Grafico Normal de \n Probabilidade dos Resíduos")
Teste para normalidade
shapiro.test(residuos)Testando a não aditividade: Primeiro vamos extrair coeficientes de tratamentos e blocos
ex02.av$coeff bl <- c(0, ex02.av$coeff[2:4]) tr <- c(0, ex02.av$coeff[5:9]) bl tr
e agora criar um novo termo e testar sua significancia na ANOVA
bltr <- rep(bl, 6) * rep(tr, rep(4,6)) ttna <- update(ex02.av, .~. + bltr) anova(ttna)
Os resultados acima indicam que os pressupostos estão obedecidos para este conjunto de dados e a análise de variância é válida. Como foi detectado efeito de tratamentos vamos proceder fazendo um teste de comparações múltiplas e encerrar as análises desanexando o objeto do caminho de procura.
ex02.tk <- TukeyHSD(ex02.av, "estag", ord=T) ex02.tk plot(ex02.tk) detach(ex02)
