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")