#-----------------------------------------------------------------------
#                                            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á
#                                       2020-abr-16 · Curitiba/PR/Brazil
#-----------------------------------------------------------------------

#-----------------------------------------------------------------------
# Pacotes.

rm(list = objects())
library(tidyverse)
library(readxl)

# knitr::opts_knit$get("rmarkdown.pandoc.to")

1 Importação dos dados

Os dados são importados do arquivo *.xlsx. Após leitura, filtra-se apenas para os genótipos de interesse imediato. São ao todo 35 genótipos para as análises subsequentes.

Neste relatório estão presentes apenas as saídas, ou seja, tabelas e gráficos. Os códigos R não serão exibidos mas podem ser fornecidos se solicitados.

#-----------------------------------------------------------------------
# Importa as tabelas de dados.

# Endereço do arquivo.
xlsx <- "./Walmes' BANCO DE DADOS TD ATUALIZADO_ABR_2020.xlsx"

# Tabela com dados dos 35 genótipos selecionados.
tb <- read_xlsx(path = xlsx,
                sheet = "4) RB08-09 CARVÃO URGENTE",
                skip = 2)
str(tb)
## Classes 'tbl_df', 'tbl' and 'data.frame':    684 obs. of  29 variables:
##  $ id            : num  1 2 3 4 5 6 7 8 9 10 ...
##  $ blc           : num  1 1 1 1 1 1 1 1 1 1 ...
##  $ lin           : num  1 1 1 1 1 2 2 2 2 2 ...
##  $ plt           : num  1 2 3 4 5 1 2 3 4 5 ...
##  $ gen           : chr  "PR093047" "PR093047" "PR093047" "PR093047" ...
##  $ PS_2019-08-09 : num  3 4 1 4 5 6 6 2 5 4 ...
##  $ PCN_2019-08-09: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-08-09 : num  3 4 1 4 5 6 6 2 5 4 ...
##  $ PCT_2019-08-09: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-09-06 : num  6 6 4 4 9 9 10 6 6 7 ...
##  $ PCN_2019-09-06: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-09-06 : num  6 6 4 4 9 9 10 6 6 7 ...
##  $ PCT_2019-09-06: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-10-23 : num  12 14 10 7 14 15 12 9 14 11 ...
##  $ PCN_2019-10-23: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-10-23 : num  12 14 10 7 14 15 12 9 14 11 ...
##  $ PCT_2019-10-23: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-11-28 : num  19 14 11 6 12 20 15 8 13 12 ...
##  $ PCN_2019-11-28: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-11-28 : num  19 14 11 6 12 20 15 8 13 12 ...
##  $ PCT_2019-11-28: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-12-26 : num  11 7 8 6 12 12 8 3 10 13 ...
##  $ PCN_2019-12-26: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-12-26 : num  11 7 8 6 12 12 8 3 10 13 ...
##  $ PCT_2019-12-26: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2020-03-04 : num  10 5 7 5 12 12 6 2 8 10 ...
##  $ PCN_2020-03-04: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2020-03-04 : num  10 5 7 5 12 12 6 2 8 10 ...
##  $ PCT_2020-03-04: num  0 0 0 0 0 0 0 0 0 0 ...
# Vetor com o nome dos 35 genótipos.
sel_gen <- c(na.omit(unique(tb$gen)))

# Tabela com todos os genótipos.
tb <- read_xlsx(path = xlsx,
                sheet = "3) RB08-09 CARVÃO",
                skip = 2)
str(tb)
## Classes 'tbl_df', 'tbl' and 'data.frame':    1140 obs. of  29 variables:
##  $ id            : num  1 2 3 4 5 6 7 8 9 10 ...
##  $ blc           : num  1 1 1 1 1 1 1 1 1 1 ...
##  $ lin           : num  1 1 1 1 1 2 2 2 2 2 ...
##  $ plt           : num  1 2 3 4 5 1 2 3 4 5 ...
##  $ gen           : chr  "PR093047" "PR093047" "PR093047" "PR093047" ...
##  $ PS_2019-08-09 : num  3 4 1 4 5 6 6 2 5 4 ...
##  $ PCN_2019-08-09: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-08-09 : num  3 4 1 4 5 6 6 2 5 4 ...
##  $ PCT_2019-08-09: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-09-06 : num  6 6 4 4 9 9 10 6 6 7 ...
##  $ PCN_2019-09-06: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-09-06 : num  6 6 4 4 9 9 10 6 6 7 ...
##  $ PCT_2019-09-06: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-10-23 : num  12 14 10 7 14 15 12 9 14 11 ...
##  $ PCN_2019-10-23: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-10-23 : num  12 14 10 7 14 15 12 9 14 11 ...
##  $ PCT_2019-10-23: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-11-28 : num  19 14 11 6 12 20 15 8 13 12 ...
##  $ PCN_2019-11-28: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-11-28 : num  19 14 11 6 12 20 15 8 13 12 ...
##  $ PCT_2019-11-28: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-12-26 : num  11 7 8 6 12 12 8 3 10 13 ...
##  $ PCN_2019-12-26: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-12-26 : num  11 7 8 6 12 12 8 3 10 13 ...
##  $ PCT_2019-12-26: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2020-03-04 : num  10 5 7 5 12 12 6 2 8 10 ...
##  $ PCN_2020-03-04: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2020-03-04 : num  10 5 7 5 12 12 6 2 8 10 ...
##  $ PCT_2020-03-04: num  0 0 0 0 0 0 0 0 0 0 ...
# Quantidade de genótipos.
length(unique(tb$gen))
## [1] 56
# Seleciona os 35 do vetor.
tb_sel <- tb %>%
    filter(gen %in% sel_gen)

# Confere a quantidade.
length(unique(tb_sel$gen))
## [1] 35
# Estrutura da tabela.
str(tb_sel)
## Classes 'tbl_df', 'tbl' and 'data.frame':    720 obs. of  29 variables:
##  $ id            : num  1 2 3 4 5 6 7 8 9 10 ...
##  $ blc           : num  1 1 1 1 1 1 1 1 1 1 ...
##  $ lin           : num  1 1 1 1 1 2 2 2 2 2 ...
##  $ plt           : num  1 2 3 4 5 1 2 3 4 5 ...
##  $ gen           : chr  "PR093047" "PR093047" "PR093047" "PR093047" ...
##  $ PS_2019-08-09 : num  3 4 1 4 5 6 6 2 5 4 ...
##  $ PCN_2019-08-09: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-08-09 : num  3 4 1 4 5 6 6 2 5 4 ...
##  $ PCT_2019-08-09: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-09-06 : num  6 6 4 4 9 9 10 6 6 7 ...
##  $ PCN_2019-09-06: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-09-06 : num  6 6 4 4 9 9 10 6 6 7 ...
##  $ PCT_2019-09-06: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-10-23 : num  12 14 10 7 14 15 12 9 14 11 ...
##  $ PCN_2019-10-23: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-10-23 : num  12 14 10 7 14 15 12 9 14 11 ...
##  $ PCT_2019-10-23: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-11-28 : num  19 14 11 6 12 20 15 8 13 12 ...
##  $ PCN_2019-11-28: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-11-28 : num  19 14 11 6 12 20 15 8 13 12 ...
##  $ PCT_2019-11-28: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2019-12-26 : num  11 7 8 6 12 12 8 3 10 13 ...
##  $ PCN_2019-12-26: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2019-12-26 : num  11 7 8 6 12 12 8 3 10 13 ...
##  $ PCT_2019-12-26: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PS_2020-03-04 : num  10 5 7 5 12 12 6 2 8 10 ...
##  $ PCN_2020-03-04: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT_2020-03-04 : num  10 5 7 5 12 12 6 2 8 10 ...
##  $ PCT_2020-03-04: num  0 0 0 0 0 0 0 0 0 0 ...
cap <- "Primeiros 10 registros da tabela de dados com valores observados de perfilhamento e ocorrência de carvão nas datas de avaliação para 35 genótipos selecionados."
knitr::kable(head(tb_sel, n = 10), caption = cap)
Primeiros 10 registros da tabela de dados com valores observados de perfilhamento e ocorrência de carvão nas datas de avaliação para 35 genótipos selecionados.
id blc lin plt gen PS_2019-08-09 PCN_2019-08-09 PT_2019-08-09 PCT_2019-08-09 PS_2019-09-06 PCN_2019-09-06 PT_2019-09-06 PCT_2019-09-06 PS_2019-10-23 PCN_2019-10-23 PT_2019-10-23 PCT_2019-10-23 PS_2019-11-28 PCN_2019-11-28 PT_2019-11-28 PCT_2019-11-28 PS_2019-12-26 PCN_2019-12-26 PT_2019-12-26 PCT_2019-12-26 PS_2020-03-04 PCN_2020-03-04 PT_2020-03-04 PCT_2020-03-04
1 1 1 1 PR093047 3 0 3 0 6 0 6 0 12 0 12 0 19 0 19 0 11 0 11 0 10 0 10 0
2 1 1 2 PR093047 4 0 4 0 6 0 6 0 14 0 14 0 14 0 14 0 7 0 7 0 5 0 5 0
3 1 1 3 PR093047 1 0 1 0 4 0 4 0 10 0 10 0 11 0 11 0 8 0 8 0 7 0 7 0
4 1 1 4 PR093047 4 0 4 0 4 0 4 0 7 0 7 0 6 0 6 0 6 0 6 0 5 0 5 0
5 1 1 5 PR093047 5 0 5 0 9 0 9 0 14 0 14 0 12 0 12 0 12 0 12 0 12 0 12 0
6 1 2 1 PR093047 6 0 6 0 9 0 9 0 15 0 15 0 20 0 20 0 12 0 12 0 12 0 12 0
7 1 2 2 PR093047 6 0 6 0 10 0 10 0 12 0 12 0 15 0 15 0 8 0 8 0 6 0 6 0
8 1 2 3 PR093047 2 0 2 0 6 0 6 0 9 0 9 0 8 0 8 0 3 0 3 0 2 0 2 0
9 1 2 4 PR093047 5 0 5 0 6 0 6 0 14 0 14 0 13 0 13 0 10 0 10 0 8 0 8 0
10 1 2 5 PR093047 4 0 4 0 7 0 7 0 11 0 11 0 12 0 12 0 13 0 13 0 10 0 10 0

2 Redisposição dos dados

Para fazer uso dos dados é necessário proceder com modificação na disposição tabular dos mesmos. As observações da mesma variável foram registradas como colunas lado a lado ao longo das avaliações (avaliações estão nas colunas). Os vetores das avaliações serão empilhados e depois será desempilhado os vetores das respostas para que a data de avaliação esteja registrada em uma coluna (avaliação nas linhas).

#-----------------------------------------------------------------------
# Arrumação dos dados para análise.

# Usa nome curto pro objeto.
tb <- tb_sel
rm(tb_sel)

# Converte variáveis para fator.
tb <- tb %>%
    mutate(blc = factor(blc),
           gen = factor(gen))

# Quantidade de genótipos.
nlevels(tb$gen)
## [1] 35
# Número de avaliações.
length(na.omit(unique(str_extract(names(tb), "_.*$"))))
## [1] 6
# Nomes das variáveis para devolver a ordem certa.
v <- names(tb) %>%
    str_remove("_.*$") %>%
    unique()

# Empilha as avaliações e desempilha nas variáveis observadas.
tbu <- tb %>%
    gather(key = "var", value = "val", -(id:gen)) %>%
    separate(col = "var", into = c("resp", "avl"), sep = "_") %>%
    mutate(avl = as.Date(avl)) %>%
    spread(key = "resp", value = "val") %>%
    arrange(id, avl) %>%
    drop_na() %>%
    select(all_of(v), everything())
str(tbu)
## Classes 'tbl_df', 'tbl' and 'data.frame':    4320 obs. of  10 variables:
##  $ id : num  1 1 1 1 1 1 2 2 2 2 ...
##  $ blc: Factor w/ 2 levels "1","2": 1 1 1 1 1 1 1 1 1 1 ...
##  $ lin: num  1 1 1 1 1 1 1 1 1 1 ...
##  $ plt: num  1 1 1 1 1 1 2 2 2 2 ...
##  $ gen: Factor w/ 35 levels "PR081019","PR081027",..: 31 31 31 31 31 31 31 31 31 31 ...
##  $ PS : num  3 6 12 19 11 10 4 6 14 14 ...
##  $ PCN: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ PT : num  3 6 12 19 11 10 4 6 14 14 ...
##  $ PCT: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ avl: Date, format: "2019-08-09" ...
# Número de registros por avaliação.
tbu %>%
    count(avl)
## # A tibble: 6 x 2
##   avl            n
##   <date>     <int>
## 1 2019-08-09   720
## 2 2019-09-06   720
## 3 2019-10-23   720
## 4 2019-11-28   720
## 5 2019-12-26   720
## 6 2020-03-04   720
# Número de genótipos.
nlevels(tbu$gen)
## [1] 35
cap <- "Primeiros 15 registros da tabela preparada para análise dos dados de perfilhamento e ocorrência de carvão para os 35 genótipos."
knitr::kable(head(tbu, n = 15), caption = cap)
Primeiros 15 registros da tabela preparada para análise dos dados de perfilhamento e ocorrência de carvão para os 35 genótipos.
id blc lin plt gen PS PCN PT PCT avl
1 1 1 1 PR093047 3 0 3 0 2019-08-09
1 1 1 1 PR093047 6 0 6 0 2019-09-06
1 1 1 1 PR093047 12 0 12 0 2019-10-23
1 1 1 1 PR093047 19 0 19 0 2019-11-28
1 1 1 1 PR093047 11 0 11 0 2019-12-26
1 1 1 1 PR093047 10 0 10 0 2020-03-04
2 1 1 2 PR093047 4 0 4 0 2019-08-09
2 1 1 2 PR093047 6 0 6 0 2019-09-06
2 1 1 2 PR093047 14 0 14 0 2019-10-23
2 1 1 2 PR093047 14 0 14 0 2019-11-28
2 1 1 2 PR093047 7 0 7 0 2019-12-26
2 1 1 2 PR093047 5 0 5 0 2020-03-04
3 1 1 3 PR093047 1 0 1 0 2019-08-09
3 1 1 3 PR093047 4 0 4 0 2019-09-06
3 1 1 3 PR093047 10 0 10 0 2019-10-23 
cap <- "Número de registros em cada avaliação por genótipo em cada bloco. Cada registro corresponde à observação de uma planta por parcela. Cada parcela tem duas fileiras com 5 plantas em cada."

# Exibição gráfica da ocorrência.
ggplot(data = filter(tbu, avl == avl[1]),
       mapping = aes(x = reorder(gen, gen, length), fill = blc)) +
    geom_bar(position = "dodge") +
    coord_flip() +
    labs(x = "Genótipos", y = "Número de registros", fill = "Bloco")
Número de registros em cada avaliação por genótipo em cada bloco. Cada registro corresponde à observação de uma planta por parcela. Cada parcela tem duas fileiras com 5 plantas em cada.

Número de registros em cada avaliação por genótipo em cada bloco. Cada registro corresponde à observação de uma planta por parcela. Cada parcela tem duas fileiras com 5 plantas em cada.

3 Agregação por unidade experimental

As variáveis respostas foram registradas ao nível de planta (5 plantas) dentro de linha (duas linhas) dentro de parcela, perfazendo assim, 10 observações por unidade experimental em cada avaliação. As plantas que foram perdidas apresentam número de perfilhos igual a 0. Então os registros de valores ausentes precisam ser apropriadamente tratados para não implicar em viéses de tamanho amostral. Após isso, valores são agregados por unidade experimental.

#-----------------------------------------------------------------------
# Define as variáveis resposta para a análise.

# NOTE: Legenda das variáveis resposta.
# PS  : Número de perfilhos sadios (numerador).
# PT  : Número de perfilhos total (denominador).
# PCN : Número de perfilhos novos com carvão.
# PCT : Número de perfilhos com carvão.

# As parcelas perdidas apresentam 0 para a variável PT.
sum(tbu$PT == 0)
## [1] 29
# Elimina as parcelas perdidas.
tbu <- tbu %>%
    filter(PT > 0)

# Agrega por unidade experimental.
tbua <- tbu %>%
    group_by(blc, gen, avl) %>%
    summarise(PT = mean(PT),
              PCT = mean(PCT),
              perc = 100 * PCT/PT,
              n = n()) %>%
    ungroup()
str(tbua)
## Classes 'tbl_df', 'tbl' and 'data.frame':    420 obs. of  7 variables:
##  $ blc : Factor w/ 2 levels "1","2": 1 1 1 1 1 1 1 1 1 1 ...
##  $ gen : Factor w/ 35 levels "PR081019","PR081027",..: 1 1 1 1 1 1 2 2 2 2 ...
##  $ avl : Date, format: "2019-08-09" ...
##  $ PT  : num  2.6 5.2 13.3 13.6 9.8 7.1 4.8 6.6 12.5 16 ...
##  $ PCT : num  0 0 0 0 0 0 0 0 0 0 ...
##  $ perc: num  0 0 0 0 0 0 0 0 0 0 ...
##  $ n   : int  10 10 10 10 10 10 10 10 10 10 ...
cap <- "Primeiros 15 registros da tabela preparada contendo o número de unidades experimentais por parcela (n) e o percentual de perfilhos com carvão (perc) para os 35 genótipos em cada avaliação por parcela."
knitr::kable(head(tbua, n = 15), caption = cap)
Primeiros 15 registros da tabela preparada contendo o número de unidades experimentais por parcela (n) e o percentual de perfilhos com carvão (perc) para os 35 genótipos em cada avaliação por parcela.
blc gen avl PT PCT perc n
1 PR081019 2019-08-09 2.6 0 0 10
1 PR081019 2019-09-06 5.2 0 0 10
1 PR081019 2019-10-23 13.3 0 0 10
1 PR081019 2019-11-28 13.6 0 0 10
1 PR081019 2019-12-26 9.8 0 0 10
1 PR081019 2020-03-04 7.1 0 0 10
1 PR081027 2019-08-09 4.8 0 0 10
1 PR081027 2019-09-06 6.6 0 0 10
1 PR081027 2019-10-23 12.5 0 0 10
1 PR081027 2019-11-28 16.0 0 0 10
1 PR081027 2019-12-26 11.6 0 0 10
1 PR081027 2020-03-04 11.0 0 0 10
1 PR081028 2019-08-09 5.0 0 0 10
1 PR081028 2019-09-06 7.0 0 0 10
1 PR081028 2019-10-23 14.0 0 0 10 
# Conta quantas parcelas com cada quantidade de plantas.
tbua %>%
    count(n)
## # A tibble: 6 x 2
##       n    nn
##   <int> <int>
## 1     7     3
## 2     8     2
## 3     9    15
## 4    10   388
## 5    19     1
## 6    20    11
cap <- "Número médio de perfilhos por unidade experimental em cada genótipo ao longo das avaliações. A linha contínua conecta a média amostral dos valores observadados nas parcelas avaliadas em uma mesma data."

# Séries com o número médio de perfilhos.
ggplot(data = tbua,
       mapping = aes(x = avl, y = PT)) +
    facet_wrap(facets = ~gen) +
    geom_point() +
    stat_summary(geom = "line", fun = "mean") +
    labs(x = "Avaliações", y = "Número médio de perfilhos")
Número médio de perfilhos por unidade experimental em cada genótipo ao longo das avaliações. A linha contínua conecta a média amostral dos valores observadados nas parcelas avaliadas em uma mesma data.

Número médio de perfilhos por unidade experimental em cada genótipo ao longo das avaliações. A linha contínua conecta a média amostral dos valores observadados nas parcelas avaliadas em uma mesma data. 

cap <- "Proporção de perfilhos com carão por unidade experimental em cada genótipo ao longo das avaliações. A linha contínua conecta a média amostral dos valores observadados nas parcelas avaliadas em uma mesma data. Os genótipos estão ordenados conforme a ocorrência de carvão."

# Séries com a proporção de perfilhos com carvão.
ggplot(data = mutate(tbua, gen = fct_reorder(gen, perc, sum)),
       mapping = aes(x = avl, y = perc)) +
    facet_wrap(facets = ~gen) +
    geom_point() +
    stat_summary(geom = "line", fun = "mean") +
    labs(x = "Avaliações", y = "Percentual de perfilhos com carvão",
         caption = "Genótipos ordenados pela ocorrência de carvão")
Proporção de perfilhos com carão por unidade experimental em cada genótipo ao longo das avaliações. A linha contínua conecta a média amostral dos valores observadados nas parcelas avaliadas em uma mesma data. Os genótipos estão ordenados conforme a ocorrência de carvão.

Proporção de perfilhos com carão por unidade experimental em cada genótipo ao longo das avaliações. A linha contínua conecta a média amostral dos valores observadados nas parcelas avaliadas em uma mesma data. Os genótipos estão ordenados conforme a ocorrência de carvão.

4 Definição das variáveis resposta

Para a análise com a finalidade principal de comparar os genótipos, é interessante eliminar as demais dimensões, como as avaliações ao longo do tempo. Sendo assim, os dados serão agregados para cada genótipo e bloco considerando a média das variáveis resposta.

Serão obtidas as variáveis para o nível de parcela:

  • O percentual médio de ocorrência de carvão nas avaliações (perc_med).
  • O percentual máximo de ocorrência de carvão nas avaliações (perc_max).
  • Quantidade máxima de perfilhos nas avaliações (perf_max).
  • Diferença entre a média de perfilhos das avaliações 3 e 4 e o número de perfilhos na avaliação 6 (dif_perf).
#-----------------------------------------------------------------------
# Marginal para genótipos.

# Obter a média de proporção de perfilhos com carvão nas avaliações para
# cada unidade experimental.
tbua_ue <- tbua %>%
    group_by(blc, gen) %>%
    arrange(avl) %>%
    summarise(perc_med = weighted.mean(x = perc, w = n),
              perc_max = max(perc),
              perf_med = weighted.mean(x = PT, w = n),
              perf_max = max(PT),
              dif_perf = mean(PT[3:4]) - PT[6])

# Fazendo a agregração por genótipo para apresentar gráficos.
tbua_gen <- tbua_ue %>%
    group_by(gen) %>%
    summarise(perc_med = mean(perc_med),
              perc_max = max(perc_max),
              perf_med = mean(perf_med),
              perf_max = max(perf_max),
              dif_perf = mean(dif_perf))
cap <- "Valores das variáveis percentual médio (perc_med) e (perc_max) de ocorrência de carvão, número máximo (perf_max) e médio (perf_med) de perfilhos e diferença no número de perfilhos entre as avaliações 3 e 4 contra a avaliação 6 (dif_perf) para os 35 genótipos selecionados."
knitr::kable(tbua_gen, caption = cap)
Valores das variáveis percentual médio (perc_med) e (perc_max) de ocorrência de carvão, número máximo (perf_max) e médio (perf_med) de perfilhos e diferença no número de perfilhos entre as avaliações 3 e 4 contra a avaliação 6 (dif_perf) para os 35 genótipos selecionados.
gen perc_med perc_max perf_med perf_max dif_perf
PR081019 0.0000000 0.000000 9.966667 17.4 6.475000
PR081027 0.0000000 0.000000 10.358333 16.0 3.175000
PR081028 0.0000000 0.000000 10.975000 18.2 8.025000
PR081031 5.3046339 13.432836 8.050000 13.5 3.800000
PR081032 0.0000000 0.000000 8.283333 15.6 7.550000
PR081043 1.5841053 6.172840 9.941667 17.4 7.600000
PR081047 1.5616176 7.352941 9.116667 17.6 6.625000
PR081050 0.0000000 0.000000 11.508333 18.1 7.025000
PR081185 1.2011688 7.500000 6.845300 14.0 3.826191
PR083007 0.0000000 0.000000 9.725000 16.0 6.100000
PR085003 0.0000000 0.000000 8.758333 14.9 5.250000
PR091104 0.4316337 1.904762 10.390819 20.0 5.016667
PR091108 0.0000000 0.000000 12.341667 21.4 8.875000
PR091114 0.0000000 0.000000 12.175000 23.7 8.025000
PR091116 0.0000000 0.000000 7.983333 15.3 2.600000
PR091119 7.9395823 18.292683 12.143629 21.6 7.344444
PR091141 13.5763416 39.534884 11.798164 21.0 8.822222
PR091142 0.0000000 0.000000 10.780508 18.9 8.316667
PR091143 1.5999092 7.954546 9.516667 18.6 5.975000
PR091166 0.1960152 1.030928 11.325666 20.2 7.426389
PR091169 0.0000000 0.000000 8.858333 15.1 6.975000
PR091171 0.0000000 0.000000 9.833333 20.9 8.800000
PR091179 6.2626598 14.102564 12.273539 18.8 3.978571
PR091183 2.3619396 5.785124 13.766667 25.4 11.800000
PR091184 0.0000000 0.000000 11.775000 19.3 8.300000
PR091204 0.0000000 0.000000 12.758333 24.7 9.275000
PR091218 1.9249484 6.593407 10.110714 17.1 6.594079
PR091219 0.0000000 0.000000 12.825000 24.0 9.125000
PR091220 1.4589320 7.692308 9.780932 17.3 6.044444
PR091228 0.2922529 1.388889 9.766667 19.7 10.025000
PR093047 2.0651407 10.769231 8.366667 13.9 5.400000
PR093053 1.5456017 9.230769 7.908333 14.8 3.725000
PR097052 0.0000000 0.000000 5.825000 13.6 3.500000
PR099448 0.0000000 0.000000 11.233333 20.0 9.400000
PRP08107 0.1527992 1.136364 8.050847 15.0 5.652778 
cap <- "Percentual médio de perfilhos com carvão por genótipo."

ggplot(data = tbua_gen,
       mapping = aes(y = reorder(gen, perc_med),
                     x = perc_med)) +
    geom_col() +
    geom_text(mapping = aes(label = sprintf("%0.2f", perc_med)),
              hjust = -0.1) +
    expand_limits(x = c(0, 15)) +
    labs(x = "Percentual de perfilhos com carvão",
         y = "Genótipos",
         caption = "Valores médios sobre blocos e avaliações")
Percentual médio de perfilhos com carvão por genótipo.

Percentual médio de perfilhos com carvão por genótipo. 

cap <- "Número máximo de perfilhos por genótipo."

ggplot(data = tbua_gen,
       mapping = aes(y = reorder(gen, perf_max),
                     x = perf_max)) +
    geom_col() +
    geom_text(mapping = aes(label = sprintf("%0.1f", perf_max)),
              hjust = -0.1) +
    expand_limits(x = c(0, 28)) +
    labs(x = "Número de perfilhos",
         y = "Genótipos",
         caption = "Valores agregados sobre blocos e avaliações")
Número máximo de perfilhos por genótipo.

Número máximo de perfilhos por genótipo. 

cap <- "Diferença no número de perfilhos médio entre as avaliações 3 e 4 e o número de perfilhos da avaliação 6 para cada genótipo."

ggplot(data = tbua_gen,
       mapping = aes(y = reorder(gen, dif_perf),
                     x = dif_perf)) +
    geom_col() +
    geom_text(mapping = aes(label = sprintf("%0.1f", dif_perf)),
              hjust = -0.1) +
    expand_limits(x = c(0, 13)) +
    labs(x = "Diferença no número de perfilhos",
         y = "Genótipos",
         caption = "Valores agragados sobre blocos e avaliações")
Diferença no número de perfilhos médio entre as avaliações 3 e 4 e o número de perfilhos da avaliação 6 para cada genótipo.

Diferença no número de perfilhos médio entre as avaliações 3 e 4 e o número de perfilhos da avaliação 6 para cada genótipo.

5 Análise estatística do número de perfilhos

Para comparação dos genótipos em relação a emissão de perfilhos, será feita a análise de variância seguido de comparações múltiplas de médias para a variável número máximo de perfilhos. Os valores foram agregrados por unidade experimental o que corresponde a um experimento em delineamento de blocos casualizados completos.

library(emmeans)

num2cld <- function(x) {
    u <- strsplit(trimws(x), split = "")
    n <- as.integer(unlist(u))
    stopifnot(!any(is.na(n)))
    k <- max(n)
    sapply(u,
           FUN = function(ui) {
               paste(letters[k:1][rev(as.integer(ui))],
                     collapse = "")
           })
}

# Ajusta o modelo aos dados.
m0 <- lm(perf_med ~ blc + gen, data = tbua_ue)

# Verifica os pressupostos.
par(mfrow = c(2, 2))
plot(m0)

layout(1)

# Quadro de análise de variância.
anova(m0)
## Analysis of Variance Table
## 
## Response: perf_med
##           Df  Sum Sq Mean Sq F value    Pr(>F)    
## blc        1   8.120  8.1204  4.7534   0.03626 *  
## gen       34 231.817  6.8181  3.9911 5.621e-05 ***
## Residuals 34  58.083  1.7083                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Aplica o teste de médias.
emm <- emmeans(m0, spec = ~gen)
tb_means <- multcomp::cld(emm) %>%
    as.data.frame() %>%
    mutate(cld = num2cld(.group),
           .group = NULL)
cap <- "Quadro de análise de variância para o número médio de perfilhos."

knitr::kable(anova(m0), caption = cap)
Quadro de análise de variância para o número médio de perfilhos.
Df Sum Sq Mean Sq F value Pr(>F)
blc 1 8.120411 8.120411 4.753443 0.0362569
gen 34 231.816579 6.818135 3.991130 0.0000562
Residuals 34 58.082940 1.708322 NA NA 
cap <- "Tabela com as médias (emmean) marginais para o número de perfilhos por genótipo. valores seguidos de mesma letra (cld) não diferem entre si pelo teste de Tukey a 5%."

knitr::kable(tb_means, caption = cap)
Tabela com as médias (emmean) marginais para o número de perfilhos por genótipo. valores seguidos de mesma letra (cld) não diferem entre si pelo teste de Tukey a 5%.
gen emmean SE df lower.CL upper.CL cld
PR097052 5.825000 0.9242082 34 3.946783 7.703217 d
PR081185 6.845300 0.9242082 34 4.967083 8.723518 cd
PR093053 7.908333 0.9242082 34 6.030116 9.786551 bcd
PR091116 7.983333 0.9242082 34 6.105116 9.861550 bcd
PR081031 8.050000 0.9242082 34 6.171783 9.928217 bcd
PRP08107 8.050847 0.9242082 34 6.172630 9.929065 bcd
PR081032 8.283333 0.9242082 34 6.405116 10.161551 bcd
PR093047 8.366667 0.9242082 34 6.488449 10.244884 abcd
PR085003 8.758333 0.9242082 34 6.880116 10.636551 abcd
PR091169 8.858333 0.9242082 34 6.980116 10.736550 abcd
PR081047 9.116667 0.9242082 34 7.238449 10.994884 abcd
PR091143 9.516667 0.9242082 34 7.638450 11.394884 abcd
PR083007 9.725000 0.9242082 34 7.846783 11.603217 abcd
PR091228 9.766667 0.9242082 34 7.888450 11.644884 abcd
PR091220 9.780932 0.9242082 34 7.902715 11.659149 abcd
PR091171 9.833333 0.9242082 34 7.955116 11.711550 abcd
PR081043 9.941667 0.9242082 34 8.063450 11.819884 abcd
PR081019 9.966667 0.9242082 34 8.088449 11.844884 abcd
PR091218 10.110714 0.9242082 34 8.232497 11.988931 abcd
PR081027 10.358333 0.9242082 34 8.480116 12.236550 abcd
PR091104 10.390819 0.9242082 34 8.512602 12.269036 abcd
PR091142 10.780508 0.9242082 34 8.902291 12.658726 abcd
PR081028 10.975000 0.9242082 34 9.096783 12.853217 abcd
PR099448 11.233333 0.9242082 34 9.355116 13.111550 abcd
PR091166 11.325666 0.9242082 34 9.447449 13.203883 abc
PR081050 11.508333 0.9242082 34 9.630116 13.386551 abc
PR091184 11.775000 0.9242082 34 9.896783 13.653217 abc
PR091141 11.798164 0.9242082 34 9.919947 13.676381 abc
PR091119 12.143629 0.9242082 34 10.265412 14.021847 abc
PR091114 12.175000 0.9242082 34 10.296783 14.053217 abc
PR091179 12.273539 0.9242082 34 10.395322 14.151756 abc
PR091108 12.341667 0.9242082 34 10.463449 14.219884 ab
PR091204 12.758333 0.9242082 34 10.880116 14.636551 ab
PR091219 12.825000 0.9242082 34 10.946783 14.703217 ab
PR091183 13.766667 0.9242082 34 11.888450 15.644884 a 
cap <- "Número médio de perfilhos por genótipo com intervalo de confiança de 95%. Valores que possuem ao menos uma letra igual não diferem entre si pelo teste de Tukey a 5% de significância."

ggplot(data = tb_means,
       mapping = aes(y = reorder(gen, emmean), x = emmean)) +
    geom_errorbarh(mapping = aes(xmin = lower.CL, xmax = upper.CL),
                   height = 0.4) +
    geom_point() +
    geom_label(mapping = aes(label = sprintf("%0.2f %s", emmean, cld)),
               hjust = 0,
               nudge_x = 0.2,
               size = 3) +
    labs(x = "Número de perfilhos",
         y = "Genótipos")
Número médio de perfilhos por genótipo com intervalo de confiança de 95%. Valores que possuem ao menos uma letra igual não diferem entre si pelo teste de Tukey a 5% de significância.

Número médio de perfilhos por genótipo com intervalo de confiança de 95%. Valores que possuem ao menos uma letra igual não diferem entre si pelo teste de Tukey a 5% de significância.

6 Análise estatística da diferença em número de perfilhos

O mesmo procedimento acima será repetido para a variável que é a diferença no número de perfilhos entre as avaliações 3 e 4 e a avaliação 6.

# Ajusta o modelo aos dados.
m0 <- lm(dif_perf ~ blc + gen, data = tbua_ue)

# Verifica os pressupostos.
par(mfrow = c(2, 2))
plot(m0)

layout(1)

# Quadro de análise de variância.
anova(m0)
## Analysis of Variance Table
## 
## Response: dif_perf
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## blc        1  14.43 14.4294  6.0563  0.01909 *  
## gen       34 326.50  9.6030  4.0306 5.06e-05 ***
## Residuals 34  81.01  2.3826                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# Aplica o teste de médias.
emm <- emmeans(m0, spec = ~gen)
tb_means <- multcomp::cld(emm) %>%
    as.data.frame() %>%
    mutate(cld = num2cld(.group),
           .group = NULL)
cap <- "Quadro de análise de variância para a diferença no número de perfilhos."

knitr::kable(anova(m0), caption = cap)
Quadro de análise de variância para a diferença no número de perfilhos.
Df Sum Sq Mean Sq F value Pr(>F)
blc 1 14.42940 14.429402 6.056268 0.0190905
gen 34 326.50273 9.603022 4.030553 0.0000506
Residuals 34 81.00693 2.382557 NA NA 
cap <- "Tabela com as médias (emmean) marginais para a diferença no número de perfilhos por genótipo. valores seguidos de mesma letra (cld) não diferem entre si pelo teste de Tukey a 5%."

knitr::kable(tb_means, caption = cap)
Tabela com as médias (emmean) marginais para a diferença no número de perfilhos por genótipo. valores seguidos de mesma letra (cld) não diferem entre si pelo teste de Tukey a 5%.
gen emmean SE df lower.CL upper.CL cld
PR091116 2.600000 1.091457 34 0.3818925 4.818107 d
PR081027 3.175000 1.091457 34 0.9568925 5.393108 cd
PR097052 3.500000 1.091457 34 1.2818925 5.718108 cd
PR093053 3.725000 1.091457 34 1.5068925 5.943107 bcd
PR081031 3.800000 1.091457 34 1.5818925 6.018108 bcd
PR081185 3.826191 1.091457 34 1.6080830 6.044298 bcd
PR091179 3.978571 1.091457 34 1.7604639 6.196679 bcd
PR091104 5.016667 1.091457 34 2.7985592 7.234774 bcd
PR085003 5.250000 1.091457 34 3.0318925 7.468108 bcd
PR093047 5.400000 1.091457 34 3.1818925 7.618107 abcd
PRP08107 5.652778 1.091457 34 3.4346703 7.870885 abcd
PR091143 5.975000 1.091457 34 3.7568925 8.193107 abcd
PR091220 6.044444 1.091457 34 3.8263370 8.262552 abcd
PR083007 6.100000 1.091457 34 3.8818925 8.318107 abcd
PR081019 6.475000 1.091457 34 4.2568925 8.693107 abcd
PR091218 6.594079 1.091457 34 4.3759715 8.812186 abcd
PR081047 6.625000 1.091457 34 4.4068925 8.843108 abcd
PR091169 6.975000 1.091457 34 4.7568925 9.193107 abcd
PR081050 7.025000 1.091457 34 4.8068925 9.243108 abcd
PR091119 7.344444 1.091457 34 5.1263370 9.562552 abcd
PR091166 7.426389 1.091457 34 5.2082814 9.644496 abcd
PR081032 7.550000 1.091457 34 5.3318925 9.768107 abcd
PR081043 7.600000 1.091457 34 5.3818925 9.818107 abcd
PR081028 8.025000 1.091457 34 5.8068925 10.243108 abcd
PR091114 8.025000 1.091457 34 5.8068925 10.243108 abcd
PR091184 8.300000 1.091457 34 6.0818925 10.518107 abcd
PR091142 8.316667 1.091457 34 6.0985592 10.534774 abcd
PR091171 8.800000 1.091457 34 6.5818925 11.018107 abcd
PR091141 8.822222 1.091457 34 6.6041147 11.040330 abcd
PR091108 8.875000 1.091457 34 6.6568925 11.093108 abcd
PR091219 9.125000 1.091457 34 6.9068925 11.343108 abc
PR091204 9.275000 1.091457 34 7.0568925 11.493108 abc
PR099448 9.400000 1.091457 34 7.1818925 11.618108 abc
PR091228 10.025000 1.091457 34 7.8068925 12.243108 ab
PR091183 11.800000 1.091457 34 9.5818925 14.018107 a 
cap <- "Diferença no número de perfilhos perfilhos por genótipo com intervalo de confiança de 95%. Valores que possuem ao menos uma letra igual não diferem entre si pelo teste de Tukey a 5% de significância."

ggplot(data = tb_means,
       mapping = aes(y = reorder(gen, emmean), x = emmean)) +
    geom_errorbarh(mapping = aes(xmin = lower.CL, xmax = upper.CL),
                   height = 0.4) +
    geom_point() +
    geom_label(mapping = aes(label = sprintf("%0.2f %s", emmean, cld)),
               hjust = 0,
               nudge_x = 0.2,
               size = 3) +
    labs(x = "Diferença no número de perfilhos",
         y = "Genótipos")
Diferença no número de perfilhos perfilhos por genótipo com intervalo de confiança de 95%. Valores que possuem ao menos uma letra igual não diferem entre si pelo teste de Tukey a 5% de significância.

Diferença no número de perfilhos perfilhos por genótipo com intervalo de confiança de 95%. Valores que possuem ao menos uma letra igual não diferem entre si pelo teste de Tukey a 5% de significância.

7 Análise estatística da ocorrência de carvão

Por não haver considerável ocorrência de carvão devido a maior porção dos genótipos ser resistente, optou-se por realizar a análise estatística. Caso ela tivesse que ser feita, a abordagem mais razoável seria considerar penas os genótipos que apresentaram alguma ocorrência de carvão, já que os que não apresentaram, devido a isso, tem variabilidade zero o que pode subestimar consideravelmente a variância residual.

O gráfico abaixo já foi apresentado anteriormente mas repetiu-se aqui por conveniência.

cap <- "Percentual médio de perfilhos com carvão por genótipo."

ggplot(data = tbua_gen,
       mapping = aes(y = reorder(gen, perc_med),
                     x = perc_med)) +
    geom_col() +
    geom_text(mapping = aes(label = sprintf("%0.2f", perc_med)),
              hjust = -0.1) +
    expand_limits(x = c(0, 15)) +
    labs(x = "Percentual de perfilhos com carvão",
         y = "Genótipos",
         caption = "Valores médios sobre blocos e avaliações")
Percentual médio de perfilhos com carvão por genótipo.

Percentual médio de perfilhos com carvão por genótipo.