source("config/setup.R")
#----------------------------------------------------------------------- # Carrega os pacotes necessários. rm(list = ls()) library(lattice) library(latticeExtra) library(captioner) library(reshape) library(doBy) library(multcomp) library(splines) library(wzRfun) library(EACS)
# Estrutura dos dados. str(biocar_eucal) # Usa objeto de nome curto mais fácil de manipular. bio <- biocar_eucal # Tabela de frequência do desenho experimental. xtabs(~dosep + biocar, data = bio) # Verifica se alguma variável resposta tem valor perdido. sapply(bio[, -(1:3)], FUN = function(x) any(is.na(x)))
# Empilha a variável resposta para fazer gráficos. bioe <- melt(bio, id.vars = 1:3, measure.vars = 6:ncol(bio), variable_name = "resp") xyplot(value ~ dosep^0.3 | resp, groups = biocar, data = bioe, as.table = TRUE, type = c("p", "smooth"), auto.key = list(title = "Biocarvão", cex.title = 1.1, columns = 3), scales = list(y = list(relation = "free")))
# Níveis únicos de dose de P. x <- sort(unique(bio$dosep)) x # Variância das distância entre níveis em escala unitária. esp <- function(p) { u <- x^p u <- (u - min(u)) u <- u/max(u) var(diff(u)) } # Otimiza para obter o valor de potência para maior uniformidade. op <- optimize(f = esp, interval = c(0, 1)) op$minimum p <- seq(0, 1, by = 0.01) v <- sapply(p, esp) plot(log(v) ~ p, type = "o") abline(v = op$minimum) # Cria potência da dose. bio$dose <- bio$dosep^0.3 bio$P <- factor(bio$dosep)
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(altfin ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 4)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(altfin ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(diafin ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(diafin ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(mfpa ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(mfpa ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(mfr ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(mfr ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(mspa ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(mspa ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(msr ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(msr ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(ph ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(ph ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(phkcl ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(phkcl ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(pcz ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(pcz ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(log(p) ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # MASS::boxcox(m0) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(log(p) ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(k ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(k ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(ca ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(ca ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(mg ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(mg ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm((hal)^0.33 ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # MASS::boxcox(m0) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(hal^0.33 ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(are ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(are ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(sil ~ bloc + biocar * P, data = bio[bio$sil < 30, ]) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(sil ~ dose, groups = biocar, data = bio[bio$sil < 30, ], auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
#----------------------------------------------------------------------- # Modelo saturado (considera dose como qualitativo). m0 <- lm(arg ~ bloc + biocar * P, data = bio) par(mfrow = c(2, 2)) plot(m0) layout(1) # Quadro de ANOVA. anova(m0) L <- LE_matrix(m0, effect = c("biocar", "P")) pred <- equallevels(attr(L, "grid"), bio) ci <- confint(glht(m0, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) str(pred) pch <- c(1, 4, 5) key <- list(title = "Biocarvão", cex.title = 1.1, text = list(levels(bio$biocar)), lines = list(pch = pch, lty = 1), columns = 3, type = "o", divide = 1) segplot(P ~ lwr + upr, centers = Estimate, data = pred, draw = FALSE, groups = biocar, xlab = "Dose de P", ylab = "Altura final (cm)", key = key, horizontal = FALSE, pch = pch, gap = 0.1, panel = panel.groups.segplot) #----------------------------------------------------------------------- # Modelo reduzido baseado em splines. m1 <- update(m0, ~ bloc + biocar * bs(dose, df = 3)) par(mfrow = c(2, 2)) plot(m1) layout(1) anova(m1) # Teste para a falta de ajuste. anova(m1, m0) # Novos valores em grid para fazer o gráfico de bandas. pred <- with(bio, expand.grid(bloc = factor(levels(bloc)[1], levels = levels(bloc)), biocar = levels(biocar), dose = seq(min(dose), max(dose), length.out = 60))) # Extrai a especificação do tempo base splines. form <- attr(terms(formula(m1)), "term.labels")[3] B <- eval(parse(text = form), envir = bio) str(B) # Gera a matriz para os valores de predição. X <- predict(B, pred$dos) # Cria a matriz do modelo para fazer a predição. L <- model.matrix(~bloc + biocar * X, data = pred) # Dá pesos para os efeitos de blocos correspondente à LSmeans. i <- attr(L, "assign") == 1 L[, i] <- 1/(sum(i) + 1) # Obtém os valores preditos com IC. ci <- confint(glht(m1, linfct = L), calpha = univariate_calpha())$confint # Junta as colunas estimadas com os valores do grid. pred <- cbind(pred, as.data.frame(ci)) xyplot(arg ~ dose, groups = biocar, data = bio, auto.key = list(points = FALSE, lines = TRUE, divide = 1, type = "o")) + as.layer(xyplot(Estimate ~ dose, groups = biocar, data = pred, type = "l", ly = pred$lwr, uy = pred$upr, cty = "bands", alpha = 0.3, prepanel = prepanel.cbH, panel.groups = panel.cbH, panel = panel.superpose))
TODO incluir o quadro de anova com todas as respostas.
cat(format(Sys.time(), format = "Atualizado em %d de %B de %Y.\n\n")) sessionInfo()
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.