Diseno de Bloques completos al azar

Ejercicio de rendimiento de almidon de yuca

Una empresa de alimentos está estudiando el efecto de 4 temperaturas sobre el rendimiento del almidón de yuca (en gramos). Cada evaluación tarda aproximadamente una hora y media, por lo que solamente se pueden realizar 4 evaluaciones por día. La idea general consiste en realizar una repetición completa del experimento en cada uno de los niveles del factor de control (factor bloque) y se realiza la aleatorización de los tratamientos dentro de cada bloque.

Temperatura	Dia 1	Dia 2	Dia 3	Dia 4
30	14	14.1	14.5	14
35	13.9	13.8	14.2	14
40	14.1	14.2	14.4	13.9
45	13	13.1	13.5	13.2

Carga de datos

#install.packages("tidyverse")
library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.4     ✔ readr     2.1.5
## ✔ forcats   1.0.0     ✔ stringr   1.5.1
## ✔ ggplot2   3.5.1     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

Yuca_desordenada <- read_csv("data/Yuca_desordenada.csv")

## Rows: 4 Columns: 5
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## dbl (5): Temperatura, Dia 1, Dia 2, Dia 3, Dia 4
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

Estructurando los datos

yuca <- Yuca_desordenada %>% 
  pivot_longer(cols = `Dia 1`:`Dia 4`, 
               names_to = "Dia",
               values_to = "Rendimiento")
knitr::kable(yuca)

Temperatura	Dia	Rendimiento
30	Dia 1	14.0
30	Dia 2	14.1
30	Dia 3	14.5
30	Dia 4	14.0
35	Dia 1	13.9
35	Dia 2	13.8
35	Dia 3	14.2
35	Dia 4	14.0
40	Dia 1	14.1
40	Dia 2	14.2
40	Dia 3	14.4
40	Dia 4	13.9
45	Dia 1	13.0
45	Dia 2	13.1
45	Dia 3	13.5
45	Dia 4	13.2

Analisis exploratorio de datos

str(yuca) #estructura de los datos

## tibble [16 × 3] (S3: tbl_df/tbl/data.frame)
##  $ Temperatura: num [1:16] 30 30 30 30 35 35 35 35 40 40 ...
##  $ Dia        : chr [1:16] "Dia 1" "Dia 2" "Dia 3" "Dia 4" ...
##  $ Rendimiento: num [1:16] 14 14.1 14.5 14 13.9 13.8 14.2 14 14.1 14.2 ...

yuca$Temperatura <- as.factor(yuca$Temperatura)
yuca$Dia <- as.factor(yuca$Dia)

#calculado rento total
sum(yuca$Rendimiento)

## [1] 221.9

# Creacion de nueva columnas
yuca <- yuca %>% 
  mutate(rento_perc=(Rendimiento/sum(yuca$Rendimiento))*100)


# Resumen de los datos
str(yuca)

## tibble [16 × 4] (S3: tbl_df/tbl/data.frame)
##  $ Temperatura: Factor w/ 4 levels "30","35","40",..: 1 1 1 1 2 2 2 2 3 3 ...
##  $ Dia        : Factor w/ 4 levels "Dia 1","Dia 2",..: 1 2 3 4 1 2 3 4 1 2 ...
##  $ Rendimiento: num [1:16] 14 14.1 14.5 14 13.9 13.8 14.2 14 14.1 14.2 ...
##  $ rento_perc : num [1:16] 6.31 6.35 6.53 6.31 6.26 ...

summary(yuca)

##  Temperatura    Dia     Rendimiento      rento_perc   
##  30:4        Dia 1:4   Min.   :13.00   Min.   :5.858  
##  35:4        Dia 2:4   1st Qu.:13.72   1st Qu.:6.185  
##  40:4        Dia 3:4   Median :14.00   Median :6.309  
##  45:4        Dia 4:4   Mean   :13.87   Mean   :6.250  
##                        3rd Qu.:14.12   3rd Qu.:6.365  
##                        Max.   :14.50   Max.   :6.534

Analisis grafico - Temperatura

bxp_temp <- yuca %>% ggplot(aes(x = Temperatura, y = Rendimiento, fill = Temperatura)) +
  geom_boxplot() +
  labs(x = "Temperatura (°C)", y  = "Rendimiento (%)", title = "Boxplot - Temperatura")


ggsave(plot = bxp_temp, "bxp_temp.png", dpi = 300, units = "in", height = 6, width = 8)

bxp_temp

Analisis grafico - Dias

bxp_dias <- yuca %>% ggplot(aes(x = Dia, y = Rendimiento, fill = Dia)) +
  geom_boxplot() +
  labs(x = "Dias", y  = "Rendimiento (%)", title = "Boxplot - Dias")

ggsave(plot = bxp_dias, "bxp_dias.png", dpi = 300, units = "in", height = 6, width = 8)

bxp_dias

Analisis grafico - Dias:Temperatura

bxp_dias_temp <- yuca %>% ggplot(aes(x = Dia, y = Rendimiento, col = Dia)) +
  facet_wrap(~Temperatura) +
  geom_point(size = 4) +
  theme(axis.text.x = element_text(angle=45, vjust=1, hjust=1)) +
  labs(x = NULL, y  = "Rendimiento (%)", title = "Boxplot - Dias:Temperatura")

ggsave(plot = bxp_dias_temp, "bxp_dias_temp.png", dpi = 300, units = "in", height = 6, width = 8)

bxp_dias_temp

# interaccion entre dias y temperatura
yuca %>% ggplot(aes(x = Dia:Temperatura, y = Rendimiento, col = Dia:Temperatura)) +
  geom_jitter(show.legend = F) +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) +
  labs(x = "Dias", y  = "Rendimiento (%)", title = "Interaccion - Dias:Temperatura")

Analisis de correlacion entre el rendimiento y el rendimiento promedio

correlacion <- yuca %>% ggplot(aes(x = Rendimiento, y = rento_perc)) +
  geom_jitter(width = 0.25, height = 0.1, color = "red")+
  theme_bw() +
  labs(x = "Rendimiento de almidon", y = "Rendimiento porcentual", title = "Rento vs rento porcentual")

ggsave(plot = correlacion, "correlacion.png", dpi = 300, units = "in", height = 4, width = 6)

correlacion

Modelacion bajo un DBCA - ANOVA

modelo1 <- aov(Rendimiento ~ Temperatura + Dia, data = yuca)
summary(modelo1)

##             Df Sum Sq Mean Sq F value   Pr(>F)    
## Temperatura  3 2.4669  0.8223   73.55 1.19e-06 ***
## Dia          3 0.4269  0.1423   12.73  0.00137 ** 
## Residuals    9 0.1006  0.0112                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Prueba post - ANOVA (Duncan)

Vamos a hacer la prueba para temperatura

#install.packages("agricolae")
library(agricolae)
compara_temp <- duncan.test(modelo1, "Temperatura")

knitr::kable(compara_temp$groups %>% rownames_to_column("Temperatura"))

Temperatura	Rendimiento	groups
30	14.150	a
40	14.150	a
35	13.975	a
45	13.200	b

Vamos a hacer la prueba para los dias (BLOQUES)

compara_dias <- duncan.test(modelo1, "Dia")
knitr::kable(compara_dias$groups %>% rownames_to_column("Dia"))

Dia	Rendimiento	groups
Dia 3	14.150	a
Dia 2	13.800	b
Dia 4	13.775	b
Dia 1	13.750	b

Prueba post - ANOVA (TUKEY)