
Wenceslau, Na formula do modelo você deve remover o intercepto (-1), deixar os blocos como último termo e usar contraste tipo soma zero para que seus coeficientes se anulem e não precisem ser representados, declarar os demais fatores com estrutura aninhada, tempo dentro de dose, veja da <- expand.grid(B=gl(4,1,la="b"), D=gl(3,1,la="d"), T=1:10) da$y <- rnorm(nrow(da)) m0 <- lm(y~-1+D/(T+I(T^2))+B, da, contrast=list(B=contr.sum)) summary(m0) Se o seu experimento você montou em blocos, então bloco deve ser termo de todos os modelos, ok, no último você só declara D e T, deixou o bloco fora. À disposição. Walmes. ========================================================================== Walmes Marques Zeviani LEG (Laboratório de Estatística e Geoinformação, 25.450418 S, 49.231759 W) Departamento de Estatística - Universidade Federal do Paraná fone: (+55) 41 3361 3573 VoIP: (3361 3600) 1053 1173 e-mail: walmes@ufpr.br twitter: @walmeszeviani homepage: http://www.leg.ufpr.br/~walmes linux user number: 531218 ========================================================================== 2011/5/4 Wenceslau <wgt007@gmail.com>
Caros colegas,
estou tentando fazer uma analise comparando modelos lineares da evolução do pH vs tempo, num experimento com diferentes doses de carbonato de cálcio (7 doses).
Tenho este Tenho este script que faz a analise e também os gráficos. Mas não consegui obter os coeficientes de cada modelo individualmente (dose 0 vs tempo; dose 4 vs tempo, etc).
Agradeço se puderem me ajudar
## Experimento pH vs tempo
dados <- read.table("dosagem.txt", h=T) names(dados) attach(dados)
B <- factor(Rep) D <- factor(Dose) T <- factor(Tempo)
modelo <- lm(Leitura ~ D*T+B:D) anova(modelo)
# Pode-se considerar o modelo fatorial diretamente
fatorial <- lm(Leitura ~ D*(I(Tempo)+I(Tempo^2))) anova(fatorial) summary(fatorial) #plot(fatorial)
#SELECIONANDO O MODELO QUADRATICO E GERANDO INTERVALOS DE CONFIANÇA DE 0.95
preditos <- predict(fatorial,interval="confidence")
media <- preditos[,1] LI <- preditos[,2] LS <- preditos[,3]
interaction.plot(Tempo,D,media,ylab="pH",xlab="Dias", ylim=c(4,8), pch=0:7, legend=T)
#Plotando os pontos e as curvas ajustadas - Grafico 1 par(new=TRUE) interaction.plot(Tempo,D,Leitura,ylab="pH",xlab="Dias", type="p",pch=0:7,col=3:9, ylim=c(4,7.5),legend=T) par(new=TRUE) interaction.plot(Tempo,D,media,ylab="pH",xlab="Dias", ylim=c(4,7.5), pch=0:7, legend=F)
#Plotando a curva ajustada e os IC
par(new=TRUE) interaction.plot(Tempo,D,ylab="pH", media,ylim=c(4,7.5),legend=F) par(new=TRUE)
interaction.plot(Tempo,D,ylab="pH",xlab="Dias",LI,col=2,ylim=c(4,7.5),legend=F) par(new=TRUE) interaction.plot(Tempo,D,ylab="pH",xlab="Dias", LS,col=3,ylim=c(4,7.5),legend=F) legend("topleft", c("IC - Superior", "IC - Inferior"), col=c("green", "red"), title="Legenda", lty=1, lwd=2)
_______________________________________________ R-br mailing list R-br@listas.c3sl.ufpr.br https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br