Definições da sessão
#-----------------------------------------------------------------------
# Prof. Dr. Walmes M. Zeviani
# leg.ufpr.br/~walmes · github.com/walmes
# walmes@ufpr.br · @walmeszeviani
# Laboratory of Statistics and Geoinformation (LEG)
# Department of Statistics · Federal University of Paraná
# 2021-fev-25 · Curitiba/PR/Brazil
#-----------------------------------------------------------------------
#-----------------------------------------------------------------------
# Pacotes.
rm(list = objects())
library(emmeans)
library(tidyverse)
library(readxl)
Importação e preparo
#-----------------------------------------------------------------------
# Importação.
# Importação de tipagem das variáveis.
tb <- read_excel("teca raiz.xlsx",
skip = 1) %>%
mutate(gen = factor(gen),
block = factor(block),
enraz = factor(enraz, levels = c("I", "II", "III", "IV")))
str(tb)
## tibble [111 × 13] (S3: tbl_df/tbl/data.frame)
## $ index : num [1:111] 1 2 3 4 5 6 7 8 9 10 ...
## $ gen : Factor w/ 37 levels "1101","1102",..: 12 7 22 37 6 35 32 8 34 24 ...
## $ block : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ altini : num [1:111] 6.5 6.5 11 9.5 4 7.5 7.5 10.5 17 9 ...
## $ alt : num [1:111] 85 109 105 52 59 42 149 59 101 155 ...
## $ diamini: num [1:111] 6.95 3.9 4.58 7.33 4.97 ...
## $ diam : num [1:111] 27.4 20.8 25.2 11.9 17.1 ...
## $ mffol : num [1:111] 456 231 272 161 172 ...
## $ mfcau : num [1:111] 434.8 199.8 276.1 56.5 126.2 ...
## $ msfol : num [1:111] 193 96 120.8 61 55.5 ...
## $ mscau : num [1:111] 206 73.6 131.6 26.3 54.9 ...
## $ msraiz : num [1:111] NA 87.3 104.8 37.9 82.3 ...
## $ enraz : Factor w/ 4 levels "I","II","III",..: 4 2 3 1 2 1 2 3 3 2 ...
# Cruzamento com nÃveis da variável resposta qualitativa.
xtabs(~gen + enraz, data = tb) %>%
addmargins()
## enraz
## gen I II III IV Sum
## 1101 0 1 1 1 3
## 1102 1 2 0 0 3
## 1106 3 0 0 0 3
## 1114 0 2 1 0 3
## 1328 1 1 1 0 3
## 1507 0 2 0 1 3
## 2011 2 1 0 0 3
## 3105 0 2 1 0 3
## 3408 3 0 0 0 3
## 4008 3 0 0 0 3
## 4010 0 1 1 0 2
## 4013 0 1 1 1 3
## 4015 0 0 1 2 3
## 4016 2 0 1 0 3
## 4017 0 2 0 1 3
## 4018 2 1 0 0 3
## BT 86 0 2 1 0 3
## BT 93 0 3 0 0 3
## BT21 1 0 1 1 3
## BT27 3 0 0 0 3
## BT61 0 1 0 1 2
## Bt62 0 1 1 0 2
## BT63 0 0 0 2 2
## BT64 0 3 0 0 3
## BT68 0 0 3 0 3
## BT70 3 0 0 0 3
## BT72 0 2 1 0 3
## BT82 2 1 0 0 3
## BT83 2 0 1 0 3
## BT85 0 1 1 1 3
## BT85 PLUS 0 3 0 0 3
## BT90 0 2 0 1 3
## BT96 0 3 0 0 3
## BT98 0 0 2 1 3
## FV03 2 1 0 0 3
## FV05 0 1 1 1 3
## FV10 3 0 0 0 3
## Sum 33 40 20 14 107
tb_means <- list()
Altura de planta (cm)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(alt ~ altini + block + gen, data = tb)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: alt
## Df Sum Sq Mean Sq F value Pr(>F)
## altini 1 160 159.64 0.2107 0.647633
## block 2 3246 1622.83 2.1417 0.124979
## gen 36 60314 1675.39 2.2111 0.002205 **
## Residuals 71 53798 757.72
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## BT82 54.3 15.9 71 1.37 107 a
## 3105 55.6 15.9 71 2.69 109 a
## 1114 59.0 15.9 71 6.04 112 ab
## FV03 68.1 16.5 71 13.11 123 ab
## 3408 71.1 16.7 71 15.42 127 ab
## 1328 72.9 16.1 71 19.31 126 ab
## FV10 73.9 15.9 71 21.01 127 ab
## BT61 75.7 16.5 71 20.98 131 ab
## BT21 78.3 15.9 71 25.26 131 ab
## 4018 78.9 15.9 71 26.04 132 ab
## 4015 79.0 16.0 71 25.60 132 ab
## 1507 80.7 17.4 71 22.87 139 ab
## 4017 81.3 16.0 71 28.18 134 ab
## 1102 81.6 15.9 71 28.66 135 ab
## 4016 83.9 16.0 71 30.72 137 ab
## FV05 85.0 16.3 71 30.76 139 ab
## BT68 87.5 16.3 71 33.39 142 ab
## 4008 90.5 16.4 71 36.04 145 ab
## 4010 91.0 16.1 71 37.50 145 ab
## BT72 91.6 15.9 71 38.64 145 ab
## BT85 94.0 16.0 71 40.60 147 ab
## BT85 PLUS 96.0 16.0 71 42.60 149 ab
## BT70 99.6 16.0 71 46.30 153 ab
## BT98 100.5 17.4 71 42.73 158 ab
## 4013 104.9 15.9 71 52.01 158 ab
## BT83 108.3 15.9 71 55.34 161 ab
## BT96 109.3 15.9 71 56.38 162 ab
## BT27 111.9 15.9 71 58.97 165 ab
## 1101 112.7 16.4 71 58.16 167 ab
## Bt62 122.7 16.0 71 69.44 176 ab
## 1106 124.4 16.3 71 70.00 179 ab
## BT63 125.1 16.4 71 70.49 180 ab
## BT90 126.5 16.5 71 71.65 181 ab
## BT 93 133.6 16.8 71 77.88 189 ab
## 2011 134.1 16.5 71 79.21 189 ab
## BT64 141.6 15.9 71 88.70 195 b
## BT 86 142.7 16.2 71 88.79 197 b
##
## Results are averaged over the levels of: block
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$alt <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "Altura de plantas (cm)")
#-----------------------------------------------------------------------
Diâmetro do caule (mm)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(diam ~ diamini + block + gen, data = tb)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: diam
## Df Sum Sq Mean Sq F value Pr(>F)
## diamini 1 72.91 72.908 4.3637 0.040298 *
## block 2 170.30 85.151 5.0965 0.008546 **
## gen 36 896.93 24.915 1.4912 0.076076 .
## Residuals 71 1186.26 16.708
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## FV10 13.3 2.41 71 5.32 21.4 a
## 3105 13.6 2.46 71 5.38 21.8 a
## BT82 13.6 2.48 71 5.38 21.9 a
## 4016 15.6 2.43 71 7.49 23.7 a
## 4015 15.8 2.54 71 7.37 24.3 a
## BT70 16.9 2.36 71 9.05 24.8 a
## 1101 17.6 2.36 71 9.71 25.4 a
## 4008 17.9 2.37 71 9.96 25.8 a
## BT72 18.0 2.36 71 10.08 25.8 a
## 4017 18.0 2.36 71 10.12 25.8 a
## 4010 18.7 2.43 71 10.63 26.8 a
## BT96 19.2 2.41 71 11.14 27.2 a
## 3408 19.3 2.37 71 11.37 27.1 a
## 4018 19.4 2.36 71 11.53 27.2 a
## 1114 20.0 2.38 71 12.13 28.0 a
## BT21 20.4 2.37 71 12.48 28.2 a
## 1106 20.5 2.37 71 12.55 28.4 a
## BT63 20.5 2.37 71 12.58 28.3 a
## BT98 20.8 2.36 71 12.91 28.6 a
## BT61 20.9 2.36 71 13.04 28.8 a
## FV05 21.3 2.36 71 13.44 29.2 a
## BT85 PLUS 21.6 2.44 71 13.43 29.7 a
## 1102 21.8 2.38 71 13.84 29.7 a
## 1507 21.8 2.38 71 13.87 29.7 a
## BT85 21.8 2.53 71 13.41 30.2 a
## BT83 21.9 2.36 71 14.04 29.8 a
## BT68 22.0 2.36 71 14.14 29.9 a
## FV03 22.0 2.38 71 14.10 29.9 a
## BT 93 22.1 2.38 71 14.15 30.0 a
## BT27 22.3 2.46 71 14.08 30.5 a
## BT90 22.6 2.52 71 14.18 31.0 a
## BT 86 22.8 2.36 71 14.98 30.7 a
## 1328 22.9 2.37 71 15.02 30.8 a
## BT64 23.9 2.43 71 15.78 32.0 a
## 4013 24.0 2.45 71 15.84 32.2 a
## 2011 24.5 2.42 71 16.43 32.5 a
## Bt62 25.0 2.40 71 17.03 33.0 a
##
## Results are averaged over the levels of: block
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$diam <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "Diâmetro do caule (mm)")
#-----------------------------------------------------------------------
Massa fresca de folhas (g)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(mffol ~ block + gen, data = tb)
MASS::boxcox(m0)
m0 <- update(m0, formula = log(.) ~ .)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: log(mffol)
## Df Sum Sq Mean Sq F value Pr(>F)
## block 2 0.3013 0.15066 1.3401 0.2683
## gen 36 6.1748 0.17152 1.5256 0.0645 .
## Residuals 72 8.0948 0.11243
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## 4016 4.79 0.194 72 4.14 5.43 a
## BT82 4.81 0.194 72 4.16 5.45 a
## 1328 4.89 0.194 72 4.25 5.53 a
## 4010 4.93 0.194 72 4.29 5.58 a
## BT21 5.00 0.194 72 4.35 5.64 a
## 1102 5.03 0.194 72 4.39 5.67 a
## 4015 5.06 0.194 72 4.41 5.70 a
## 3105 5.07 0.194 72 4.42 5.71 a
## 4018 5.08 0.194 72 4.44 5.72 a
## 1114 5.11 0.194 72 4.47 5.76 a
## BT96 5.12 0.194 72 4.48 5.77 a
## BT70 5.13 0.194 72 4.49 5.77 a
## BT68 5.15 0.194 72 4.51 5.80 a
## FV10 5.18 0.194 72 4.54 5.83 a
## BT72 5.18 0.194 72 4.54 5.83 a
## 4017 5.21 0.194 72 4.57 5.86 a
## 1106 5.25 0.194 72 4.61 5.89 a
## FV03 5.25 0.194 72 4.61 5.90 a
## 3408 5.27 0.194 72 4.62 5.91 a
## BT 93 5.28 0.194 72 4.64 5.92 a
## 1507 5.29 0.194 72 4.64 5.93 a
## BT98 5.29 0.194 72 4.65 5.94 a
## 1101 5.31 0.194 72 4.67 5.96 a
## BT85 5.38 0.194 72 4.74 6.03 a
## BT 86 5.40 0.194 72 4.75 6.04 a
## BT64 5.40 0.194 72 4.76 6.04 a
## 4008 5.46 0.194 72 4.82 6.11 a
## Bt62 5.47 0.194 72 4.82 6.11 a
## BT27 5.50 0.194 72 4.86 6.15 a
## BT85 PLUS 5.50 0.194 72 4.86 6.15 a
## BT63 5.51 0.194 72 4.86 6.15 a
## 2011 5.51 0.194 72 4.86 6.15 a
## BT61 5.54 0.194 72 4.90 6.19 a
## BT83 5.59 0.194 72 4.94 6.23 a
## BT90 5.66 0.194 72 5.01 6.30 a
## 4013 5.66 0.194 72 5.01 6.30 a
## FV05 5.69 0.194 72 5.05 6.34 a
##
## Results are averaged over the levels of: block
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$log_mffol <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "log da massa fresca de folhas (g)")
#-----------------------------------------------------------------------
Massa fresca de caules (g)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(mfcau ~ block + gen, data = tb)
MASS::boxcox(m0)
m0 <- update(m0, formula = log(.) ~ .)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: log(mfcau)
## Df Sum Sq Mean Sq F value Pr(>F)
## block 2 0.0574 0.028707 0.1641 0.84895
## gen 36 11.2631 0.312863 1.7888 0.01832 *
## Residuals 72 12.5932 0.174906
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## BT82 4.35 0.241 72 3.55 5.15 a
## FV10 4.36 0.241 72 3.55 5.16 a
## 3105 4.41 0.241 72 3.61 5.21 a
## 4016 4.79 0.241 72 3.99 5.59 a
## 4017 4.86 0.241 72 4.05 5.66 a
## 4015 4.89 0.241 72 4.09 5.70 a
## 3408 4.90 0.241 72 4.09 5.70 a
## 1114 4.92 0.241 72 4.11 5.72 a
## 4008 4.95 0.241 72 4.14 5.75 a
## 4018 4.95 0.241 72 4.14 5.75 a
## 1101 4.96 0.241 72 4.16 5.77 a
## BT21 4.98 0.241 72 4.18 5.78 a
## 4010 5.02 0.241 72 4.21 5.82 a
## BT70 5.09 0.241 72 4.29 5.90 a
## BT98 5.11 0.241 72 4.31 5.91 a
## BT68 5.12 0.241 72 4.31 5.92 a
## BT85 PLUS 5.16 0.241 72 4.36 5.97 a
## BT96 5.20 0.241 72 4.40 6.01 a
## 1102 5.21 0.241 72 4.41 6.01 a
## 1106 5.21 0.241 72 4.41 6.02 a
## 1328 5.22 0.241 72 4.41 6.02 a
## FV03 5.23 0.241 72 4.43 6.04 a
## 1507 5.26 0.241 72 4.46 6.06 a
## BT 93 5.26 0.241 72 4.46 6.07 a
## BT27 5.26 0.241 72 4.46 6.07 a
## BT72 5.29 0.241 72 4.49 6.09 a
## BT61 5.37 0.241 72 4.57 6.17 a
## FV05 5.40 0.241 72 4.60 6.20 a
## BT90 5.40 0.241 72 4.60 6.21 a
## BT83 5.40 0.241 72 4.60 6.21 a
## BT64 5.44 0.241 72 4.64 6.24 a
## BT85 5.46 0.241 72 4.65 6.26 a
## 4013 5.53 0.241 72 4.73 6.33 a
## Bt62 5.55 0.241 72 4.75 6.35 a
## BT63 5.56 0.241 72 4.75 6.36 a
## 2011 5.59 0.241 72 4.78 6.39 a
## BT 86 5.61 0.241 72 4.80 6.41 a
##
## Results are averaged over the levels of: block
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$log_mfcau <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "log da massa fresca de caules (g)")
#-----------------------------------------------------------------------
Massa seca de folhas (g)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(msfol ~ block + gen, data = tb)
MASS::boxcox(m0)
m0 <- update(m0, formula = log(.) ~ .)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: log(msfol)
## Df Sum Sq Mean Sq F value Pr(>F)
## block 2 0.4683 0.23417 1.8659 0.16216
## gen 36 8.0170 0.22270 1.7745 0.01966 *
## Residuals 72 9.0357 0.12549
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## 4016 3.78 0.205 72 3.10 4.46 a
## BT82 3.80 0.205 72 3.12 4.48 a
## 3105 3.99 0.205 72 3.31 4.67 a
## 1102 4.02 0.205 72 3.34 4.70 a
## 4015 4.16 0.205 72 3.48 4.84 a
## 4018 4.16 0.205 72 3.48 4.84 a
## FV10 4.18 0.205 72 3.50 4.86 a
## 1328 4.18 0.205 72 3.50 4.86 a
## BT21 4.19 0.205 72 3.51 4.87 a
## 4017 4.24 0.205 72 3.56 4.92 a
## 1114 4.24 0.205 72 3.56 4.92 a
## 3408 4.29 0.205 72 3.61 4.97 a
## BT96 4.29 0.205 72 3.61 4.98 a
## BT68 4.30 0.205 72 3.62 4.98 a
## 1101 4.32 0.205 72 3.64 5.00 a
## 4010 4.34 0.205 72 3.66 5.02 a
## BT85 4.37 0.205 72 3.69 5.05 a
## BT72 4.39 0.205 72 3.71 5.07 a
## 1106 4.40 0.205 72 3.72 5.08 a
## BT70 4.40 0.205 72 3.72 5.08 a
## BT98 4.42 0.205 72 3.74 5.10 a
## 1507 4.42 0.205 72 3.74 5.10 a
## FV03 4.43 0.205 72 3.74 5.11 a
## 4008 4.47 0.205 72 3.79 5.15 a
## BT 86 4.53 0.205 72 3.85 5.21 a
## BT61 4.58 0.205 72 3.90 5.26 a
## BT90 4.61 0.205 72 3.93 5.29 a
## BT 93 4.65 0.205 72 3.97 5.33 a
## BT63 4.65 0.205 72 3.97 5.33 a
## BT83 4.67 0.205 72 3.99 5.35 a
## BT85 PLUS 4.69 0.205 72 4.01 5.37 a
## Bt62 4.69 0.205 72 4.01 5.37 a
## BT27 4.74 0.205 72 4.06 5.42 a
## BT64 4.74 0.205 72 4.06 5.42 a
## 2011 4.81 0.205 72 4.13 5.49 a
## 4013 4.84 0.205 72 4.16 5.52 a
## FV05 4.84 0.205 72 4.16 5.52 a
##
## Results are averaged over the levels of: block
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$log_msfol <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "log da seca de folhas (g)")
#-----------------------------------------------------------------------
Massa seca de caules (g)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(mscau ~ block + gen, data = tb)
MASS::boxcox(m0)
m0 <- update(m0, formula = log(.) ~ .)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: log(mscau)
## Df Sum Sq Mean Sq F value Pr(>F)
## block 2 0.1292 0.06462 0.3698 0.692179
## gen 36 12.4827 0.34674 1.9843 0.006857 **
## Residuals 72 12.5814 0.17474
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## FV10 3.48 0.241 72 2.68 4.28 a
## BT82 3.50 0.241 72 2.69 4.30 ab
## 3105 3.50 0.241 72 2.70 4.31 ab
## 4016 3.91 0.241 72 3.11 4.71 ab
## 1114 3.93 0.241 72 3.13 4.73 ab
## 4017 4.01 0.241 72 3.21 4.81 ab
## 4015 4.02 0.241 72 3.22 4.83 ab
## 1106 4.04 0.241 72 3.23 4.84 ab
## 4018 4.07 0.241 72 3.27 4.87 ab
## 3408 4.09 0.241 72 3.29 4.89 ab
## 4008 4.12 0.241 72 3.32 4.93 ab
## BT21 4.15 0.241 72 3.35 4.95 ab
## BT70 4.20 0.241 72 3.39 5.00 ab
## 4010 4.21 0.241 72 3.40 5.01 ab
## 1102 4.21 0.241 72 3.40 5.01 ab
## 1328 4.21 0.241 72 3.41 5.01 ab
## 1101 4.23 0.241 72 3.43 5.03 ab
## BT98 4.26 0.241 72 3.46 5.06 ab
## BT85 PLUS 4.28 0.241 72 3.47 5.08 ab
## FV03 4.29 0.241 72 3.48 5.09 ab
## BT72 4.31 0.241 72 3.51 5.12 ab
## BT68 4.36 0.241 72 3.56 5.16 ab
## BT96 4.39 0.241 72 3.59 5.20 ab
## BT27 4.40 0.241 72 3.60 5.21 ab
## 1507 4.45 0.241 72 3.64 5.25 ab
## BT61 4.51 0.241 72 3.71 5.31 ab
## BT 93 4.53 0.241 72 3.73 5.33 ab
## BT83 4.53 0.241 72 3.73 5.33 ab
## BT85 4.54 0.241 72 3.74 5.34 ab
## FV05 4.54 0.241 72 3.74 5.34 ab
## BT90 4.55 0.241 72 3.74 5.35 ab
## 2011 4.68 0.241 72 3.87 5.48 ab
## BT63 4.70 0.241 72 3.90 5.51 ab
## 4013 4.74 0.241 72 3.93 5.54 ab
## BT 86 4.74 0.241 72 3.94 5.54 ab
## BT64 4.78 0.241 72 3.97 5.58 b
## Bt62 4.78 0.241 72 3.97 5.58 b
##
## Results are averaged over the levels of: block
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$log_mscau <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "log da seca de caules (g)")
#-----------------------------------------------------------------------
Massa seca de raizes (g)
#-----------------------------------------------------------------------
# Declaração do modelo de ANCOVA.
m0 <- lm(msraiz ~ block + gen, data = tb)
MASS::boxcox(m0)
m0 <- update(m0, formula = log(.) ~ .)
# Avaliação dos pressupostos.
par(mfrow = c(2, 2))
plot(m0)
layout(1)
# Quadro de anova.
anova(m0)
## Analysis of Variance Table
##
## Response: log(msraiz)
## Df Sum Sq Mean Sq F value Pr(>F)
## block 2 0.3177 0.158873 1.7371 0.1834354
## gen 36 7.8000 0.216667 2.3689 0.0009788 ***
## Residuals 71 6.4938 0.091461
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## FV10 3.70 0.175 71 3.12 4.29 a
## 3105 3.81 0.175 71 3.23 4.39 a
## BT82 3.93 0.175 71 3.35 4.51 ab
## 3408 4.04 0.175 71 3.46 4.63 abc
## 4016 4.05 0.175 71 3.47 4.63 abc
## 4008 4.07 0.175 71 3.48 4.65 abc
## BT90 4.11 0.175 71 3.53 4.69 abc
## BT70 4.14 0.175 71 3.56 4.73 abc
## 1101 4.20 0.175 71 3.62 4.78 abc
## 4015 4.24 0.175 71 3.66 4.82 abc
## 4017 4.25 0.175 71 3.67 4.83 abc
## 4013 4.30 0.215 71 3.59 5.02 abc
## BT85 PLUS 4.31 0.175 71 3.73 4.90 abc
## BT72 4.32 0.175 71 3.74 4.90 abc
## 1507 4.32 0.175 71 3.74 4.90 abc
## 1102 4.33 0.175 71 3.74 4.91 abc
## BT61 4.33 0.175 71 3.75 4.91 abc
## 1114 4.35 0.175 71 3.77 4.94 abc
## BT98 4.38 0.175 71 3.79 4.96 abc
## 4018 4.38 0.175 71 3.80 4.96 abc
## BT96 4.40 0.175 71 3.82 4.98 abc
## 1106 4.40 0.175 71 3.82 4.98 abc
## BT21 4.42 0.175 71 3.84 5.00 abc
## 1328 4.42 0.175 71 3.84 5.01 abc
## 2011 4.43 0.175 71 3.85 5.01 abc
## BT27 4.52 0.175 71 3.93 5.10 abc
## BT83 4.54 0.175 71 3.96 5.12 abc
## BT 93 4.56 0.175 71 3.98 5.14 abc
## Bt62 4.56 0.175 71 3.98 5.15 abc
## BT68 4.57 0.175 71 3.99 5.15 abc
## FV03 4.59 0.175 71 4.01 5.17 abc
## BT 86 4.61 0.175 71 4.03 5.19 abc
## BT85 4.62 0.175 71 4.04 5.20 abc
## 4010 4.62 0.175 71 4.04 5.21 abc
## FV05 4.84 0.175 71 4.26 5.42 bc
## BT63 4.87 0.175 71 4.29 5.45 c
## BT64 4.89 0.175 71 4.31 5.47 c
##
## Results are averaged over the levels of: block
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means$log_msraiz <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.0f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "log da massa fresca de raizes (g)")
#-----------------------------------------------------------------------
Grau de enraizamento
Análise com ordinal logistic regression
#-----------------------------------------------------------------------
library(MASS)
# Determina os genótipos que não tem variação.
g <- tb %>%
group_by(gen) %>%
summarise(a = diff(range(as.integer(enraz))),
b = all(enraz == "I") | all(enraz == "IV")) %>%
ungroup() %>%
filter(!b) %>%
pull(gen)
# filter(tb, !is.element(gen, g)) %>%
# arrange(gen, block, enraz) %>%
# dplyr::select(gen, block, enraz)
# Genótipos que todas as repetições e uma das classes extremas.
filter(tb, !is.element(gen, g)) %>%
distinct(gen, enraz)
## # A tibble: 8 x 2
## gen enraz
## <fct> <fct>
## 1 FV10 I
## 2 BT70 I
## 3 BT63 IV
## 4 BT27 I
## 5 3408 I
## 6 4008 I
## 7 1106 I
## 8 BT63 <NA>
# Regressão logÃstica ordinal.
m0 <- polr(enraz ~ block + gen,
data = droplevels(filter(tb, is.element(gen, g))),
# data = tb,
Hess = TRUE)
m_null <- update(m0, formula = . ~ block)
# class(m0)
# methods(class = "polr")
# Teste para o efeito dos genótipos.
anova(m0, m_null)
## Likelihood ratio tests of ordinal regression models
##
## Response: enraz
## Model Resid. df Resid. Dev Test Df LR stat. Pr(Chi)
## 1 block 82 220.5416
## 2 block + gen 53 162.4829 1 vs 2 29 58.0587 0.001070333
# Resumo do modelo.
summary(m0)
## Call:
## polr(formula = enraz ~ block + gen, data = droplevels(filter(tb,
## is.element(gen, g))), Hess = TRUE)
##
## Coefficients:
## Value Std. Error t value
## block2 1.079e-01 0.5233 2.061e-01
## block3 6.413e-01 0.5601 1.145e+00
## gen1102 -4.085e+00 1.7098 -2.389e+00
## gen1114 -1.662e+00 1.5864 -1.047e+00
## gen1328 -2.939e+00 1.8526 -1.587e+00
## gen1507 -1.019e+00 1.6602 -6.136e-01
## gen2011 -5.364e+00 1.7972 -2.985e+00
## gen3105 -1.638e+00 1.6133 -1.015e+00
## gen4010 -9.468e-01 1.7200 -5.505e-01
## gen4013 5.022e-05 1.6207 3.099e-05
## gen4015 1.701e+00 1.6642 1.022e+00
## gen4016 -4.542e+00 1.8673 -2.432e+00
## gen4017 -1.019e+00 1.6602 -6.136e-01
## gen4018 -5.366e+00 1.7962 -2.988e+00
## genBT 86 -1.662e+00 1.5864 -1.047e+00
## genBT 93 -2.809e+00 1.6852 -1.667e+00
## genBT21 -5.272e-01 1.7243 -3.058e-01
## genBT61 1.668e-01 1.9460 8.571e-02
## genBt62 -9.468e-01 1.7200 -5.505e-01
## genBT64 -2.809e+00 1.6852 -1.667e+00
## genBT68 -5.182e-02 1.4653 -3.536e-02
## genBT72 -1.662e+00 1.5864 -1.047e+00
## genBT82 -5.364e+00 1.7972 -2.985e+00
## genBT83 -4.542e+00 1.8673 -2.432e+00
## genBT85 -1.379e-01 1.6028 -8.606e-02
## genBT85 PLUS -2.809e+00 1.6852 -1.667e+00
## genBT90 -1.140e+00 1.7158 -6.646e-01
## genBT96 -2.809e+00 1.6852 -1.667e+00
## genBT98 6.302e-01 1.5308 4.117e-01
## genFV03 -5.366e+00 1.7962 -2.988e+00
## genFV05 3.026e-02 1.6227 1.865e-02
##
## Intercepts:
## Value Std. Error t value
## I|II -4.3947 1.3457 -3.2658
## II|III -0.7227 1.2169 -0.5939
## III|IV 1.1175 1.2134 0.9209
##
## Residual Deviance: 162.4829
## AIC: 230.4829
## (3 observations deleted due to missingness)
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df asymp.LCL asymp.UCL .group
## FV03 -3.7831 1.297 Inf -7.853 0.287 a
## 4018 -3.7831 1.297 Inf -7.853 0.287 a
## 2011 -3.7815 1.298 Inf -7.854 0.291 a
## BT82 -3.7815 1.298 Inf -7.854 0.291 a
## BT83 -2.9588 1.408 Inf -7.375 1.457 ab
## 4016 -2.9588 1.408 Inf -7.375 1.457 ab
## 1102 -2.5019 1.187 Inf -6.225 1.221 ab
## 1328 -1.3565 1.413 Inf -5.790 3.077 ab
## BT 93 -1.2261 1.190 Inf -4.958 2.506 ab
## BT64 -1.2261 1.190 Inf -4.958 2.506 ab
## BT85 PLUS -1.2261 1.190 Inf -4.958 2.506 ab
## BT96 -1.2261 1.190 Inf -4.958 2.506 ab
## BT 86 -0.0785 1.088 Inf -3.492 3.335 ab
## BT72 -0.0785 1.088 Inf -3.492 3.335 ab
## 1114 -0.0785 1.088 Inf -3.492 3.335 ab
## 3105 -0.0547 1.131 Inf -3.601 3.492 ab
## BT90 0.4427 1.269 Inf -3.537 4.423 ab
## 1507 0.5644 1.205 Inf -3.217 4.345 ab
## 4017 0.5644 1.205 Inf -3.217 4.345 ab
## 4010 0.6362 1.294 Inf -3.424 4.696 ab
## Bt62 0.6362 1.294 Inf -3.424 4.696 ab
## BT21 1.0558 1.309 Inf -3.052 5.163 ab
## BT85 1.4451 1.130 Inf -2.100 4.990 ab
## BT68 1.5312 0.933 Inf -1.394 4.456 ab
## 1101 1.5830 1.169 Inf -2.084 5.250 ab
## 4013 1.5831 1.169 Inf -2.084 5.250 ab
## FV05 1.6133 1.163 Inf -2.035 5.262 ab
## BT61 1.7498 1.591 Inf -3.240 6.739 ab
## BT98 2.2132 1.034 Inf -1.030 5.457 ab
## 4015 3.2844 1.230 Inf -0.574 7.143 b
##
## Results are averaged over the levels of: block
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 30 estimates
## P value adjustment: tukey method for comparing a family of 30 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean),
.group = str_trim(.group))
tb_emm
## gen emmean SE df asymp.LCL asymp.UCL .group
## 29 FV03 -3.78309031 1.2973931 Inf -7.852688 0.2865077 a
## 13 4018 -3.78309031 1.2973931 Inf -7.852688 0.2865077 a
## 6 2011 -3.78146496 1.2982471 Inf -7.853742 0.2908120 a
## 22 BT82 -3.78146496 1.2982471 Inf -7.853742 0.2908120 a
## 23 BT83 -2.95881903 1.4078387 Inf -7.374857 1.4572194 ab
## 11 4016 -2.95881903 1.4078387 Inf -7.374857 1.4572194 ab
## 2 1102 -2.50189317 1.1870032 Inf -6.225226 1.2214396 ab
## 4 1328 -1.35646685 1.4133808 Inf -5.789889 3.0769558 ab
## 15 BT 93 -1.22606139 1.1898310 Inf -4.958264 2.5061414 ab
## 19 BT64 -1.22606139 1.1898310 Inf -4.958264 2.5061414 ab
## 25 BT85 PLUS -1.22606139 1.1898310 Inf -4.958264 2.5061414 ab
## 27 BT96 -1.22606139 1.1898310 Inf -4.958264 2.5061414 ab
## 14 BT 86 -0.07853765 1.0881929 Inf -3.491927 3.3348519 ab
## 21 BT72 -0.07853765 1.0881929 Inf -3.491927 3.3348519 ab
## 3 1114 -0.07853765 1.0881929 Inf -3.491927 3.3348519 ab
## 7 3105 -0.05467011 1.1305275 Inf -3.600852 3.4915123 ab
## 26 BT90 0.44267769 1.2687762 Inf -3.537156 4.4225119 ab
## 5 1507 0.56435548 1.2054098 Inf -3.216714 4.3454250 ab
## 12 4017 0.56435548 1.2054098 Inf -3.216714 4.3454250 ab
## 8 4010 0.63617901 1.2943932 Inf -3.424009 4.6963671 ab
## 18 Bt62 0.63617901 1.2943932 Inf -3.424009 4.6963671 ab
## 16 BT21 1.05577442 1.3094974 Inf -3.051792 5.1633407 ab
## 24 BT85 1.44507699 1.1300572 Inf -2.099630 4.9897843 ab
## 20 BT68 1.53119237 0.9325331 Inf -1.393931 4.4563157 ab
## 1 1101 1.58301295 1.1691678 Inf -2.084374 5.2504003 ab
## 9 4013 1.58306317 1.1691682 Inf -2.084326 5.2504519 ab
## 30 FV05 1.61327295 1.1630598 Inf -2.034955 5.2615010 ab
## 17 BT61 1.74979613 1.5906748 Inf -3.239753 6.7393457 ab
## 28 BT98 2.21319217 1.0339823 Inf -1.030152 5.4565365 ab
## 10 4015 3.28440738 1.2301079 Inf -0.574134 7.1429487 b
tb_means$enraz <- tb_emm
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = asymp.LCL,
xmax = asymp.UCL,
label = sprintf("%0.2f %s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "Grau de enraizamento")
#-----------------------------------------------------------------------
Análise de componentes principais com valores individuais
#-----------------------------------------------------------------------
str(tb)
## tibble [111 × 13] (S3: tbl_df/tbl/data.frame)
## $ index : num [1:111] 1 2 3 4 5 6 7 8 9 10 ...
## $ gen : Factor w/ 37 levels "1101","1102",..: 12 7 22 37 6 35 32 8 34 24 ...
## $ block : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
## $ altini : num [1:111] 6.5 6.5 11 9.5 4 7.5 7.5 10.5 17 9 ...
## $ alt : num [1:111] 85 109 105 52 59 42 149 59 101 155 ...
## $ diamini: num [1:111] 6.95 3.9 4.58 7.33 4.97 ...
## $ diam : num [1:111] 27.4 20.8 25.2 11.9 17.1 ...
## $ mffol : num [1:111] 456 231 272 161 172 ...
## $ mfcau : num [1:111] 434.8 199.8 276.1 56.5 126.2 ...
## $ msfol : num [1:111] 193 96 120.8 61 55.5 ...
## $ mscau : num [1:111] 206 73.6 131.6 26.3 54.9 ...
## $ msraiz : num [1:111] NA 87.3 104.8 37.9 82.3 ...
## $ enraz : Factor w/ 4 levels "I","II","III",..: 4 2 3 1 2 1 2 3 3 2 ...
# Fórmulas para ajustar os modelos.
formulas <- list(alt = alt ~ altini + block + gen,
diam = diam ~ diamini + block + gen,
mffol = log(mffol) ~ block + gen,
mfcau = log(mfcau) ~ block + gen,
msfol = log(msfol) ~ block + gen,
mscau = log(mscau) ~ block + gen,
msraiz = log(msraiz) ~ block + gen)
# m0 <- lm(formulas[[7]],
# data = tb,
# contrasts = list(block = "contr.sum"))
# m0
# Ajusta o modelo para ajustar para o efeito de covariáveis e blocos.
tb_fit <-
imap(formulas,
function(f, nm) {
# Ajuste do modelo aos dados.
m0 <- lm(f,
data = tb,
contrasts = list(block = "contr.sum"))
# Matriz do modelo para predição.
X <- model.matrix(formula(m0)[-2], data = tb)
# Faz média de termos para obter médias marginais.
a <- range(m0$assign)
i <- !(m0$assign %in% a)
X[, i] <- apply(X[, i],
MARGIN = 2,
FUN = function(x) rep(mean(x), length(x)))
# Médias ajustadas para covariáveis e blocos.
fit <- X %*% coef(m0)
tb_fit <- cbind(tb[, c("gen", "block"), drop = FALSE],
fit = fit)
# Soma resÃduos e médias para retornas observações
# ajustadas.
tb_fit <- left_join(tb_fit,
cbind(m0$model[, c("gen", "block")],
res = residuals(m0)),
by = c("gen", "block")) %>%
replace_na(list(res = 0)) %>%
mutate(fit = fit + res,
res = NULL)
names(tb_fit)[3] <- nm
return(tb_fit)
}) %>%
reduce(left_join)
# Confere.
# anova(lm(alt ~ block + gen, data = tb_fit))
# anova(lm(alt ~ altini + block + gen, data = tb))
#------------------------------------------------------------------------
# Análise de componentes principais.
X <- dplyr::select(tb_fit, -gen, -block)
rownames(X) <- with(tb_fit, sprintf("%s-%s", gen, block))
lattice::splom(X,
auto.key = TRUE,
cex = 0.4,
col = "black",
as.matrix = TRUE)
# Análise de componentes principais.
acp <- princomp(X, cor = TRUE)
# Importância dos componentes e variância explicada.
summary(acp)
## Importance of components:
## Comp.1 Comp.2 Comp.3 Comp.4
## Standard deviation 2.2399318 0.82526012 0.7953390 0.58202883
## Proportion of Variance 0.7167564 0.09729347 0.0903663 0.04839394
## Cumulative Proportion 0.7167564 0.81404983 0.9044161 0.95281007
## Comp.5 Comp.6 Comp.7
## Standard deviation 0.42055733 0.31813428 0.228586185
## Proportion of Variance 0.02526692 0.01445849 0.007464521
## Cumulative Proportion 0.97807699 0.99253548 1.000000000
# Screeplot.
plot(acp, type = "l")
# screeplot(acp, type = "l")
# Carregamentos.
print(acp)
## Call:
## princomp(x = X, cor = TRUE)
##
## Standard deviations:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## 2.2399318 0.8252601 0.7953390 0.5820288 0.4205573 0.3181343 0.2285862
##
## 7 variables and 111 observations.
# Gráfico de biplot.
biplot(acp)
abline(v = 0, h = 0, lty = 2)
# Os componentes ou escores.
head(acp$scores)
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5
## 4013-1 3.7273531 1.19379413 -1.505833339 -0.8713683 -0.88814969
## 2011-1 0.4055854 0.22430132 0.193520858 0.3153508 -0.02234785
## Bt62-1 2.2211729 -0.01744098 -0.003946862 -0.1924202 -0.15583537
## FV10-1 -4.4692922 1.49406059 -0.423757110 0.5674485 0.45085666
## 1507-1 -2.0475050 -0.57181187 -0.489079673 0.1162954 -0.51007943
## FV03-1 -0.2655243 -0.44851195 -1.639958487 -0.4071202 0.19929455
## Comp.6 Comp.7
## 4013-1 -0.19563829 -0.11662529
## 2011-1 0.07174615 0.23074318
## Bt62-1 -0.05739745 -0.14044731
## FV10-1 -0.19972343 -0.03951118
## 1507-1 0.24096648 -0.18219640
## FV03-1 -0.02193009 -0.03212463
# lattice::splom(acp$scores,
# auto.key = TRUE,
# cex = 0.4,
# col = "black",
# as.matrix = TRUE)
# Carregamentos.
acp$loadings
##
## Loadings:
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## alt 0.315 0.228 0.817 0.273 0.306 0.111
## diam 0.390 -0.193 -0.149 -0.611 0.609 0.207
## mffol 0.366 0.431 -0.441 0.362 0.565 -0.182
## mfcau 0.421 0.112 -0.261 -0.545 0.180 0.641
## msfol 0.403 0.284 -0.285 0.193 0.216 -0.721 0.269
## mscau 0.424 0.125 -0.244 -0.435 -0.271 -0.695
## msraiz 0.307 -0.801 0.506
##
## Comp.1 Comp.2 Comp.3 Comp.4 Comp.5 Comp.6 Comp.7
## SS loadings 1.000 1.000 1.000 1.000 1.000 1.000 1.000
## Proportion Var 0.143 0.143 0.143 0.143 0.143 0.143 0.143
## Cumulative Var 0.143 0.286 0.429 0.571 0.714 0.857 1.000
# print(acp$loadings, cutoff = 0.4)
# print(acp$loadings, cutoff = 0.3)
# print(acp$loadings, cutoff = 0.2)
# head(acp$scores) # Escores dos indivÃduos usados no ajuste.
# head(predict(acp)) # Idem.
# predict(acp, newdata = tail(X)) # Pode usar em uma nova amostra.
# NOTE: os resultados da `princomp()` são melhor organizados. Ela
# permite passar matrizes de covariância. Ou seja, pode-se usar outros
# estimadores (mais robustos) da matriz de covariância.
tb_means <- list()
#-----------------------------------------------------------------------
# Componente do tamanho.
m0 <- lm(acp$scores[, 1] ~ tb_fit$gen)
anova(m0)
## Analysis of Variance Table
##
## Response: acp$scores[, 1]
## Df Sum Sq Mean Sq F value Pr(>F)
## tb_fit$gen 36 280.50 7.7918 2.0859 0.003874 **
## Residuals 74 276.42 3.7354
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## BT82 -3.8462 1.12 74 -7.55 -0.138 a
## 3105 -3.4203 1.12 74 -7.13 0.288 ab
## FV10 -3.1134 1.12 74 -6.82 0.595 ab
## 4016 -2.5600 1.12 74 -6.27 1.148 ab
## 4015 -1.5629 1.12 74 -5.27 2.145 ab
## 1114 -1.2036 1.12 74 -4.91 2.505 ab
## 4017 -1.1441 1.12 74 -4.85 2.564 ab
## 3408 -1.0836 1.12 74 -4.79 2.625 ab
## 4018 -1.0109 1.12 74 -4.72 2.697 ab
## BT21 -0.8450 1.12 74 -4.55 2.863 ab
## BT70 -0.6840 1.12 74 -4.39 3.024 ab
## 1102 -0.6603 1.12 74 -4.37 3.048 ab
## 4008 -0.5456 1.12 74 -4.25 3.163 ab
## 4010 -0.5284 1.12 74 -4.24 3.180 ab
## 1328 -0.5167 1.12 74 -4.22 3.192 ab
## 1101 -0.4371 1.12 74 -4.15 3.271 ab
## BT72 -0.2002 1.12 74 -3.91 3.508 ab
## BT96 -0.0133 1.12 74 -3.72 3.695 ab
## BT68 0.0975 1.12 74 -3.61 3.806 ab
## BT98 0.1184 1.12 74 -3.59 3.827 ab
## 1106 0.1699 1.12 74 -3.54 3.878 ab
## FV03 0.2007 1.12 74 -3.51 3.909 ab
## 1507 0.2706 1.12 74 -3.44 3.979 ab
## BT85 PLUS 0.6461 1.12 74 -3.06 4.354 ab
## BT61 0.7321 1.12 74 -2.98 4.440 ab
## BT85 0.9593 1.12 74 -2.75 4.668 ab
## BT27 1.2884 1.12 74 -2.42 4.997 ab
## BT 93 1.3202 1.12 74 -2.39 5.028 ab
## BT90 1.3919 1.12 74 -2.32 5.100 ab
## BT83 1.5100 1.12 74 -2.20 5.218 ab
## FV05 1.7760 1.12 74 -1.93 5.484 ab
## 4013 2.0005 1.12 74 -1.71 5.709 ab
## BT63 2.0106 1.12 74 -1.70 5.719 ab
## BT 86 2.0147 1.12 74 -1.69 5.723 ab
## Bt62 2.1903 1.12 74 -1.52 5.899 b
## 2011 2.2496 1.12 74 -1.46 5.958 b
## BT64 2.4286 1.12 74 -1.28 6.137 b
##
## Results are given on the [ (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## Note: contrasts are still on the [ scale
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means[[1]] <- tb_emm %>%
dplyr::select(gen, emmean, ends_with("CL")) %>%
rename_if(is.numeric, ~str_c(., "_1"))
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.2f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "Componente de tamanho")
#-----------------------------------------------------------------------
# Componente de relação parte áerea e radicular.
m0 <- lm(acp$scores[, 2] ~ tb_fit$gen)
# Comparações múltiplas.
emm <- emmeans(m0, specs = ~gen) %>%
multcomp::cld(Letters = letters, adjust = "tukey", alpha = 0.1)
emm
## gen emmean SE df lower.CL upper.CL .group
## 1328 -1.0364 0.395 74 -2.3497 0.277 a
## 4010 -0.9901 0.395 74 -2.3034 0.323 a
## BT68 -0.8187 0.395 74 -2.1320 0.495 ab
## FV03 -0.7975 0.395 74 -2.1108 0.516 ab
## BT21 -0.7149 0.395 74 -2.0282 0.598 ab
## 1102 -0.6553 0.395 74 -1.9686 0.658 abc
## BT64 -0.6433 0.395 74 -1.9566 0.670 abc
## BT85 -0.5725 0.395 74 -1.8857 0.741 abc
## BT63 -0.5156 0.395 74 -1.8288 0.798 abc
## 4018 -0.5070 0.395 74 -1.8203 0.806 abc
## 1114 -0.5036 0.395 74 -1.8169 0.810 abc
## FV05 -0.3983 0.395 74 -1.7116 0.915 abc
## BT96 -0.1996 0.395 74 -1.5128 1.114 abc
## 4016 -0.1729 0.395 74 -1.4862 1.140 abc
## BT 86 -0.1454 0.395 74 -1.4586 1.168 abc
## BT 93 -0.0831 0.395 74 -1.3964 1.230 abc
## Bt62 -0.0714 0.395 74 -1.3846 1.242 abc
## 1507 -0.0571 0.395 74 -1.3704 1.256 abc
## 4015 -0.0560 0.395 74 -1.3693 1.257 abc
## BT98 0.0268 0.395 74 -1.2865 1.340 abc
## BT72 0.0425 0.395 74 -1.2708 1.356 abc
## BT82 0.0702 0.395 74 -1.2431 1.383 abc
## 1106 0.0846 0.395 74 -1.2287 1.398 abc
## 4017 0.0928 0.395 74 -1.2205 1.406 abc
## BT83 0.1637 0.395 74 -1.1496 1.477 abc
## BT27 0.1904 0.395 74 -1.1228 1.504 abc
## BT61 0.3378 0.395 74 -0.9755 1.651 abc
## 2011 0.4584 0.395 74 -0.8548 1.772 abc
## BT70 0.4985 0.395 74 -0.8147 1.812 abc
## 3408 0.5277 0.395 74 -0.7856 1.841 abc
## BT85 PLUS 0.5392 0.395 74 -0.7740 1.853 abc
## 1101 0.6089 0.395 74 -0.7044 1.922 abc
## 4013 0.7707 0.395 74 -0.5426 2.084 abc
## 3105 0.7897 0.395 74 -0.5236 2.103 abc
## 4008 1.0319 0.395 74 -0.2814 2.345 abc
## BT90 1.2639 0.395 74 -0.0494 2.577 bc
## FV10 1.4411 0.395 74 0.1278 2.754 c
##
## Results are given on the [ (not the response) scale.
## Confidence level used: 0.95
## Conf-level adjustment: sidak method for 37 estimates
## Note: contrasts are still on the [ scale
## P value adjustment: tukey method for comparing a family of 37 estimates
## significance level used: alpha = 0.1
tb_emm <- emm %>%
as.data.frame() %>%
mutate(gen = fct_reorder(gen, emmean))
tb_means[[2]] <- tb_emm %>%
dplyr::select(gen, emmean, ends_with("CL")) %>%
rename_if(is.numeric, ~str_c(., "_2"))
ggplot(data = tb_emm,
mapping = aes(y = gen,
x = emmean,
xmin = lower.CL,
xmax = upper.CL,
label = sprintf("%0.2f%s", emmean, .group))) +
geom_errorbarh() +
geom_point() +
geom_text(nudge_y = 0.25, vjust = 0, size = 3) +
expand_limits(y = c(NA, nlevels(tb_emm$gen) + 1)) +
labs(y = "Genótipos", x = "Componente de relação")
#-----------------------------------------------------------------------
# Juntando.
# tb_means <- tb_means %>%
# reduce(left_join)
#
# ggplot(data = tb_means,
# mapping = aes(x = emmean_1,
# y = emmean_2,
# xmin = lower.CL_1,
# xmax = upper.CL_1,
# ymin = lower.CL_2,
# ymax = upper.CL_2,
# label = gen)) +
# geom_vline(xintercept = 0, color = "red") +
# geom_hline(yintercept = 0, color = "red") +
# geom_errorbar() +
# geom_errorbarh() +
# geom_point() +
# geom_label() +
# labs(x = "Componente 1", y = "Componente 2") +
# coord_equal()
#-----------------------------------------------------------------------