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: [email protected] twitter: @walmeszeviani homepage: http://www.leg.ufpr.br/~walmes linux user number: 531218 ========================================================================== 2011/5/4 Wenceslau <[email protected]> > 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 > [email protected] > https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br > >
_______________________________________________ R-br mailing list [email protected] https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br
