Análise da qualidade de ovos
Maria Fernanda Burbarelli, Leonardo de Oliveira Seno & Walmes Marques Zeviani
Preparação e análise exploratória
library(nlme) # Modelos de efeitos mistos.
# Para não exibir o bloco de correlação entre estimativas dos efeitos
# fixos na saída do `summary()`. ATTENTION: tem que vir antes da
# `lmerTest`.
assignInNamespace("print.correlation",
function(x, title) return(),
ns = "nlme")
library(lmerTest) # Testes de hipótese para classe lme e similares.
library(emmeans) # Médias ajustadas.
library(tidyverse) # Manipulação e visualização de dados.
# Para apresentar letras ordenadas com o "a" vínculado a maior média.
numbers_to_letters <- function(.group) {
list_num <- trimws(.group) %>%
strsplit(split = "")
n <- max(as.integer(unlist(list_num)))
list_num %>%
lapply(FUN = function(x) sort(letters[n:1][as.integer(x)])) %>%
sapply(FUN = paste, collapse = "")
}
# Função para fazer os gráficos com anotações das médias.
gg_means <- function(data = grid,
mapping = list(x = "Mate"),
facets = ~Amb,
fmt = "%0.2f %s",
labs = list(x = "Mate",
y = "???",
subtitle = "Ambiente")) {
ggplot(data = data,
mapping = aes(x = eval(parse(text = mapping$x)),
y = emmean)) +
facet_grid(facets = facets) +
geom_point() +
geom_errorbar(mapping = aes(ymin = lower.CL, ymax = upper.CL),
width = 0.1) +
geom_text(mapping = aes(label = sprintf(fmt, emmean, cld)),
angle = 90,
vjust = 0,
nudge_x = -0.075) +
labs(x = labs$x, y = labs$y) +
ggtitle(label = NULL, subtitle = labs$subtitle)
}
# Texto para a legenda conforme cada fator.
legends <- list(Mate = "Concentração de mate",
Per = "Período de postura",
Amb = "Ambiente de armazenamento")
## Leitura da tabela de dados.
url <- "Dados níveis de erva-mate - PARA ENVIAR ESTATISTICA.xlsx"
tb <- gdata::read.xls(xls = url,
sheet = 2,
na.strings = ".",
check.names = FALSE)
tb <- tb[, nchar(names(tb)) > 1]
# Adequação dos nomes das variáveis.
names(tb) <- names(tb) %>%
str_to_lower() %>%
str_replace_all("[[:punct:]]", "")
str(tb)## 'data.frame': 360 obs. of 19 variables:
## $ mate : num 0 0 0 0 0 0 0 0 0 0 ...
## $ rep : int 1 2 3 4 5 6 1 2 3 4 ...
## $ per : int 1 1 1 1 1 1 2 2 2 2 ...
## $ amb : Factor w/ 3 levels "AM","FR","REF": 1 1 1 1 1 1 1 1 1 1 ...
## $ pi : num 60.6 50.6 55.9 62.6 55.5 ...
## $ gravesp : int 1050 1050 1060 1050 1050 1055 1050 1060 1055 1050 ...
## $ altalb : num 6.25 6.87 6.95 7.54 8.13 8.21 5.44 5.19 4.36 5.43 ...
## $ uh : num 79.1 86.5 85.5 87 92.2 ...
## $ altge : num 15.1 15.2 14.4 15.7 15.6 ...
## $ diamge : num 48.4 45.4 50 48.5 48 ...
## $ pesoge : num 17.2 14.8 14.8 16.8 15.8 ...
## $ pecasca : num 7.53 5.6 5.39 5.64 5.34 4.43 5.8 5.65 5.93 6.36 ...
## $ pealb : num 33.3 27.4 32.1 36.9 31 ...
## $ masovo : num 50.5 42.2 47 53.6 46.8 ...
## $ espcasca: num 0.87 0.93 0.94 0.82 0.83 0.82 1.06 0.92 0.72 1.05 ...
## $ l : num 62.1 56.9 57.5 63.7 60.1 ...
## $ a : num 1.87 2.62 2.76 1.51 1.82 1.71 2.17 1.29 1.71 2.26 ...
## $ b : num 50.8 58.3 57.3 57.8 52.8 ...
## $ ph : num 6.8 6.5 6.7 6.3 6.3 6.3 6.3 6.3 6.3 6.3 ...# Faz ordenação dos registros.
tb <- tb %>%
arrange(mate, rep, per, amb)
# Número de repetições de cada ponto experimental.
ftable(xtabs(~amb + per + mate, data = tb))## mate 0 1.5 3 4.5 6
## amb per
## AM 1 6 6 6 6 6
## 2 6 6 6 6 6
## 3 6 6 6 6 6
## 4 6 6 6 6 6
## FR 1 6 6 6 6 6
## 2 6 6 6 6 6
## 3 6 6 6 6 6
## 4 6 6 6 6 6
## REF 1 6 6 6 6 6
## 2 6 6 6 6 6
## 3 6 6 6 6 6
## 4 6 6 6 6 6# Versão empilhada ou verticalizada das respostas.
tb_long <- tb %>%
gather(key = "variable", value = "value", pi:ph)
str(tb_long)## 'data.frame': 5400 obs. of 6 variables:
## $ mate : num 0 0 0 0 0 0 0 0 0 0 ...
## $ rep : int 1 1 1 1 1 1 1 1 1 1 ...
## $ per : int 1 1 1 2 2 2 3 3 3 4 ...
## $ amb : Factor w/ 3 levels "AM","FR","REF": 1 2 3 1 2 3 1 2 3 1 ...
## $ variable: chr "pi" "pi" "pi" "pi" ...
## $ value : num 60.6 60.3 60.8 58.6 64.7 ...ggplot(data = tb_long,
mapping = aes(x = mate, y = value, color = factor(per))) +
facet_wrap(facets = ~variable, scale = "free_y") +
# geom_point() +
geom_jitter(width = 0.1, height = 0) +
stat_summary(mapping = aes(group = per),
fun.y = "mean",
geom = "line") +
labs(x = "Concentração de mate",
y = "Valor de cada resposta")
ggplot(data = tb_long,
mapping = aes(x = per, y = value, color = factor(mate))) +
facet_wrap(facets = ~variable, scale = "free_y") +
# geom_point() +
geom_jitter(width = 0.1, height = 0) +
stat_summary(mapping = aes(group = mate),
fun.y = "mean",
geom = "line") +
labs(x = "Período de postura",
y = "Valor de cada resposta")
# ATTENTION: valores atípicos para algumas variáveis.
tb$a[tb$a > 7.5] <- NA
tb$espcasca[tb$espcasca > 2] <- NA
tb$gravesp[tb$gravesp > 1150] <- NA
tb$diamge[tb$diamge < 25] <- NA
# Conta os valores presentes (tamanho de amostra).
tb %>%
group_by(mate, per) %>%
summarise_at(vars(-(1:4)),
function(x) sum(is.finite(x)))## # A tibble: 20 x 15
## # Groups: mate [5]
## mate per altalb uh altge diamge pesoge pecasca pealb masovo
## <dbl> <int> <int> <int> <int> <int> <int> <int> <int> <int>
## 1 0 1 18 18 18 18 18 18 18 18
## 2 0 2 18 18 18 18 18 18 18 18
## 3 0 3 18 18 18 18 18 18 18 18
## 4 0 4 18 18 18 18 18 18 18 18
## 5 1.5 1 18 18 18 18 18 18 18 18
## 6 1.5 2 18 18 18 18 18 18 18 18
## 7 1.5 3 18 18 18 18 18 18 18 18
## 8 1.5 4 18 18 18 18 18 18 18 18
## 9 3 1 18 18 18 18 18 18 18 18
## 10 3 2 18 18 18 18 18 18 18 18
## 11 3 3 18 18 18 18 18 18 18 18
## 12 3 4 18 18 18 18 18 18 18 18
## 13 4.5 1 18 18 18 18 18 18 18 18
## 14 4.5 2 18 18 18 18 18 18 18 18
## 15 4.5 3 18 18 18 18 18 18 18 18
## 16 4.5 4 18 18 18 18 18 18 18 18
## 17 6 1 18 18 18 18 18 18 18 18
## 18 6 2 18 18 18 18 18 18 18 18
## 19 6 3 18 18 18 17 18 18 18 18
## 20 6 4 18 18 18 18 18 18 18 18
## # … with 5 more variables: espcasca <int>, l <int>, a <int>, b <int>,
## # ph <int># Conta os valores ausentes.
tb %>%
group_by(mate, per) %>%
summarise_at(vars(-(1:4)),
function(x) sum(is.na(x)))## # A tibble: 20 x 15
## # Groups: mate [5]
## mate per altalb uh altge diamge pesoge pecasca pealb masovo
## <dbl> <int> <int> <int> <int> <int> <int> <int> <int> <int>
## 1 0 1 0 0 0 0 0 0 0 0
## 2 0 2 0 0 0 0 0 0 0 0
## 3 0 3 0 0 0 0 0 0 0 0
## 4 0 4 0 0 0 0 0 0 0 0
## 5 1.5 1 0 0 0 0 0 0 0 0
## 6 1.5 2 0 0 0 0 0 0 0 0
## 7 1.5 3 0 0 0 0 0 0 0 0
## 8 1.5 4 0 0 0 0 0 0 0 0
## 9 3 1 0 0 0 0 0 0 0 0
## 10 3 2 0 0 0 0 0 0 0 0
## 11 3 3 0 0 0 0 0 0 0 0
## 12 3 4 0 0 0 0 0 0 0 0
## 13 4.5 1 0 0 0 0 0 0 0 0
## 14 4.5 2 0 0 0 0 0 0 0 0
## 15 4.5 3 0 0 0 0 0 0 0 0
## 16 4.5 4 0 0 0 0 0 0 0 0
## 17 6 1 0 0 0 0 0 0 0 0
## 18 6 2 0 0 0 0 0 0 0 0
## 19 6 3 0 0 0 1 0 0 0 0
## 20 6 4 0 0 0 0 0 0 0 0
## # … with 5 more variables: espcasca <int>, l <int>, a <int>, b <int>,
## # ph <int># Cria variável para representar a parcela e cria cópias como fator.
tb <- tb %>%
mutate(parcela = interaction(mate, rep, drop = TRUE),
subparc = interaction(mate, rep, per, drop = TRUE),
Mate = factor(mate),
Per = factor(per),
Amb = factor(amb))
#Especificação do modelo
O modelo correspondente ao experimento em parcelas subsubdivididas é \[ \begin{aligned} y_{ijkl}| ijkl &\sim \text{Normal}(\mu_{ijk}, \sigma^2) \\ \mu_{ijk} &= \mu + \alpha_i + e_{il} + \beta_j + \gamma_{ij} + e_{ijl} + \psi_k + \phi_{ik} + \kappa_{jk} + \lambda_{ijk}\\ e_{il} &\sim \text{Normal}(0, \sigma_{p}^2)\\ e_{ijl} &\sim \text{Normal}(0, \sigma_{s}^2)\\ \end{aligned} \]
em que * \(y_{ijkl}\) é o valor observado na repetição \(l\) (\(l = 1, \ldots, 6\)) ponto experimental \(ijk\) (\(i = 1, \ldots, 6\); \(j = 1, \ldots, 4\); \(k = 1, \ldots, 3\)), * \(\alpha_i\) é o conjunto de parâmetros que acomoda o efeito principal de concentração de mate, * \(\beta_j\) é o conjunto de parâmetros que acomoda o efeito principal de período de postura, * \(\psi_k\) é o conjunto de parâmetros que acomoda o efeito principal de ambiente de armazenamento, * \(\gamma_{ij}\) é o conjunto de parâmetros que acomoda o efeito de interação entre mate e período, * \(\phi_{ik}\) é o conjunto de parâmetros que acomoda o efeito de interação entre mate e ambiente, * \(\kappa_{jk}\) é o conjunto de parâmetros que acomoda o efeito de interação entre período e ambiente, * \(\lambda_{ijk}\) é o conjunto de parâmetros que acomoda o efeito de interação entre mate, período e ambiente, * \(e_{il}\) acomoda o erro experimental ao nível de parcelas e * \(e_{ijl}\) acomoda o erro experimental ao nível de subparcelas.
O modelo é ajustado aos dados pelo método da máxima verossimilhança restrita (REML).
Havendo interação entre duplas ou tripla, a prioridade para estudo dos efeitos será Mate > Ambiente > Período.
Análise da variável altalb
# Modelo de 1 extrato para avaliar os pressupostos.
m0 <- lm(terms(altalb ~
Mate + Mate:factor(rep) +
Per + Mate:Per + Mate:Per:factor(rep) +
Amb + Mate:Amb + Per:Amb + Mate:Per:Amb,
keep.order = TRUE),
data = tb)
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Apenas para verificação. Estatísticas F erradas pros estratos
# superiores.
anova(m0)## Analysis of Variance Table
##
## Response: altalb
## Df Sum Sq Mean Sq F value Pr(>F)
## Mate 4 10.131 2.533 2.7792 0.02803 *
## Mate:factor(rep) 25 24.609 0.984 1.0802 0.36813
## Per 3 131.225 43.742 48.0009 < 2.2e-16 ***
## Mate:Per 12 22.670 1.889 2.0732 0.02018 *
## Mate:factor(rep):Per 75 78.097 1.041 1.1427 0.23236
## Amb 2 5.703 2.851 3.1291 0.04591 *
## Mate:Amb 8 14.928 1.866 2.0477 0.04268 *
## Per:Amb 6 68.506 11.418 12.5294 5.618e-12 ***
## Mate:Per:Amb 24 26.833 1.118 1.2269 0.22210
## Residuals 200 182.253 0.911
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Ajuste do modelo de efeitos aleatórios.
mm0 <- lme(altalb ~ Mate * Per * Amb,
random = ~1 | parcela/subparc,
method = "ML",
na.action = na.omit,
data = tb)
# Quadro de testes para o termos de efeito fixo.
anova(mm0)## numDF denDF F-value p-value
## (Intercept) 1 200 17625.197 <.0001
## Mate 4 25 2.466 0.0710
## Per 3 75 42.589 <.0001
## Amb 2 200 3.129 0.0459
## Mate:Per 12 75 1.839 0.0567
## Mate:Amb 8 200 2.048 0.0427
## Per:Amb 6 200 12.529 <.0001
## Mate:Per:Amb 24 200 1.227 0.2221# Atualiza mantendo até os termos de dupla ordem.
# mm1 <- update(mm0, . ~ (Mate + Per + Amb)^2)
mm1 <- update(mm0, . ~ (Mate + Per) * Amb)
# Teste da razão de verossimilhanças (é anticonservador para teste de
# termos de efeitos fixos).
anova(mm1, mm0)## Model df AIC BIC logLik Test L.Ratio p-value
## mm1 1 27 1047.801 1152.726 -496.9006
## mm0 2 63 1062.903 1307.728 -468.4517 1 vs 2 56.89776 0.0147# Estimativas dos parâmetros.
summary(mm1)## Linear mixed-effects model fit by maximum likelihood
## Data: tb
## AIC BIC logLik
## 1047.801 1152.726 -496.9006
##
## Random effects:
## Formula: ~1 | parcela
## (Intercept)
## StdDev: 4.854454e-05
##
## Formula: ~1 | subparc %in% parcela
## (Intercept) Residual
## StdDev: 0.2405639 0.9333769
##
## Fixed effects: altalb ~ Mate + Per + Amb + Mate:Amb + Per:Amb
## Value Std.Error DF t-value p-value
## (Intercept) 7.037250 0.2576076 224 27.317712 0.0000
## Mate1.5 0.660417 0.2880140 25 2.293002 0.0305
## Mate3 0.536250 0.2880140 25 1.861888 0.0744
## Mate4.5 0.898333 0.2880140 25 3.119061 0.0045
## Mate6 0.298750 0.2880140 25 1.037276 0.3095
## Per2 -1.209333 0.2576076 87 -4.694479 0.0000
## Per3 0.686667 0.2576076 87 2.665553 0.0092
## Per4 -1.008000 0.2576076 87 -3.912928 0.0002
## AmbFR -0.015083 0.3527833 224 -0.042755 0.9659
## AmbREF -0.099083 0.3527833 224 -0.280862 0.7791
## Mate1.5:AmbFR -1.044583 0.3944237 224 -2.648379 0.0087
## Mate3:AmbFR -0.565417 0.3944237 224 -1.433526 0.1531
## Mate4.5:AmbFR -0.993750 0.3944237 224 -2.519499 0.0124
## Mate6:AmbFR -0.340833 0.3944237 224 -0.864130 0.3884
## Mate1.5:AmbREF -0.480833 0.3944237 224 -1.219078 0.2241
## Mate3:AmbREF -0.407083 0.3944237 224 -1.032096 0.3031
## Mate4.5:AmbREF -0.213333 0.3944237 224 -0.540873 0.5891
## Mate6:AmbREF 0.450000 0.3944237 224 1.140905 0.2551
## Per2:AmbFR 0.283000 0.3527833 224 0.802192 0.4233
## Per3:AmbFR 0.670000 0.3527833 224 1.899183 0.0588
## Per4:AmbFR 0.610667 0.3527833 224 1.730997 0.0848
## Per2:AmbREF 0.259333 0.3527833 224 0.735107 0.4630
## Per3:AmbREF -0.745333 0.3527833 224 -2.112723 0.0357
## Per4:AmbREF 1.749000 0.3527833 224 4.957718 0.0000
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -5.171890245 -0.625263160 0.002635295 0.619253133 4.158987211
##
## Number of Observations: 360
## Number of Groups:
## parcela subparc %in% parcela
## 30 120# Médias ajustadas.
emm <- emmeans(mm1, specs = ~Mate + Per + Amb)
ggplot(data = as.data.frame(emm),
mapping = aes(x = Amb,
y = emmean)) +
facet_grid(facets = Mate ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Ambiente de armazenamento",
y = "altalb") +
ggtitle(label = NULL,
subtitle = "Concentração de mate x período de postura")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Per,
y = emmean)) +
facet_grid(facets = Mate ~ Amb) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Período de postura",
y = "altalb") +
ggtitle(label = NULL,
subtitle = "Concentração de mate x ambiente de armazenamento")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Mate,
y = emmean)) +
facet_grid(facets = Amb ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Concentração de mate",
y = "altalb") +
ggtitle(label = NULL,
subtitle = "Ambiente de armazenamento x período de postura")
# Estudo da interação Per:Amb via Amb | Per.
emm <- emmeans(mm1, specs = ~Amb | Per)
grid <- multcomp::cld(emm)
grid## Per = 1:
## Amb emmean SE df lower.CL upper.CL .group
## FR 6.91 0.182 25 6.54 7.29 1
## REF 7.29 0.182 25 6.91 7.66 12
## AM 7.52 0.182 25 7.14 7.89 2
##
## Per = 2:
## Amb emmean SE df lower.CL upper.CL .group
## FR 5.99 0.182 25 5.61 6.36 1
## AM 6.31 0.182 25 5.93 6.68 1
## REF 6.34 0.182 25 5.96 6.71 1
##
## Per = 3:
## Amb emmean SE df lower.CL upper.CL .group
## REF 7.23 0.182 25 6.85 7.60 1
## AM 8.20 0.182 25 7.83 8.58 2
## FR 8.27 0.182 25 7.89 8.64 2
##
## Per = 4:
## Amb emmean SE df lower.CL upper.CL .group
## AM 6.51 0.182 25 6.13 6.88 1
## FR 6.51 0.182 25 6.14 6.89 1
## REF 8.03 0.182 25 7.65 8.40 2
##
## Results are averaged over the levels of: Mate
## d.f. method: containment
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 3 estimates
## significance level used: alpha = 0.05grid <- grid %>%
group_by(Per) %>%
mutate(cld = numbers_to_letters(.group)) %>%
ungroup()
gg_means(grid,
mapping = list(x = "Amb"),
facets = ~Per,
labs = list(x = legends$Amb,
y = "altalb",
subtitle = legends$Per))
# Estudo da interação Mate:Amb via Mate | Amb.
emm <- emmeans(mm1, specs = ~Mate | Amb)
grid <- multcomp::cld(emm)
grid## Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 0 6.65 0.204 29 6.24 7.07 1
## 6 6.95 0.204 25 6.53 7.37 12
## 3 7.19 0.204 25 6.77 7.61 12
## 1.5 7.32 0.204 25 6.90 7.73 12
## 4.5 7.55 0.204 25 7.13 7.97 2
##
## Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 6.65 0.204 25 6.23 7.07 1
## 4.5 6.93 0.204 25 6.52 7.35 1
## 6 6.99 0.204 25 6.57 7.41 1
## 3 7.00 0.204 25 6.58 7.42 1
## 0 7.03 0.204 29 6.61 7.45 1
##
## Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 0 6.87 0.204 29 6.45 7.29 1
## 3 7.00 0.204 25 6.58 7.42 1
## 1.5 7.05 0.204 25 6.63 7.47 1
## 4.5 7.56 0.204 25 7.14 7.98 1
## 6 7.62 0.204 25 7.20 8.04 1
##
## Results are averaged over the levels of: Per
## d.f. method: containment
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 5 estimates
## significance level used: alpha = 0.05grid <- grid %>%
group_by(Amb) %>%
mutate(cld = numbers_to_letters(.group)) %>%
ungroup()
gg_means(grid,
mapping = list(x = "Mate"),
facets = ~Amb,
labs = list(x = legends$Mate,
y = "altalb",
subtitle = legends$Amb))
#Análise da variável altge
# Modelo de 1 extrato para avaliar os pressupostos.
m0 <- lm(terms(altge ~
Mate + Mate:factor(rep) +
Per + Mate:Per + Mate:Per:factor(rep) +
Amb + Mate:Amb + Per:Amb + Mate:Per:Amb,
keep.order = TRUE),
data = tb)
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Apenas para verificação. Estatísticas F erradas pros estratos
# superiores.
anova(m0)## Analysis of Variance Table
##
## Response: altge
## Df Sum Sq Mean Sq F value Pr(>F)
## Mate 4 6.95 1.736 1.6516 0.16279
## Mate:factor(rep) 25 46.11 1.844 1.7543 0.01855 *
## Per 3 39.45 13.149 12.5077 1.573e-07 ***
## Mate:Per 12 53.33 4.444 4.2277 6.118e-06 ***
## Mate:factor(rep):Per 75 92.57 1.234 1.1741 0.19069
## Amb 2 157.80 78.900 75.0518 < 2.2e-16 ***
## Mate:Amb 8 43.12 5.390 5.1274 7.989e-06 ***
## Per:Amb 6 524.52 87.420 83.1564 < 2.2e-16 ***
## Mate:Per:Amb 24 69.82 2.909 2.7673 5.618e-05 ***
## Residuals 200 210.25 1.051
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Ajuste do modelo de efeitos aleatórios.
mm0 <- lme(altge ~ Mate * Per * Amb,
random = ~1 | parcela/subparc,
method = "ML",
na.action = na.omit,
data = tb)
# Quadro de testes para o termos de efeito fixo.
anova(mm0)## numDF denDF F-value p-value
## (Intercept) 1 200 46032.19 <.0001
## Mate 4 25 0.94 0.4564
## Per 3 75 10.65 <.0001
## Amb 2 200 75.05 <.0001
## Mate:Per 12 75 3.60 0.0003
## Mate:Amb 8 200 5.13 <.0001
## Per:Amb 6 200 83.16 <.0001
## Mate:Per:Amb 24 200 2.77 0.0001mm1 <- mm0
# Estimativas dos parâmetros.
summary(mm1)## Linear mixed-effects model fit by maximum likelihood
## Data: tb
## AIC BIC logLik
## 1131.306 1376.13 -502.6529
##
## Random effects:
## Formula: ~1 | parcela
## (Intercept)
## StdDev: 0.2058209
##
## Formula: ~1 | subparc %in% parcela
## (Intercept) Residual
## StdDev: 0.2254629 0.93598
##
## Fixed effects: altge ~ Mate * Per * Amb
## Value Std.Error DF t-value p-value
## (Intercept) 15.191667 0.4402850 200 34.50417 0.0000
## Mate1.5 0.806667 0.6226570 25 1.29552 0.2070
## Mate3 -0.053333 0.6226570 25 -0.08565 0.9324
## Mate4.5 0.331667 0.6226570 25 0.53266 0.5990
## Mate6 -1.991667 0.6226570 25 -3.19866 0.0037
## Per2 -2.578333 0.6088981 75 -4.23443 0.0001
## Per3 0.390000 0.6088981 75 0.64050 0.5238
## Per4 -2.061667 0.6088981 75 -3.38590 0.0011
## AmbFR 2.178333 0.5919657 200 3.67983 0.0003
## AmbREF 0.415000 0.5919657 200 0.70105 0.4841
## Mate1.5:Per2 0.015000 0.8611119 75 0.01742 0.9861
## Mate3:Per2 1.623333 0.8611119 75 1.88516 0.0633
## Mate4.5:Per2 0.275000 0.8611119 75 0.31935 0.7503
## Mate6:Per2 2.310000 0.8611119 75 2.68258 0.0090
## Mate1.5:Per3 0.121667 0.8611119 75 0.14129 0.8880
## Mate3:Per3 0.458333 0.8611119 75 0.53226 0.5961
## Mate4.5:Per3 -0.523333 0.8611119 75 -0.60774 0.5452
## Mate6:Per3 2.671667 0.8611119 75 3.10258 0.0027
## Mate1.5:Per4 -0.823333 0.8611119 75 -0.95613 0.3421
## Mate3:Per4 0.498333 0.8611119 75 0.57871 0.5645
## Mate4.5:Per4 1.795000 0.8611119 75 2.08451 0.0405
## Mate6:Per4 3.161667 0.8611119 75 3.67161 0.0004
## Mate1.5:AmbFR -2.570000 0.8371660 200 -3.06988 0.0024
## Mate3:AmbFR -2.553333 0.8371660 200 -3.04997 0.0026
## Mate4.5:AmbFR -2.623333 0.8371660 200 -3.13359 0.0020
## Mate6:AmbFR 1.356667 0.8371660 200 1.62055 0.1067
## Mate1.5:AmbREF -0.463333 0.8371660 200 -0.55345 0.5806
## Mate3:AmbREF -0.500000 0.8371660 200 -0.59725 0.5510
## Mate4.5:AmbREF 0.875000 0.8371660 200 1.04519 0.2972
## Mate6:AmbREF -0.821667 0.8371660 200 -0.98149 0.3275
## Per2:AmbFR -1.043333 0.8371660 200 -1.24627 0.2141
## Per3:AmbFR -0.386667 0.8371660 200 -0.46188 0.6447
## Per4:AmbFR -0.015000 0.8371660 200 -0.01792 0.9857
## Per2:AmbREF 5.305000 0.8371660 200 6.33686 0.0000
## Per3:AmbREF -2.223333 0.8371660 200 -2.65579 0.0086
## Per4:AmbREF 3.136667 0.8371660 200 3.74677 0.0002
## Mate1.5:Per2:AmbFR 0.676667 1.1839315 200 0.57154 0.5683
## Mate3:Per2:AmbFR 0.986667 1.1839315 200 0.83338 0.4056
## Mate4.5:Per2:AmbFR 2.026667 1.1839315 200 1.71181 0.0885
## Mate6:Per2:AmbFR -1.478333 1.1839315 200 -1.24866 0.2132
## Mate1.5:Per3:AmbFR 0.953333 1.1839315 200 0.80523 0.4216
## Mate3:Per3:AmbFR 2.058333 1.1839315 200 1.73856 0.0837
## Mate4.5:Per3:AmbFR 2.170000 1.1839315 200 1.83288 0.0683
## Mate6:Per3:AmbFR -2.013333 1.1839315 200 -1.70055 0.0906
## Mate1.5:Per4:AmbFR 2.121667 1.1839315 200 1.79205 0.0746
## Mate3:Per4:AmbFR 1.845000 1.1839315 200 1.55837 0.1207
## Mate4.5:Per4:AmbFR 0.020000 1.1839315 200 0.01689 0.9865
## Mate6:Per4:AmbFR -2.781667 1.1839315 200 -2.34952 0.0198
## Mate1.5:Per2:AmbREF -0.318333 1.1839315 200 -0.26888 0.7883
## Mate3:Per2:AmbREF -0.515000 1.1839315 200 -0.43499 0.6640
## Mate4.5:Per2:AmbREF -2.045000 1.1839315 200 -1.72730 0.0857
## Mate6:Per2:AmbREF 0.405000 1.1839315 200 0.34208 0.7326
## Mate1.5:Per3:AmbREF 0.013333 1.1839315 200 0.01126 0.9910
## Mate3:Per3:AmbREF 0.981667 1.1839315 200 0.82916 0.4080
## Mate4.5:Per3:AmbREF 1.158333 1.1839315 200 0.97838 0.3291
## Mate6:Per3:AmbREF 1.151667 1.1839315 200 0.97275 0.3319
## Mate1.5:Per4:AmbREF -0.168333 1.1839315 200 -0.14218 0.8871
## Mate3:Per4:AmbREF -0.098333 1.1839315 200 -0.08306 0.9339
## Mate4.5:Per4:AmbREF -2.935000 1.1839315 200 -2.47903 0.0140
## Mate6:Per4:AmbREF -0.975000 1.1839315 200 -0.82353 0.4112
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -3.50051768 -0.55275335 -0.01985377 0.51530953 3.69723642
##
## Number of Observations: 360
## Number of Groups:
## parcela subparc %in% parcela
## 30 120# Médias ajustadas.
emm <- emmeans(mm1, specs = ~Mate + Per + Amb)
ggplot(data = as.data.frame(emm),
mapping = aes(x = Amb,
y = emmean)) +
facet_grid(facets = Mate ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Ambiente de armazenamento",
y = "altge") +
ggtitle(label = NULL,
subtitle = "Concentração de mate x período de postura")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Per,
y = emmean)) +
facet_grid(facets = Mate ~ Amb) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Período de postura",
y = "altge") +
ggtitle(label = NULL,
subtitle = "Concentração de mate x ambiente de armazenamento")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Mate,
y = emmean)) +
facet_grid(facets = Amb ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Concentração de mate",
y = "altge") +
ggtitle(label = NULL,
subtitle = "Ambiente de armazenamento x período de postura")
# Estudo da interação Mate:Per:Amb via Mate | Per + Amb.
emm <- emmeans(mm1, specs = ~Mate | Per + Amb)
grid <- multcomp::cld(emm)
grid## Per = 1, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 6 13.2 0.44 25 12.3 14.1 1
## 3 15.1 0.44 25 14.2 16.0 2
## 0 15.2 0.44 29 14.3 16.1 2
## 4.5 15.5 0.44 25 14.6 16.4 2
## 1.5 16.0 0.44 25 15.1 16.9 2
##
## Per = 1, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 3 14.8 0.44 25 13.9 15.7 1
## 4.5 15.1 0.44 25 14.2 16.0 12
## 1.5 15.6 0.44 25 14.7 16.5 123
## 6 16.7 0.44 25 15.8 17.6 23
## 0 17.4 0.44 29 16.5 18.3 3
##
## Per = 1, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 6 12.8 0.44 25 11.9 13.7 1
## 3 15.1 0.44 25 14.1 16.0 2
## 0 15.6 0.44 29 14.7 16.5 2
## 1.5 15.9 0.44 25 15.0 16.9 2
## 4.5 16.8 0.44 25 15.9 17.7 2
##
## Per = 2, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 0 12.6 0.44 29 11.7 13.5 1
## 6 12.9 0.44 25 12.0 13.8 1
## 4.5 13.2 0.44 25 12.3 14.1 1
## 1.5 13.4 0.44 25 12.5 14.3 1
## 3 14.2 0.44 25 13.3 15.1 1
##
## Per = 2, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 12.7 0.44 25 11.8 13.6 1
## 0 13.7 0.44 29 12.8 14.6 1
## 3 13.8 0.44 25 12.8 14.7 1
## 4.5 13.8 0.44 25 12.9 14.7 1
## 6 13.9 0.44 25 13.0 14.9 1
##
## Per = 2, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 4.5 17.8 0.44 25 16.9 18.7 1
## 6 18.2 0.44 25 17.3 19.1 1
## 0 18.3 0.44 29 17.4 19.2 1
## 1.5 18.4 0.44 25 17.5 19.3 1
## 3 18.9 0.44 25 18.0 19.8 1
##
## Per = 3, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 4.5 15.4 0.44 25 14.5 16.3 1
## 0 15.6 0.44 29 14.7 16.5 1
## 3 16.0 0.44 25 15.1 16.9 1
## 6 16.3 0.44 25 15.4 17.2 1
## 1.5 16.5 0.44 25 15.6 17.4 1
##
## Per = 3, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 16.7 0.44 25 15.8 17.6 1
## 4.5 16.7 0.44 25 15.8 17.6 1
## 3 17.3 0.44 25 16.4 18.2 1
## 0 17.4 0.44 29 16.5 18.3 1
## 6 17.4 0.44 25 16.5 18.3 1
##
## Per = 3, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 0 13.8 0.44 29 12.9 14.7 1
## 1.5 14.3 0.44 25 13.3 15.2 12
## 3 14.7 0.44 25 13.8 15.6 12
## 6 14.8 0.44 25 13.9 15.7 12
## 4.5 15.6 0.44 25 14.7 16.5 2
##
## Per = 4, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 13.1 0.44 25 12.2 14.0 1
## 0 13.1 0.44 29 12.2 14.0 1
## 3 13.6 0.44 25 12.7 14.5 12
## 6 14.3 0.44 25 13.4 15.2 12
## 4.5 15.3 0.44 25 14.3 16.2 2
##
## Per = 4, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 4.5 14.8 0.44 25 13.9 15.7 1
## 1.5 14.8 0.44 25 13.9 15.7 1
## 3 15.0 0.44 25 14.1 15.9 1
## 6 15.0 0.44 25 14.1 15.9 1
## 0 15.3 0.44 29 14.4 16.2 1
##
## Per = 4, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 16.0 0.44 25 15.1 16.9 1
## 6 16.1 0.44 25 15.1 17.0 1
## 3 16.5 0.44 25 15.6 17.4 1
## 0 16.7 0.44 29 15.8 17.6 1
## 4.5 16.7 0.44 25 15.8 17.7 1
##
## d.f. method: containment
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 5 estimates
## significance level used: alpha = 0.05grid <- grid %>%
group_by(Per, Amb) %>%
mutate(cld = numbers_to_letters(.group)) %>%
ungroup()
gg_means(grid,
mapping = list(x = "Mate"),
facets = Amb ~ Per,
labs = list(x = legends$Mate,
y = "altge",
subtitle = paste(legends$Amb, "x", legends$Per)))
#Análise da variável diamge
# Modelo de 1 extrato para avaliar os pressupostos.
m0 <- lm(terms(diamge ~
Mate + Mate:factor(rep) +
Per + Mate:Per + Mate:Per:factor(rep) +
Amb + Mate:Amb + Per:Amb + Mate:Per:Amb,
keep.order = TRUE),
data = tb)
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Apenas para verificação. Estatísticas F erradas pros estratos
# superiores.
anova(m0)## Analysis of Variance Table
##
## Response: diamge
## Df Sum Sq Mean Sq F value Pr(>F)
## Mate 4 37.5 9.4 3.3293 0.0114982 *
## Mate:factor(rep) 25 209.6 8.4 2.9752 1.171e-05 ***
## Per 3 460.5 153.5 54.4714 < 2.2e-16 ***
## Mate:Per 12 58.2 4.9 1.7223 0.0642844 .
## Mate:factor(rep):Per 75 273.6 3.6 1.2944 0.0810644 .
## Amb 2 7577.5 3788.7 1344.5511 < 2.2e-16 ***
## Mate:Amb 8 18.2 2.3 0.8069 0.5972109
## Per:Amb 6 202.1 33.7 11.9512 1.905e-11 ***
## Mate:Per:Amb 24 160.7 6.7 2.3766 0.0006076 ***
## Residuals 199 560.8 2.8
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Ajuste do modelo de efeitos aleatórios.
mm0 <- lme(diamge ~ Mate * Per * Amb,
random = ~1 | parcela/subparc,
method = "ML",
na.action = na.omit,
data = tb)
# Quadro de testes para o termos de efeito fixo.
anova(mm0)## numDF denDF F-value p-value
## (Intercept) 1 199 74132.11 <.0001
## Mate 4 25 1.19 0.3385
## Per 3 75 43.94 <.0001
## Amb 2 199 1347.59 <.0001
## Mate:Per 12 75 1.41 0.1819
## Mate:Amb 8 199 0.82 0.5833
## Per:Amb 6 199 12.07 <.0001
## Mate:Per:Amb 24 199 2.39 0.0006mm1 <- mm0
# Estimativas dos parâmetros.
summary(mm1)## Linear mixed-effects model fit by maximum likelihood
## Data: tb
## AIC BIC logLik
## 1501.408 1746.058 -687.7042
##
## Random effects:
## Formula: ~1 | parcela
## (Intercept)
## StdDev: 0.5592441
##
## Formula: ~1 | subparc %in% parcela
## (Intercept) Residual
## StdDev: 0.4339038 1.531334
##
## Fixed effects: diamge ~ Mate * Per * Amb
## Value Std.Error DF t-value p-value
## (Intercept) 48.44667 0.7546651 199 64.19625 0.0000
## Mate1.5 -1.96500 1.0672576 25 -1.84117 0.0775
## Mate3 -0.55833 1.0672576 25 -0.52315 0.6055
## Mate4.5 1.08000 1.0672576 25 1.01194 0.3213
## Mate6 -1.23833 1.0672576 25 -1.16029 0.2569
## Per2 -4.92167 1.0069098 75 -4.88789 0.0000
## Per3 0.10167 1.0069098 75 0.10097 0.9198
## Per4 -2.03833 1.0069098 75 -2.02435 0.0465
## AmbFR -7.48667 0.9687706 199 -7.72801 0.0000
## AmbREF -12.22000 0.9687706 199 -12.61393 0.0000
## Mate1.5:Per2 4.26833 1.4239855 75 2.99746 0.0037
## Mate3:Per2 2.40500 1.4239855 75 1.68892 0.0954
## Mate4.5:Per2 2.21333 1.4239855 75 1.55432 0.1243
## Mate6:Per2 4.87833 1.4239855 75 3.42583 0.0010
## Mate1.5:Per3 0.38500 1.4239855 75 0.27037 0.7876
## Mate3:Per3 -0.42167 1.4239855 75 -0.29612 0.7680
## Mate4.5:Per3 -2.52833 1.4239855 75 -1.77553 0.0799
## Mate6:Per3 0.76167 1.4239855 75 0.53488 0.5943
## Mate1.5:Per4 2.68333 1.4239855 75 1.88438 0.0634
## Mate3:Per4 1.50667 1.4239855 75 1.05806 0.2934
## Mate4.5:Per4 -1.43000 1.4239855 75 -1.00422 0.3185
## Mate6:Per4 3.27000 1.4239855 75 2.29637 0.0244
## Mate1.5:AmbFR 2.36667 1.3700485 199 1.72743 0.0856
## Mate3:AmbFR 0.08000 1.3700485 199 0.05839 0.9535
## Mate4.5:AmbFR -1.57667 1.3700485 199 -1.15081 0.2512
## Mate6:AmbFR 0.38833 1.3700485 199 0.28344 0.7771
## Mate1.5:AmbREF 3.47833 1.3700485 199 2.53884 0.0119
## Mate3:AmbREF 3.37000 1.3700485 199 2.45977 0.0148
## Mate4.5:AmbREF 2.81000 1.3700485 199 2.05102 0.0416
## Mate6:AmbREF 3.49167 1.3700485 199 2.54857 0.0116
## Per2:AmbFR -1.16333 1.3700485 199 -0.84912 0.3968
## Per3:AmbFR -4.88000 1.3700485 199 -3.56192 0.0005
## Per4:AmbFR -3.00333 1.3700485 199 -2.19214 0.0295
## Per2:AmbREF 5.34333 1.3700485 199 3.90011 0.0001
## Per3:AmbREF 0.82000 1.3700485 199 0.59852 0.5502
## Per4:AmbREF 2.59167 1.3700485 199 1.89166 0.0600
## Mate1.5:Per2:AmbFR -6.34500 1.9375412 199 -3.27477 0.0012
## Mate3:Per2:AmbFR -1.19833 1.9375412 199 -0.61848 0.5370
## Mate4.5:Per2:AmbFR 0.43333 1.9375412 199 0.22365 0.8233
## Mate6:Per2:AmbFR -2.38333 1.9375412 199 -1.23008 0.2201
## Mate1.5:Per3:AmbFR -0.31167 1.9375412 199 -0.16086 0.8724
## Mate3:Per3:AmbFR 1.21500 1.9375412 199 0.62708 0.5313
## Mate4.5:Per3:AmbFR 2.58833 1.9375412 199 1.33589 0.1831
## Mate6:Per3:AmbFR -0.48994 1.9643420 199 -0.24941 0.8033
## Mate1.5:Per4:AmbFR -1.72333 1.9375412 199 -0.88944 0.3748
## Mate3:Per4:AmbFR -0.61167 1.9375412 199 -0.31569 0.7526
## Mate4.5:Per4:AmbFR 2.21000 1.9375412 199 1.14062 0.2554
## Mate6:Per4:AmbFR -1.72833 1.9375412 199 -0.89202 0.3735
## Mate1.5:Per2:AmbREF -6.60000 1.9375412 199 -3.40638 0.0008
## Mate3:Per2:AmbREF -4.44833 1.9375412 199 -2.29587 0.0227
## Mate4.5:Per2:AmbREF -7.13500 1.9375412 199 -3.68250 0.0003
## Mate6:Per2:AmbREF -8.46500 1.9375412 199 -4.36894 0.0000
## Mate1.5:Per3:AmbREF -1.36333 1.9375412 199 -0.70364 0.4825
## Mate3:Per3:AmbREF -1.24000 1.9375412 199 -0.63999 0.5229
## Mate4.5:Per3:AmbREF -0.96167 1.9375412 199 -0.49633 0.6202
## Mate6:Per3:AmbREF -2.30167 1.9375412 199 -1.18793 0.2363
## Mate1.5:Per4:AmbREF -3.48500 1.9375412 199 -1.79867 0.0736
## Mate3:Per4:AmbREF -3.65333 1.9375412 199 -1.88555 0.0608
## Mate4.5:Per4:AmbREF -0.35000 1.9375412 199 -0.18064 0.8568
## Mate6:Per4:AmbREF -3.77333 1.9375412 199 -1.94749 0.0529
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.603676510 -0.609075414 -0.003396124 0.537618690 3.351788401
##
## Number of Observations: 359
## Number of Groups:
## parcela subparc %in% parcela
## 30 120# Médias ajustadas.
emm <- emmeans(mm1, specs = ~Mate + Per + Amb)
ggplot(data = as.data.frame(emm),
mapping = aes(x = Amb,
y = emmean)) +
facet_grid(facets = Mate ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Ambiente de armazenamento",
y = "diamge") +
ggtitle(label = NULL,
subtitle = "Concentração de mate x período de postura")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Per,
y = emmean)) +
facet_grid(facets = Mate ~ Amb) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Período de postura",
y = "diamge") +
ggtitle(label = NULL,
subtitle = "Concentração de mate x ambiente de armazenamento")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Mate,
y = emmean)) +
facet_grid(facets = Amb ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Concentração de mate",
y = "diamge") +
ggtitle(label = NULL,
subtitle = "Ambiente de armazenamento x período de postura")
# Estudo da interação Mate:Per:Amb via Mate | Per + Amb.
emm <- emmeans(mm1, specs = ~Mate | Per + Amb)
grid <- multcomp::cld(emm)
grid## Per = 1, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 46.5 0.755 25 44.9 48.0 1
## 6 47.2 0.755 25 45.7 48.8 1
## 3 47.9 0.755 25 46.3 49.4 1
## 0 48.4 0.755 29 46.9 50.0 1
## 4.5 49.5 0.755 25 48.0 51.1 1
##
## Per = 1, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 6 40.1 0.755 25 38.6 41.7 1
## 4.5 40.5 0.755 25 38.9 42.0 1
## 3 40.5 0.755 25 38.9 42.0 1
## 0 41.0 0.755 29 39.4 42.5 1
## 1.5 41.4 0.755 25 39.8 42.9 1
##
## Per = 1, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 0 36.2 0.755 29 34.7 37.8 1
## 1.5 37.7 0.755 25 36.2 39.3 12
## 6 38.5 0.755 25 36.9 40.0 12
## 3 39.0 0.755 25 37.5 40.6 12
## 4.5 40.1 0.755 25 38.6 41.7 2
##
## Per = 2, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 0 43.5 0.755 29 42.0 45.1 1
## 3 45.4 0.755 25 43.8 46.9 12
## 1.5 45.8 0.755 25 44.3 47.4 12
## 4.5 46.8 0.755 25 45.3 48.4 2
## 6 47.2 0.755 25 45.6 48.7 2
##
## Per = 2, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 33.2 0.755 25 31.6 34.8 1
## 0 34.9 0.755 29 33.3 36.4 12
## 3 35.6 0.755 25 34.0 37.2 12
## 6 36.5 0.755 25 35.0 38.1 2
## 4.5 37.0 0.755 25 35.5 38.6 2
##
## Per = 2, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 6 35.3 0.755 25 33.8 36.9 1
## 4.5 35.6 0.755 25 34.1 37.2 1
## 1.5 35.8 0.755 25 34.3 37.4 1
## 0 36.6 0.755 29 35.1 38.2 1
## 3 37.4 0.755 25 35.9 39.0 1
##
## Per = 3, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 1.5 47.0 0.755 25 45.4 48.5 1
## 4.5 47.1 0.755 25 45.5 48.7 1
## 3 47.6 0.755 25 46.0 49.1 1
## 6 48.1 0.755 25 46.5 49.6 1
## 0 48.5 0.755 29 47.0 50.1 1
##
## Per = 3, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 6 35.6 0.821 25 33.9 37.3 1
## 4.5 35.7 0.755 25 34.2 37.3 1
## 0 36.2 0.755 29 34.6 37.7 1
## 3 36.5 0.755 25 34.9 38.1 1
## 1.5 36.7 0.755 25 35.1 38.2 1
##
## Per = 3, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 0 37.1 0.755 29 35.6 38.7 1
## 4.5 37.5 0.755 25 36.0 39.1 1
## 1.5 37.7 0.755 25 36.1 39.2 1
## 6 37.9 0.755 25 36.3 39.4 1
## 3 38.3 0.755 25 36.7 39.9 1
##
## Per = 4, Amb = AM:
## Mate emmean SE df lower.CL upper.CL .group
## 4.5 46.1 0.755 25 44.5 47.6 1
## 0 46.4 0.755 29 44.9 48.0 1
## 1.5 47.1 0.755 25 45.6 48.7 1
## 3 47.4 0.755 25 45.8 48.9 1
## 6 48.4 0.755 25 46.9 50.0 1
##
## Per = 4, Amb = FR:
## Mate emmean SE df lower.CL upper.CL .group
## 0 35.9 0.755 29 34.4 37.5 1
## 4.5 36.2 0.755 25 34.6 37.8 1
## 3 36.3 0.755 25 34.8 37.9 1
## 6 36.6 0.755 25 35.1 38.2 1
## 1.5 37.3 0.755 25 35.7 38.8 1
##
## Per = 4, Amb = REF:
## Mate emmean SE df lower.CL upper.CL .group
## 0 36.8 0.755 29 35.2 38.3 1
## 3 37.4 0.755 25 35.9 39.0 1
## 1.5 37.5 0.755 25 35.9 39.0 1
## 6 38.5 0.755 25 37.0 40.1 1
## 4.5 38.9 0.755 25 37.3 40.4 1
##
## d.f. method: containment
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 5 estimates
## significance level used: alpha = 0.05grid <- grid %>%
group_by(Per, Amb) %>%
mutate(cld = numbers_to_letters(.group)) %>%
ungroup()
gg_means(grid,
mapping = list(x = "Mate"),
facets = Amb ~ Per,
labs = list(x = legends$Mate,
y = "diamge",
subtitle = paste(legends$Amb, "x", legends$Per)))
#Análise da variável pesoge
# Modelo de 1 extrato para avaliar os pressupostos.
m0 <- lm(terms(pesoge ~
Mate + Mate:factor(rep) +
Per + Mate:Per + Mate:Per:factor(rep) +
Amb + Mate:Amb + Per:Amb + Mate:Per:Amb,
keep.order = TRUE),
data = tb)
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Apenas para verificação. Estatísticas F erradas pros estratos
# superiores.
anova(m0)## Analysis of Variance Table
##
## Response: pesoge
## Df Sum Sq Mean Sq F value Pr(>F)
## Mate 4 24.596 6.1489 5.9698 0.0001474 ***
## Mate:factor(rep) 25 69.100 2.7640 2.6835 7.451e-05 ***
## Per 3 21.338 7.1128 6.9056 0.0001902 ***
## Mate:Per 12 13.196 1.0996 1.0676 0.3893575
## Mate:factor(rep):Per 75 71.381 0.9518 0.9240 0.6481299
## Amb 2 9.443 4.7217 4.5842 0.0113085 *
## Mate:Amb 8 6.792 0.8490 0.8243 0.5820649
## Per:Amb 6 59.726 9.9543 9.6643 2.446e-09 ***
## Mate:Per:Amb 24 29.773 1.2406 1.2044 0.2414432
## Residuals 200 206.001 1.0300
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Ajuste do modelo de efeitos aleatórios.
mm0 <- lme(pesoge ~ Mate * Per * Amb,
random = ~1 | parcela/subparc,
method = "ML",
na.action = na.omit,
data = tb)
# Quadro de testes para o termos de efeito fixo.
anova(mm0)## numDF denDF F-value p-value
## (Intercept) 1 200 34411.10 <.0001
## Mate 4 25 2.22 0.0952
## Per 3 75 7.05 0.0003
## Amb 2 200 4.68 0.0103
## Mate:Per 12 75 1.09 0.3808
## Mate:Amb 8 200 0.84 0.5670
## Per:Amb 6 200 9.87 <.0001
## Mate:Per:Amb 24 200 1.23 0.2196# Atualiza mantendo os termos relevantes.
mm1 <- update(mm0, . ~ Per * Amb)
# Teste da razão de verossimilhanças (é anticonservador para teste de
# termos de efeitos fixos).
anova(mm1, mm0)## Model df AIC BIC logLik Test L.Ratio p-value
## mm1 1 15 1082.933 1141.224 -526.4663
## mm0 2 63 1115.347 1360.172 -494.6736 1 vs 2 63.5855 0.0653# Estimativas dos parâmetros.
summary(mm1)## Linear mixed-effects model fit by maximum likelihood
## Data: tb
## AIC BIC logLik
## 1082.933 1141.224 -526.4663
##
## Random effects:
## Formula: ~1 | parcela
## (Intercept)
## StdDev: 0.4214914
##
## Formula: ~1 | subparc %in% parcela
## (Intercept) Residual
## StdDev: 0.0001005047 0.9956627
##
## Fixed effects: pesoge ~ Per + Amb + Per:Amb
## Value Std.Error DF t-value p-value
## (Intercept) 16.543000 0.2007743 232 82.39600 0.0000
## Per2 0.235667 0.2614738 87 0.90130 0.3699
## Per3 0.051000 0.2614738 87 0.19505 0.8458
## Per4 -0.533000 0.2614738 87 -2.03844 0.0445
## AmbFR 0.491667 0.2614738 232 1.88037 0.0613
## AmbREF -0.378667 0.2614738 232 -1.44820 0.1489
## Per2:AmbFR -2.203000 0.3697798 232 -5.95760 0.0000
## Per3:AmbFR -1.256000 0.3697798 232 -3.39662 0.0008
## Per4:AmbFR 0.026667 0.3697798 232 0.07211 0.9426
## Per2:AmbREF -0.564000 0.3697798 232 -1.52523 0.1286
## Per3:AmbREF 0.203667 0.3697798 232 0.55078 0.5823
## Per4:AmbREF 0.615333 0.3697798 232 1.66405 0.0975
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -3.91543657 -0.62607997 0.06669838 0.61290354 3.11232986
##
## Number of Observations: 360
## Number of Groups:
## parcela subparc %in% parcela
## 30 120# Médias ajustadas.
emm <- emmeans(mm1, specs = ~Per + Amb)
ggplot(data = as.data.frame(emm),
mapping = aes(x = Amb,
y = emmean)) +
facet_grid(facets = ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Ambiente de armazenamento",
y = "pesoge") +
ggtitle(label = NULL,
subtitle = "Período de postura")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Per,
y = emmean)) +
facet_grid(facets = ~Amb) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Período de postura",
y = "pesoge") +
ggtitle(label = NULL,
subtitle = "Ambiente de armazenamento")
# Estudo da interação Per:Amb via Amb | Per.
emm <- emmeans(mm1, specs = ~Amb | Per)
grid <- multcomp::cld(emm)
grid## Per = 1:
## Amb emmean SE df lower.CL upper.CL .group
## REF 16.2 0.201 29 15.8 16.6 1
## AM 16.5 0.201 29 16.1 17.0 12
## FR 17.0 0.201 29 16.6 17.4 2
##
## Per = 2:
## Amb emmean SE df lower.CL upper.CL .group
## FR 15.1 0.201 29 14.7 15.5 1
## REF 15.8 0.201 29 15.4 16.2 2
## AM 16.8 0.201 29 16.4 17.2 3
##
## Per = 3:
## Amb emmean SE df lower.CL upper.CL .group
## FR 15.8 0.201 29 15.4 16.2 1
## REF 16.4 0.201 29 16.0 16.8 12
## AM 16.6 0.201 29 16.2 17.0 2
##
## Per = 4:
## Amb emmean SE df lower.CL upper.CL .group
## AM 16.0 0.201 29 15.6 16.4 1
## REF 16.2 0.201 29 15.8 16.7 1
## FR 16.5 0.201 29 16.1 16.9 1
##
## d.f. method: containment
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 3 estimates
## significance level used: alpha = 0.05grid <- grid %>%
group_by(Per) %>%
mutate(cld = numbers_to_letters(.group)) %>%
ungroup()
gg_means(grid,
mapping = list(x = "Amb"),
facets = ~Per,
labs = list(x = legends$Amb,
y = "pesoge",
subtitle = legends$Per))
#Análise da variável pH
# Modelo de 1 extrato para avaliar os pressupostos.
m0 <- lm(terms(ph ~
Mate + Mate:factor(rep) +
Per + Mate:Per + Mate:Per:factor(rep) +
Amb + Mate:Amb + Per:Amb + Mate:Per:Amb,
keep.order = TRUE),
data = tb)
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Apenas para verificação. Estatísticas F erradas pros estratos
# superiores.
anova(m0)## Analysis of Variance Table
##
## Response: ph
## Df Sum Sq Mean Sq F value Pr(>F)
## Mate 4 0.2707 0.06768 2.6076 0.03691 *
## Mate:factor(rep) 25 0.7875 0.03150 1.2136 0.23049
## Per 3 0.0552 0.01841 0.7092 0.54755
## Mate:Per 12 0.5184 0.04320 1.6643 0.07702 .
## Mate:factor(rep):Per 75 1.8814 0.02509 0.9665 0.55893
## Amb 2 1.8169 0.90844 35.0000 9.260e-14 ***
## Mate:Amb 8 0.3323 0.04153 1.6002 0.12664
## Per:Amb 6 1.4484 0.24141 9.3008 5.399e-09 ***
## Mate:Per:Amb 24 0.7713 0.03214 1.2381 0.21287
## Residuals 200 5.1911 0.02596
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1# Ajuste do modelo de efeitos aleatórios.
mm0 <- lme(ph ~ Mate * Per * Amb,
random = ~1 | parcela/subparc,
method = "ML",
na.action = na.omit,
data = tb)
# Quadro de testes para o termos de efeito fixo.
anova(mm0)## numDF denDF F-value p-value
## (Intercept) 1 200 443180.5 <.0001
## Mate 4 25 2.1 0.1045
## Per 3 75 0.7 0.5457
## Amb 2 200 35.3 <.0001
## Mate:Per 12 75 1.7 0.0885
## Mate:Amb 8 200 1.6 0.1224
## Per:Amb 6 200 9.4 <.0001
## Mate:Per:Amb 24 200 1.2 0.2038# Atualiza mantendo os termos relevantes.
mm1 <- update(mm0, . ~ Per * Amb)
# Teste da razão de verossimilhanças (é anticonservador para teste de
# termos de efeitos fixos).
anova(mm1, mm0)## Model df AIC BIC logLik Test L.Ratio p-value
## mm1 1 15 -248.7262 -190.43461 139.3631
## mm0 2 63 -229.7169 15.10767 177.8584 1 vs 2 76.99071 0.005# Estimativas dos parâmetros.
summary(mm1)## Linear mixed-effects model fit by maximum likelihood
## Data: tb
## AIC BIC logLik
## -248.7262 -190.4346 139.3631
##
## Random effects:
## Formula: ~1 | parcela
## (Intercept)
## StdDev: 0.026786
##
## Formula: ~1 | subparc %in% parcela
## (Intercept) Residual
## StdDev: 0.01206122 0.1619499
##
## Fixed effects: ph ~ Per + Amb + Per:Amb
## Value Std.Error DF t-value p-value
## (Intercept) 6.330000 0.03056409 232 207.10580 0.0000
## Per2 0.010000 0.04264792 87 0.23448 0.8152
## Per3 0.080000 0.04264792 87 1.87582 0.0640
## Per4 -0.136667 0.04264792 87 -3.20453 0.0019
## AmbFR -0.183333 0.04253013 232 -4.31067 0.0000
## AmbREF -0.100000 0.04253013 232 -2.35127 0.0195
## Per2:AmbFR 0.026667 0.06014669 232 0.44336 0.6579
## Per3:AmbFR -0.206667 0.06014669 232 -3.43604 0.0007
## Per4:AmbFR 0.220000 0.06014669 232 3.65772 0.0003
## Per2:AmbREF -0.040000 0.06014669 232 -0.66504 0.5067
## Per3:AmbREF -0.070000 0.06014669 232 -1.16382 0.2457
## Per4:AmbREF 0.110000 0.06014669 232 1.82886 0.0687
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.4417630 -0.4940944 -0.1146826 0.2385605 5.9209929
##
## Number of Observations: 360
## Number of Groups:
## parcela subparc %in% parcela
## 30 120# Médias ajustadas.
emm <- emmeans(mm1, specs = ~Per + Amb)
ggplot(data = as.data.frame(emm),
mapping = aes(x = Amb,
y = emmean)) +
facet_grid(facets = ~ Per) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Ambiente de armazenamento",
y = "ph") +
ggtitle(label = NULL,
subtitle = "Período de postura")
ggplot(data = as.data.frame(emm),
mapping = aes(x = Per,
y = emmean)) +
facet_grid(facets = ~Amb) +
geom_point() +
geom_line(mapping = aes(group = 1)) +
geom_errorbar(mapping = aes(ymin = lower.CL,
ymax = upper.CL),
width = 0.1) +
labs(x = "Período de postura",
y = "ph") +
ggtitle(label = NULL,
subtitle = "Ambiente de armazenamento")
# Estudo da interação Per:Amb via Amb | Per.
emm <- emmeans(mm1, specs = ~Amb | Per)
grid <- multcomp::cld(emm)
grid## Per = 1:
## Amb emmean SE df lower.CL upper.CL .group
## FR 6.15 0.0306 29 6.08 6.21 1
## REF 6.23 0.0306 29 6.17 6.29 12
## AM 6.33 0.0306 29 6.27 6.39 2
##
## Per = 2:
## Amb emmean SE df lower.CL upper.CL .group
## FR 6.18 0.0306 29 6.12 6.25 1
## REF 6.20 0.0306 29 6.14 6.26 1
## AM 6.34 0.0306 29 6.28 6.40 2
##
## Per = 3:
## Amb emmean SE df lower.CL upper.CL .group
## FR 6.02 0.0306 29 5.96 6.08 1
## REF 6.24 0.0306 29 6.18 6.30 2
## AM 6.41 0.0306 29 6.35 6.47 3
##
## Per = 4:
## Amb emmean SE df lower.CL upper.CL .group
## AM 6.19 0.0306 29 6.13 6.26 1
## REF 6.20 0.0306 29 6.14 6.27 1
## FR 6.23 0.0306 29 6.17 6.29 1
##
## d.f. method: containment
## Confidence level used: 0.95
## P value adjustment: tukey method for comparing a family of 3 estimates
## significance level used: alpha = 0.05grid <- grid %>%
group_by(Per) %>%
mutate(cld = numbers_to_letters(.group)) %>%
ungroup()
gg_means(grid,
mapping = list(x = "Amb"),
facets = ~Per,
labs = list(x = legends$Amb,
y = "ph",
subtitle = legends$Per))
#