ANOVA PDP. Actividad
Luis Fernando Rosero Males
2024-01-13
Crear un vector de paquetes
pacotes <- c("agricolae", "performance", "ggplot2", "doebioresearch",
"outliers","car", "readxl","tidyverse",
"RColorBrewer","nortest", "stats",
"agridat", "performance")Script para instalar y cargar librerias y dependencias
if(sum(as.numeric(!pacotes %in% installed.packages())) != 0){
instalador <- pacotes[!pacotes %in% installed.packages()]
for(i in 1:length(instalador)) {
install.packages(instalador, dependencies = T)
break()}
sapply(pacotes, require, character = T)
} else {
sapply(pacotes, require, character = T)
}## Loading required package: agricolae## Loading required package: performance## Loading required package: ggplot2## Loading required package: doebioresearch## Warning: package 'doebioresearch' was built under R version 4.3.2## Loading required package: outliers## Loading required package: car## Loading required package: carData## Loading required package: readxl## Loading required package: tidyverse## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.3 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ lubridate 1.9.3 ✔ tibble 3.2.1
## ✔ purrr 1.0.2 ✔ tidyr 1.3.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ dplyr::recode() masks car::recode()
## ✖ purrr::some() masks car::some()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
## Loading required package: RColorBrewer
##
## Loading required package: nortest
##
## Loading required package: agridat## Warning: package 'agridat' was built under R version 4.3.2## agricolae performance ggplot2 doebioresearch outliers
## TRUE TRUE TRUE TRUE TRUE
## car readxl tidyverse RColorBrewer nortest
## TRUE TRUE TRUE TRUE TRUE
## stats agridat performance
## TRUE TRUE TRUEDiseño en Parcelas Divididas con dos factores ——–
df <- read_excel("C:/Users/PC/Desktop/Ing Agronomica/semestre 9/Software R/CLase 6/ACT.xlsx")
str(df)## tibble [24 × 4] (S3: tbl_df/tbl/data.frame)
## $ MetodoSiembra: chr [1:24] "S" "S" "S" "S" ...
## $ Variedades : chr [1:24] "BB" "BB" "BB" "BB" ...
## $ Bloques : num [1:24] 1 2 3 4 1 2 3 4 1 2 ...
## $ Rendimiento : num [1:24] 8.5 8.2 8 7.6 8 7 7.3 6.9 10 10.8 ...Convertir en factores
df$MetodoSiembra <- as.factor(df$MetodoSiembra)
df$Variedades <- as.factor(df$Variedades)
df$Bloques <- as.factor(df$Bloques)
str(df)## tibble [24 × 4] (S3: tbl_df/tbl/data.frame)
## $ MetodoSiembra: Factor w/ 2 levels "I","S": 2 2 2 2 2 2 2 2 2 2 ...
## $ Variedades : Factor w/ 3 levels "BB","BP","CR": 1 1 1 1 2 2 2 2 3 3 ...
## $ Bloques : Factor w/ 4 levels "1","2","3","4": 1 2 3 4 1 2 3 4 1 2 ...
## $ Rendimiento : num [1:24] 8.5 8.2 8 7.6 8 7 7.3 6.9 10 10.8 ...Valores atípicos
df%>% select_if(is.numeric)%>% outlier()
ANOVA DPD bifactorial ———-
Variable de respuesta Rendimiento
df[3]## # A tibble: 24 × 1
## Bloques
## <fct>
## 1 1
## 2 2
## 3 3
## 4 4
## 5 1
## 6 2
## 7 3
## 8 4
## 9 1
## 10 2
## # ℹ 14 more rowsdf## # A tibble: 24 × 4
## MetodoSiembra Variedades Bloques Rendimiento
## <fct> <fct> <fct> <dbl>
## 1 S BB 1 8.5
## 2 S BB 2 8.2
## 3 S BB 3 8
## 4 S BB 4 7.6
## 5 S BP 1 8
## 6 S BP 2 7
## 7 S BP 3 7.3
## 8 S BP 4 6.9
## 9 S CR 1 10
## 10 S CR 2 10.8
## # ℹ 14 more rowsRealizar el análisis de varianza (ANOVA)
?splitplot## starting httpd help server ... donemodelo <- splitplot(df[4], df$Bloques, df$MetodoSiembra, df$Variedades, 3)
modelo## $Rendimiento
## $Rendimiento[[1]]
## $Rendimiento[[1]][[1]]
## Analysis of Variance Table
##
## Response: dependent.var
## Df Sum Sq Mean Sq F value Pr(>F)
## block 3 0.175 0.0583 0.8607 0.5476318
## main.plot 1 18.027 18.0267 265.9672 0.0005016 ***
## Ea 3 0.203 0.0678
## sub.plot 2 34.308 17.1538 75.2174 1.626e-07 ***
## main.plot:sub.plot 2 7.776 3.8879 17.0481 0.0003112 ***
## Eb 12 2.737 0.2281
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## $Rendimiento[[1]][[2]]
## [1] "CV(a): 3.306 , CV(b) : 6.064"
##
## $Rendimiento[[1]][[3]]
## [1] "R Square 0.957"
##
## $Rendimiento[[1]][[4]]
##
## Shapiro-Wilk normality test
##
## data: model$residuals
## W = 0.95513, p-value = 0.3486
##
##
## $Rendimiento[[1]][[5]]
## [1] "Normality assumption is not violated"
##
## $Rendimiento[[1]][[6]]
## [1] "The means of one or more levels of main plot factor are not same, so go for multiple comparison test"
##
## $Rendimiento[[1]][[7]]
## $Rendimiento[[1]][[7]][[1]]
## MSerror Df Mean CV MSD
## 0.06777778 3 7.875 3.305926 0.3382432
##
## $Rendimiento[[1]][[7]][[2]]
## test name.t ntr StudentizedRange alpha
## Tukey main.plot 2 4.500659 0.05
##
## $Rendimiento[[1]][[7]][[3]]
## dependent.var groups
## S 8.741667 a
## I 7.008333 b
##
##
## $Rendimiento[[1]][[8]]
## [1] "The means of one or more levels of sub plot factor are not same, so go for multiple comparison test"
##
## $Rendimiento[[1]][[9]]
## $Rendimiento[[1]][[9]][[1]]
## MSerror Df Mean CV MSD
## 0.2280556 12 7.875 6.064148 0.6370213
##
## $Rendimiento[[1]][[9]][[2]]
## test name.t ntr StudentizedRange alpha
## Tukey sub.plot 3 3.772929 0.05
##
## $Rendimiento[[1]][[9]][[3]]
## dependent.var groups
## CR 9.5500 a
## BP 7.2375 b
## BB 6.8375 b
##
##
## $Rendimiento[[1]][[10]]
## [1] "The means of one or more levels of interaction are not same, so go for multiple comparison test"
##
## $Rendimiento[[1]][[11]]
## $Rendimiento[[1]][[11]][[1]]
## MSerror Df Mean CV MSD
## 0.2280556 12 7.875 6.064148 1.13424
##
## $Rendimiento[[1]][[11]][[2]]
## test name.t ntr StudentizedRange alpha
## Tukey main.plot:sub.plot 6 4.750231 0.05
##
## $Rendimiento[[1]][[11]][[3]]
## dependent.var groups
## S:CR 10.850 a
## I:CR 8.250 b
## S:BB 8.075 b
## S:BP 7.300 b
## I:BP 7.175 b
## I:BB 5.600 crendimiento_data <- as.data.frame(modelo$Rendimiento[[1]][[11]][[3]])
rendimiento_data <- rownames_to_column(rendimiento_data, var = "Interaction")
names(rendimiento_data) <- c("Interacción", "Rendimiento", "Grupos") ANOVA DPD —————
Crear el gráfico
library(ggplot2)
ggplot(df, aes(x = MetodoSiembra, y = Rendimiento, fill = Variedades)) +
geom_bar(stat = "summary", fun = "mean", position = position_dodge(width = 0.8), width = 0.7) +
labs(title = "Rendimiento por Método de siembra y Variedad",
x = "Métodos de Siembra",
y = "Rendimiento") +
scale_fill_manual(values = c("BP" = "blue", "CR" = "green", "BB" = "yellow"))