Análise a partir dos dados brutos das Doenças

library(tidyverse)
library(dplyr)
library(agricolae)
library(psych)

## Importação e limpeza dos dados
df = read.csv("Dataset.csv", stringsAsFactors = FALSE)


################### Filtering treatments ###################
target <- c("T1R1", "T1R2", "T1R3", "T1R4", "T4R1", "T4R2", "T4R3", "T4R4",
            "T5R1", "T5R2", "T5R3", "T5R4", "T6R1", "T6R2", "T6R3", "T6R4",
            "T11R1", "T11R2", "T11R3", "T11R4", 
            "T12R1", "T12R2", "T12R3", "T12R4")
df <- df %>%
    filter(TRATAMENTO %in% target)

#### Função para renomear os tratamentos
###################################################################################
my_fun <- function(x) { 
  x %>%       
    mutate(TRATAMENTO = case_when(TRATAMENTO == 'T1R1' ~'T1',
                                  TRATAMENTO == 'T1R2' ~'T1',
                                  TRATAMENTO == 'T1R3' ~'T1',
                                  TRATAMENTO == 'T1R4' ~'T1',
                                  TRATAMENTO == 'T2R1' ~'T2',
                                  TRATAMENTO == 'T2R2' ~'T2',
                                  TRATAMENTO == 'T2R3' ~'T2',
                                  TRATAMENTO == 'T2R4' ~'T2',
                                  TRATAMENTO == 'T3R1' ~'T3',
                                  TRATAMENTO == 'T3R2' ~'T3',
                                  TRATAMENTO == 'T3R3' ~'T3',
                                  TRATAMENTO == 'T3R4' ~'T3',
                                  TRATAMENTO == 'T4R1' ~'T4',
                                  TRATAMENTO == 'T4R2' ~'T4',
                                  TRATAMENTO == 'T4R3' ~'T4',
                                  TRATAMENTO == 'T4R4' ~'T4',
                                  TRATAMENTO == 'T5R1' ~'T5',
                                  TRATAMENTO == 'T5R2' ~'T5',
                                  TRATAMENTO == 'T5R3' ~'T5',
                                  TRATAMENTO == 'T5R4' ~'T5',
                                  TRATAMENTO == 'T6R1' ~'T6',
                                  TRATAMENTO == 'T6R2' ~'T6',
                                  TRATAMENTO == 'T6R3' ~'T6',
                                  TRATAMENTO == 'T6R4' ~'T6',
                                  TRATAMENTO == 'T7R1' ~'T7',
                                  TRATAMENTO == 'T7R2' ~'T7',
                                  TRATAMENTO == 'T7R3' ~'T7',
                                  TRATAMENTO == 'T7R4' ~'T7',
                                  TRATAMENTO == 'T8R1' ~'T8',
                                  TRATAMENTO == 'T8R2' ~'T8',
                                  TRATAMENTO == 'T8R3' ~'T8',
                                  TRATAMENTO == 'T8R4' ~'T8',
                                  TRATAMENTO == 'T9R1' ~'T9',
                                  TRATAMENTO == 'T9R2' ~'T9',
                                  TRATAMENTO == 'T9R3' ~'T9',
                                  TRATAMENTO == 'T9R4' ~'T9',
                                  TRATAMENTO == 'T10R1' ~'T10',
                                  TRATAMENTO == 'T10R2' ~'T10',
                                  TRATAMENTO == 'T10R3' ~'T10',
                                  TRATAMENTO == 'T10R4' ~'T10',
                                  TRATAMENTO == 'T11R1' ~'T11',
                                  TRATAMENTO == 'T11R2' ~'T11',
                                  TRATAMENTO == 'T11R3' ~'T11',
                                  TRATAMENTO == 'T11R4' ~'T11',
                                  TRATAMENTO == 'T12R1' ~'T12',
                                  TRATAMENTO == 'T12R2' ~'T12',
                                  TRATAMENTO == 'T12R3' ~'T12',
                                  TRATAMENTO == 'T12R4' ~'T12'))
}

### Separa as bases por doenças e remove outliers
D1 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(DIPLOCARPON.In = ifelse((abs(DIPLOCARPON.In - median(DIPLOCARPON.In)) > 2*sd(DIPLOCARPON.In)), mean(DIPLOCARPON.In), DIPLOCARPON.In)) %>%
  select(SEMANA, TRATAMENTO, DIPLOCARPON.In) 

D2 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(DIPLOCARPON.Sev = ifelse((abs(DIPLOCARPON.Sev - median(DIPLOCARPON.Sev)) > 2*sd(DIPLOCARPON.Sev)), mean(DIPLOCARPON.Sev), DIPLOCARPON.Sev)) %>%
  select(SEMANA, TRATAMENTO, DIPLOCARPON.Sev) 

D3 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(MICOSFAERELA.In = ifelse(!(abs(MICOSFAERELA.In - median(MICOSFAERELA.In)) > 2*sd(MICOSFAERELA.In)), mean(MICOSFAERELA.In), MICOSFAERELA.In)) %>%
  select(SEMANA, TRATAMENTO, MICOSFAERELA.In) 

D4 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(MICOSFAERELA.Sev = ifelse((abs(MICOSFAERELA.Sev - median(MICOSFAERELA.Sev)) > 2*sd(MICOSFAERELA.Sev)), mean(MICOSFAERELA.Sev), MICOSFAERELA.Sev)) %>%
  select(SEMANA, TRATAMENTO, MICOSFAERELA.Sev) 

D5 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(ANTRACNOSE.In = ifelse((abs(ANTRACNOSE.In - median(ANTRACNOSE.In)) > 2*sd(ANTRACNOSE.In)), mean(ANTRACNOSE.In), ANTRACNOSE.In)) %>%
  select(SEMANA, TRATAMENTO, ANTRACNOSE.In) 

D6 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(ANTRACNOSE.Sev = ifelse(!(abs(ANTRACNOSE.Sev - median(ANTRACNOSE.Sev)) > 2*sd(ANTRACNOSE.Sev)), mean(ANTRACNOSE.Sev), ANTRACNOSE.Sev)) %>%
  select(SEMANA, TRATAMENTO, ANTRACNOSE.Sev) 

D7 = df %>%
  group_by(TRATAMENTO) %>%
  mutate(BOTRYTIS.Inf = ifelse(!(abs(BOTRYTIS.Inf - median(BOTRYTIS.Inf)) > 2*sd(BOTRYTIS.Inf)), mean(BOTRYTIS.Inf), BOTRYTIS.Inf)) %>%
  select(SEMANA, TRATAMENTO, BOTRYTIS.Inf) 

### Renomeia os tratamentos
D1 <- my_fun(D1)
D2 <- my_fun(D2)
D3 <- my_fun(D3)
D4 <- my_fun(D4)
D5 <- my_fun(D5)
D6 <- my_fun(D6)
D7 <- my_fun(D7)

#################################################################
####################### DIPLOCARPON.In ##########################
#################################################################

### Plot Boxplot disease
D1 %>% 
  ggplot(aes(TRATAMENTO, DIPLOCARPON.In)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

###################################################
kruskal.test(DIPLOCARPON.In~TRATAMENTO,D1)

    Kruskal-Wallis rank sum test

data:  DIPLOCARPON.In by TRATAMENTO
Kruskal-Wallis chi-squared = 14.041, df = 5, p-value = 0.01535

pairwise.wilcox.test(D1$DIPLOCARPON.In, D1$TRATAMENTO)

    Pairwise comparisons using Wilcoxon rank sum test 

data:  D1$DIPLOCARPON.In and D1$TRATAMENTO 

    T1    T11   T12   T4    T5   
T11 0.354 -     -     -     -    
T12 0.049 1.000 -     -     -    
T4  1.000 1.000 0.428 -     -    
T5  0.039 1.000 1.000 0.400 -    
T6  0.684 1.000 1.000 1.000 1.000

P value adjustment method: holm 

##################################################

################ SUMMARY
saida <- describeBy(D1 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média    DesvP
T1  5.6725 4.636974
T11 3.3250 3.257162
T12 2.7000 2.823346
T4  4.2325 3.501383
T5  2.5050 2.153465
T6  3.8300 3.827646

#################################################################
####################### DIPLOCARPON.Sev ##########################
#################################################################

### Plot Boxplot disease
D2 %>% 
  ggplot(aes(TRATAMENTO, DIPLOCARPON.Sev)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

ANOVA<-aov(DIPLOCARPON.Sev~TRATAMENTO,D2)
res <- ANOVA$residuals
shapiro.test(res)

    Shapiro-Wilk normality test

data:  res
W = 0.81426, p-value = 3.158e-16

summary(ANOVA)

             Df Sum Sq Mean Sq F value Pr(>F)
TRATAMENTO    5   1.24  0.2485    0.92  0.468
Residuals   234  63.20  0.2701               

TukeyHSD(ANOVA)

  Tukey multiple comparisons of means
    95% family-wise confidence level

Fit: aov(formula = DIPLOCARPON.Sev ~ TRATAMENTO, data = D2)

$TRATAMENTO
           diff        lwr       upr     p adj
T11-T1   0.0725 -0.2614094 0.4064094 0.9891805
T12-T1   0.0350 -0.2989094 0.3689094 0.9996650
T4-T1   -0.0425 -0.3764094 0.2914094 0.9991367
T5-T1   -0.0625 -0.3964094 0.2714094 0.9945560
T6-T1   -0.1500 -0.4839094 0.1839094 0.7898663
T12-T11 -0.0375 -0.3714094 0.2964094 0.9995305
T4-T11  -0.1150 -0.4489094 0.2189094 0.9209167
T5-T11  -0.1350 -0.4689094 0.1989094 0.8544640
T6-T11  -0.2225 -0.5564094 0.1114094 0.3957297
T4-T12  -0.0775 -0.4114094 0.2564094 0.9853522
T5-T12  -0.0975 -0.4314094 0.2364094 0.9598962
T6-T12  -0.1850 -0.5189094 0.1489094 0.6046418
T5-T4   -0.0200 -0.3539094 0.3139094 0.9999789
T6-T4   -0.1075 -0.4414094 0.2264094 0.9397452
T6-T5   -0.0875 -0.4214094 0.2464094 0.9748612

################ SUMMARY
saida <- describeBy(D2 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média     DesvP
T1  0.8775 0.4054042
T11 0.9500 0.6774765
T12 0.9125 0.5979398
T4  0.8350 0.3641287
T5  0.8150 0.5484477
T6  0.7275 0.4540346

#################################################################
####################### MICOSFAERELA.In ##########################
#################################################################
### Plot Boxplot disease
D3 %>% 
  ggplot(aes(TRATAMENTO, MICOSFAERELA.In)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

kruskal.test(MICOSFAERELA.In~TRATAMENTO,D3)

    Kruskal-Wallis rank sum test

data:  MICOSFAERELA.In by TRATAMENTO
Kruskal-Wallis chi-squared = 54.403, df = 5, p-value = 1.732e-10

pairwise.wilcox.test(D3$MICOSFAERELA.In, D3$TRATAMENTO)

    Pairwise comparisons using Wilcoxon rank sum test 

data:  D3$MICOSFAERELA.In and D3$TRATAMENTO 

    T1      T11     T12    T4     T5    
T11 0.0080  -       -      -      -     
T12 1.0000  0.2311  -      -      -     
T4  0.0011  1.3e-07 0.0014 -      -     
T5  4.2e-05 1.1e-05 0.0187 1.0000 -     
T6  0.0232  5.7e-05 0.1694 1.0000 1.0000

P value adjustment method: holm 

saida <- describeBy(D3 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média     DesvP
T1  0.9075 0.5050451
T11 1.6050 1.8886910
T12 1.1575 1.0928409
T4  0.5400 0.8414883
T5  0.5325 0.7630531
T6  0.6300 0.7759857

#################################################################
####################### MICOSFAERELA.Sev ##########################
#################################################################
### Plot Boxplot disease
D4 %>% 
  ggplot(aes(TRATAMENTO, MICOSFAERELA.Sev)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

kruskal.test(MICOSFAERELA.Sev~TRATAMENTO,D4)

    Kruskal-Wallis rank sum test

data:  MICOSFAERELA.Sev by TRATAMENTO
Kruskal-Wallis chi-squared = 28.231, df = 5, p-value = 3.28e-05

pairwise.wilcox.test(D4$MICOSFAERELA.Sev, D4$TRATAMENTO)

    Pairwise comparisons using Wilcoxon rank sum test 

data:  D4$MICOSFAERELA.Sev and D4$TRATAMENTO 

    T1     T11    T12    T4     T5    
T11 1.0000 -      -      -      -     
T12 1.0000 0.6616 -      -      -     
T4  0.0140 0.0079 0.2617 -      -     
T5  0.0077 0.0051 0.1908 1.0000 -     
T6  0.0397 0.0239 0.6416 1.0000 1.0000

P value adjustment method: holm 

saida <- describeBy(D4 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média     DesvP
T1  0.4725 0.5392243
T11 0.6325 0.7353309
T12 0.3525 0.5218250
T4  0.1225 0.3100765
T5  0.0850 0.2348486
T6  0.1400 0.3326660

#################################################################
####################### ANTRACNOSE.In ##########################
#################################################################
### Plot Boxplot disease
D5 %>% 
  ggplot(aes(TRATAMENTO, ANTRACNOSE.In)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

kruskal.test(ANTRACNOSE.In~TRATAMENTO,D5)

    Kruskal-Wallis rank sum test

data:  ANTRACNOSE.In by TRATAMENTO
Kruskal-Wallis chi-squared = 17.353, df = 5, p-value = 0.003877

pairwise.wilcox.test(D5$ANTRACNOSE.In, D5$TRATAMENTO)

    Pairwise comparisons using Wilcoxon rank sum test 

data:  D5$ANTRACNOSE.In and D5$TRATAMENTO 

    T1    T11   T12   T4    T5   
T11 0.027 -     -     -     -    
T12 0.050 1.000 -     -     -    
T4  0.024 1.000 1.000 -     -    
T5  0.012 1.000 1.000 1.000 -    
T6  0.017 1.000 1.000 1.000 1.000

P value adjustment method: holm 

saida <- describeBy(D5 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média    DesvP
T1  4.4050 3.143895
T11 2.2850 1.948379
T12 2.4425 1.854155
T4  2.3225 2.336169
T5  2.1100 1.896123
T6  2.1650 1.866788

#################################################################
####################### ANTRACNOSE.Sev ##########################
#################################################################
### Plot Boxplot disease
D6 %>% 
  ggplot(aes(TRATAMENTO, ANTRACNOSE.Sev)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

kruskal.test(ANTRACNOSE.Sev~TRATAMENTO,D6)

    Kruskal-Wallis rank sum test

data:  ANTRACNOSE.Sev by TRATAMENTO
Kruskal-Wallis chi-squared = 27.756, df = 5, p-value = 4.062e-05

pairwise.wilcox.test(D6$ANTRACNOSE.Sev, D6$TRATAMENTO)

    Pairwise comparisons using Wilcoxon rank sum test 

data:  D6$ANTRACNOSE.Sev and D6$TRATAMENTO 

    T1      T11    T12    T4     T5    
T11 0.6004  -      -      -      -     
T12 0.0086  1.0000 -      -      -     
T4  1.3e-07 1.0000 1.0000 -      -     
T5  0.0114  1.0000 1.0000 1.0000 -     
T6  1.0e-05 1.0000 1.0000 1.0000 1.0000

P value adjustment method: holm 

saida <- describeBy(D6 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média     DesvP
T1  1.9175 0.6126393
T11 1.4550 0.7214480
T12 1.5525 0.6151871
T4  1.3525 0.3862957
T5  1.4675 0.8483430
T6  1.5100 0.7989737

#################################################################
####################### BOTRYTIS.Inf ##########################
#################################################################
### Plot Boxplot disease
D7 %>% 
  ggplot(aes(TRATAMENTO, BOTRYTIS.Inf)) +
  geom_jitter(size=3, width = .1)+
  geom_boxplot(size=1)

kruskal.test(BOTRYTIS.Inf~TRATAMENTO,D7)

    Kruskal-Wallis rank sum test

data:  BOTRYTIS.Inf by TRATAMENTO
Kruskal-Wallis chi-squared = 35.737, df = 5, p-value = 1.072e-06

pairwise.wilcox.test(D7$BOTRYTIS.Inf, D7$TRATAMENTO)

    Pairwise comparisons using Wilcoxon rank sum test 

data:  D7$BOTRYTIS.Inf and D7$TRATAMENTO 

    T1      T11     T12     T4      T5     
T11 1.00000 -       -       -       -      
T12 1.00000 1.00000 -       -       -      
T4  1.00000 1.00000 1.00000 -       -      
T5  1.00000 0.24717 0.47524 0.71947 -      
T6  0.00231 9.5e-06 0.00161 8.2e-06 0.00038

P value adjustment method: holm 

saida <- describeBy(D7 ~ TRATAMENTO, skew=FALSE,ranges=FALSE,mat=TRUE)
mat <- c(saida[13,5], saida[13,6], saida[14,5], saida[14,6],
         saida[15,5], saida[15,6], saida[16,5], saida[16,6],
         saida[17,5], saida[17,6], saida[18,5], saida[18,6])
rnames <- c(toString(saida[13,2]),toString(saida[14,2]),toString(saida[15,2]),
            toString(saida[16,2]),toString(saida[17,2]),toString(saida[18,2]))
cnames <- c("Média","DesvP")
mat <- matrix(mat,nrow=6, ncol=2,byrow=TRUE,dimnames=list(rnames,cnames))
print("Valores de Média e Desvio Padrão")

[1] "Valores de Média e Desvio Padrão"

print(mat)

     Média     DesvP
T1  1.7275 1.5207939
T11 1.6900 1.4840692
T12 1.4050 0.6898346
T4  1.8650 1.8259595
T5  1.3500 1.1502508
T6  0.9175 1.3157600