La idea es comparar 2 variedades de arroz y tres plántulas (2 semanas de edad, 3 semanas de edad y 4 semanas de edad).
Se realizó un experimento en un diseño completamente aleatorio con 3 repeticiones, mientras que las variedades se mantuvieron en la parcela principal y la edad de las plántulas en subparcelas.
La variable de respuesta se registró como tiempo en días desde el inicio hasta la madurez de las plantas de arroz.
Split Plot Analysis
se limpian los valores en en entorno global, terminal y cierra los plots
rm(list = ls(all = TRUE))
graphics.off()
shell("cls") #shell("cls") = Limpiar consola
Se sube la data para el estudio, el cual es una evaluación del estado de desarrollo de la planta respecto a su
library(readxl)
data <- read.csv(file = "data_split.csv",
header = TRUE)
data
## block var age heading
## 1 1 Super 2W 89
## 2 1 Super 3W 94
## 3 1 Super 4W 96
## 4 1 Shaheen 2W 75
## 5 1 Shaheen 3W 80
## 6 1 Shaheen 4W 84
## 7 2 Super 2W 85
## 8 2 Super 3W 99
## 9 2 Super 4W 101
## 10 2 Shaheen 2W 72
## 11 2 Shaheen 3W 74
## 12 2 Shaheen 4W 82
## 13 3 Super 2W 85
## 14 3 Super 3W 93
## 15 3 Super 4W 99
## 16 3 Shaheen 2W 67
## 17 3 Shaheen 3W 79
## 18 3 Shaheen 4W 79
## 19 4 Super 2W 81
## 20 4 Super 3W 93
## 21 4 Super 4W 101
## 22 4 Shaheen 2W 75
## 23 4 Shaheen 3W 81
## 24 4 Shaheen 4W 87
Se establece que el modelo experimental será el siguiente:
library(collapsibleTree)
collapsibleTree(data,hierarchy = c("var","age","heading"))
Observar Variables
Ver la estructura de las variables:
var y age como factores:
data$var = as.factor(data$var)
data$age = as.factor(data$age)
str(data)
## 'data.frame': 24 obs. of 4 variables:
## $ block : int 1 1 1 1 1 1 2 2 2 2 ...
## $ var : Factor w/ 2 levels "Shaheen","Super": 2 2 2 1 1 1 2 2 2 1 ...
## $ age : Factor w/ 3 levels "2W","3W","4W": 1 2 3 1 2 3 1 2 3 1 ...
## $ heading: int 89 94 96 75 80 84 85 99 101 72 ...
Análisis de ajuste del modelo de varianza
Al usar la librería agricolae podemos ajustar el modelo de analisis de varianza para el diseño de parcelas divididas asignando una función para el modelo con sp.plot para la variable de bloqueo(block) y las sub- parcelas
age, block, heading, var
library(agricolae)
## Warning: package 'agricolae' was built under R version 4.2.2
attach(data)
model <- sp.plot(block = block,
pplot = var,
splot = age,
Y = heading)
##
## ANALYSIS SPLIT PLOT: heading
## Class level information
##
## var : Super Shaheen
## age : 2W 3W 4W
## block : 1 2 3 4
##
## Number of observations: 24
##
## Analysis of Variance Table
##
## Response: heading
## Df Sum Sq Mean Sq F value Pr(>F)
## block 3 28.46 9.49 NaN NaN
## var 1 1365.04 1365.04 63.5314 0.00412 **
## Ea 3 64.46 21.49 NaN NaN
## age 2 641.33 320.67 44.5714 2.789e-06 ***
## var:age 2 16.33 8.17 1.1351 0.35359
## Eb 12 86.33 7.19 NaN NaN
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## cv(a) = 5.4 %, cv(b) = 3.1 %, Mean = 85.45833
El análisis de la tabla de varianza mostró que ambos efectos principales afectaron significativamente el tiempo que transcurre desde la salida hasta la madurez del arroz (var, age).
Si bien la interacción no fue significativa con respecto a los días desde el inicio hasta el vencimiento.
Sin embargo, procederé con la prueba de comparación media tanto para el tratamiento individual como para el efecto de interacción en caso de que obtenga el efecto significativo del término de interacción.
Prueba de separación media o LSD.test
Primero se establecen los errores por así decirlo para cada factor e interacción.
Edf_a = error del grado de libertad por la parcela principal
Edf_b = para el error del grado de libertad por la subparcela
# obteneoms primer error df
Edf_a <- model$gl.a
Edf_a
## [1] 3
# segundo error df
Edf_b <- model$gl.b
Edf_b
## [1] 12
# Primer error MS
EMS_a <- model$Ea
EMS_a
## [1] 21.48611
# segundo error MS
EMS_b <- model$Eb
EMS_b
## [1] 7.194444
Para LSD.test se requiere algunos argumentos especificos.
Var respuesta = heading (días totales en completar su desarrollo) en Y
El trt que va en X y el DFerror y MSerror con un nivel del 5%.
p.adj = método
out1 <- LSD.test(y = heading,
trt = var,
DFerror = Edf_a,
MSerror = EMS_a,
alpha = 0.05,
p.adj = "bonferroni",
group = TRUE,
console = TRUE)
##
## Study: heading ~ var
##
## LSD t Test for heading
## P value adjustment method: bonferroni
##
## Mean Square Error: 21.48611
##
## var, means and individual ( 95 %) CI
##
## heading std r LCL UCL Min Max
## Shaheen 77.91667 5.550730 12 73.65824 82.17510 67 87
## Super 93.00000 6.728501 12 88.74157 97.25843 81 101
##
## Alpha: 0.05 ; DF Error: 3
## Critical Value of t: 3.182446
##
## Minimum Significant Difference: 6.022327
##
## Treatments with the same letter are not significantly different.
##
## heading groups
## Super 93.00000 a
## Shaheen 77.91667 b
Es similar para la subparcela solo se cambian los valores de trt, DFerror y MSerror.
out2 <- LSD.test(y = heading,
trt = age,
DFerror = Edf_b,
MSerror = EMS_b,
alpha = 0.05,
p.adj = "bonferroni",
group = TRUE,
console = TRUE)
##
## Study: heading ~ age
##
## LSD t Test for heading
## P value adjustment method: bonferroni
##
## Mean Square Error: 7.194444
##
## age, means and individual ( 95 %) CI
##
## heading std r LCL UCL Min Max
## 2W 78.625 7.558108 8 76.55879 80.69121 67 89
## 3W 86.625 9.117291 8 84.55879 88.69121 74 99
## 4W 91.125 9.093758 8 89.05879 93.19121 79 101
##
## Alpha: 0.05 ; DF Error: 12
## Critical Value of t: 2.779473
##
## Minimum Significant Difference: 3.727616
##
## Treatments with the same letter are not significantly different.
##
## heading groups
## 4W 91.125 a
## 3W 86.625 b
## 2W 78.625 c
LSD_test en términos de interacción:
en el trt se concadenan los factores (age, var) concadenandolos
out3 <- LSD.test(y = heading,
trt = var:age,
DFerror = Edf_b,
MSerror = EMS_b,
alpha = 0.05,
p.adj = "bonferroni",
group = TRUE,
console = TRUE)
##
## Study: heading ~ var:age
##
## LSD t Test for heading
## P value adjustment method: bonferroni
##
## Mean Square Error: 7.194444
##
## var:age, means and individual ( 95 %) CI
##
## heading std r LCL UCL Min Max
## Shaheen:2W 72.25 3.774917 4 69.32794 75.17206 67 75
## Shaheen:3W 78.50 3.109126 4 75.57794 81.42206 74 81
## Shaheen:4W 83.00 3.366502 4 80.07794 85.92206 79 87
## Super:2W 85.00 3.265986 4 82.07794 87.92206 81 89
## Super:3W 94.75 2.872281 4 91.82794 97.67206 93 99
## Super:4W 99.25 2.362908 4 96.32794 102.17206 96 101
##
## Alpha: 0.05 ; DF Error: 12
## Critical Value of t: 3.648889
##
## Minimum Significant Difference: 6.920609
##
## Treatments with the same letter are not significantly different.
##
## heading groups
## Super:4W 99.25 a
## Super:3W 94.75 a
## Super:2W 85.00 b
## Shaheen:4W 83.00 b
## Shaheen:3W 78.50 bc
## Shaheen:2W 72.25 c
Graficamos los resultados del LSD.test
plot(out1,
xlab = "variedades",
ylab = "Madurez (días)",
las = 1,
variation = "IQR")
plot(out2,
xlab = "edad de plantulas",
ylab = "Madurez (días)",
las = 1,
variation = "IQR")
plot(out3,
xlab = "variedad:edad de plantulas",
ylab = "madurez (días)",
las = 1,
variation = "IQR")
# Main plot factor (varieties)
bar.err(out1$means,
variation = "SE",
ylim = c(0, 100),
names.arg = c("Shaheen rice","Super rice"))
#adicionanado titulo y XY labels
title(main = "Variedades vs madurez",
cex.main = 0.8,
xlab = "Variedades",
ylab = "madurez (días)")
# Sub-plot factor (seedling age)
bar.err(out2$means,
variation = "SE",
ylim = c(0, 100),
names.arg = c("2-semanas", "3-semanas", "4-semanas"))
#adicionanado titulo y XY labels
title(main = "Edad de plantula (semanas) vs madurez",
cex.main = 0.8,
xlab = "Edad de las plantulas",
ylab = "Madurez (días)")
# parcela principal , factor (variedades)
bar.err(out3$means,
variation = "SE",
ylim = c(0, 110),
names.arg = c("V1W1","V1W2", "V1W3", "V2W1", "V2W2", "V2W3"))
#adicionanado titulo y XY labels
title(main = "grafico de barra con error estándar",
cex.main = 0.8,
xlab = "variedad:edad",
ylab = "Madurez (días)")
Medidas Repetidas
boxplot(heading ~ var, data = data,
lwd = 2,
xlab = "Variedad",
ylab="Madurez (días)", font.lab = 2,
cex.lab = 2,
main = "Diagrama de caja")
boxplot(heading ~ age, data = data,
lwd = 2,
xlab = "edad (semanas)",
ylab="Madurez (días)", font.lab = 2,
cex.lab = 2,
main = "Diagrama de caja")
boxplot(heading ~ var:age, data = data,
lwd = 1,
xlab = "Variedad",
ylab="Madurez (días)", font.lab = 1,
cex.lab = 2,
main = "Diagrama de caja")