Example 8.1 from Experimental Design and Analysis for Tree Improvement
Source:R/Exam8.1.R
Exam8.1.RdExam8.1 presents the Diameter at breast height (Dbh) of 60 SeedLots under layout of row column design with 6 rows and 10 columns in 18 countries and 59 provinces of 18 selected countries.
References
E.R. Williams, C.E. Harwood and A.C. Matheson (2023). Experimental Design and Analysis for Tree Improvement. CSIRO Publishing (https://www.publish.csiro.au/book/3145/).
Author
Muhammad Yaseen (myaseen208@gmail.com)
Sami Ullah (samiullahuos@gmail.com)
Examples
library(car)
library(dae)
library(dplyr)
library(emmeans)
library(ggplot2)
library(lmerTest)
library(magrittr)
library(predictmeans)
data(DataExam8.1)
# Pg. 141
fm8.4 <-
aov(
formula = dbh ~ inoc + Error(repl/inoc) +
inoc*country*prov
, data = DataExam8.1
)
#> Warning: Error() model is singular
# Pg. 150
summary(fm8.4)
#>
#> Error: repl
#> Df Sum Sq Mean Sq F value Pr(>F)
#> inoc 1 11.542 11.542 11.46 0.0773 .
#> Residuals 2 2.014 1.007
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#>
#> Error: Within
#> Df Sum Sq Mean Sq F value Pr(>F)
#> country 17 54.62 3.213 5.305 0.0000000159 ***
#> prov 41 18.61 0.454 0.749 0.854
#> inoc:country 17 10.07 0.592 0.978 0.487
#> inoc:prov 41 21.46 0.523 0.864 0.698
#> Residuals 116 70.26 0.606
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
# Pg. 150
model.tables(x = fm8.4, type = "means")
#> Tables of means
#> Grand mean
#>
#> 3.40411
#>
#> inoc
#> 7 weeks 1 week
#> 3.625 3.183
#> rep 118.000 118.000
#>
#> country
#> India Vietnam Egypt Kenya Fiji Thailand Malaysia Philippines Australia
#> 3.575 3.276 2.498 3.491 2.612 3.841 4.031 3.612 2.631
#> rep 24.000 20.000 12.000 32.000 12.000 16.000 36.000 12.000 16.000
#> PNG Solomon Is. Mauritius Sri Lanka Guam China Puerto Rico Vanuatu Benin
#> 3.65 3.699 3.122 3.243 2.342 3.686 3.345 2.762 3.342
#> rep 4.00 8.000 4.000 12.000 4.000 12.000 4.000 4.000 4.000
#>
#> prov
#> 1 2 3 4 5 6 7 8 10 11 12 13 14
#> 2.623 4.013 3.71 3.27 3.404 3.093 3.701 3.541 3.371 3.301 3.18 3.37 3.404
#> rep 4.000 4.000 4.00 4.00 4.000 4.000 4.000 4.000 4.000 4.000 4.00 4.00 4.000
#> 15 16 17 18 19 20 21 22 23 24 25 26 27
#> 3.595 3.43 3.275 3.085 3.66 3.382 3.235 3.46 3.08 3.555 3.918 3.648 3.43
#> rep 4.000 4.00 4.000 4.000 4.00 4.000 4.000 4.00 4.00 4.000 4.000 4.000 4.00
#> 28 29 30 31 32 33 34 35 36 37 38 39 40
#> 2.905 3.708 3.196 3.761 3.416 3.178 2.958 3.636 3.376 3.404 3.252 3.15 3.81
#> rep 4.000 4.000 4.000 4.000 4.000 4.000 4.000 4.000 4.000 4.000 4.000 4.00 4.00
#> 41 42 45 46 47 48 50 51 52 53 54 55 56
#> 3.195 3.613 3.518 2.76 3.733 3.605 3.404 3.685 3.235 3.755 3.605 2.74 3.662
#> rep 4.000 4.000 4.000 4.00 4.000 4.000 4.000 4.000 4.000 4.000 4.000 4.00 4.000
#> 57 58 59 60 61 62 63
#> 3.408 3.404 3.404 3.528 3.178 3.506 3.418
#> rep 4.000 4.000 4.000 4.000 4.000 4.000 4.000
#>
#> inoc:country
#> country
#> inoc India Vietnam Egypt Kenya Fiji Thailand Malaysia Philippines
#> 7 weeks 3.672 3.443 2.747 3.609 2.955 3.611 4.502 3.558
#> rep 12.000 10.000 6.000 16.000 6.000 8.000 18.000 6.000
#> 1 week 3.477 3.110 2.250 3.373 2.268 4.071 3.559 3.665
#> rep 12.000 10.000 6.000 16.000 6.000 8.000 18.000 6.000
#> country
#> inoc Australia PNG Solomon Is. Mauritius Sri Lanka Guam China
#> 7 weeks 2.959 3.850 4.200 3.390 3.695 2.245 4.030
#> rep 8.000 2.000 4.000 2.000 6.000 2.000 6.000
#> 1 week 2.304 3.450 3.197 2.855 2.792 2.440 3.342
#> rep 8.000 2.000 4.000 2.000 6.000 2.000 6.000
#> country
#> inoc Puerto Rico Vanuatu Benin
#> 7 weeks 3.540 2.720 3.870
#> rep 2.000 2.000 2.000
#> 1 week 3.150 2.805 2.815
#> rep 2.000 2.000 2.000
#>
#> inoc:prov
#> prov
#> inoc 1 2 3 4 5 6 7 8 10 11 12
#> 7 weeks 2.427 4.682 3.757 3.637 3.625 3.540 4.100 3.774 3.559 3.544 3.135
#> rep 2.000 2.819 3.344 3.664 2.904 3.183 2.646 3.301 3.308 3.183 3.058
#> 1 week 2.819 3.344 3.664 2.904 3.183 2.646 3.301 3.308 3.183 3.058 3.225
#> rep 2.427 4.682 3.757 3.637 3.625 3.540 4.100 3.774 3.559 3.544 3.135
#> prov
#> inoc 13 14 15 16 17 18 19 20 21 22 23
#> 7 weeks 3.845 3.625 4.104 4.119 3.604 3.114 3.509 3.304 3.591 3.801 3.276
#> rep 3.225 2.895 3.183 3.085 2.740 2.945 3.055 3.810 3.460 2.880 3.120
#> 1 week 2.895 3.183 3.085 2.740 2.945 3.055 3.810 3.460 2.880 3.120 2.885
#> rep 3.845 3.625 4.104 4.119 3.604 3.114 3.509 3.304 3.591 3.801 3.276
#> prov
#> inoc 24 25 26 27 28 29 30 31 32 33 34
#> 7 weeks 3.021 3.976 4.286 4.186 2.866 3.663 2.988 3.793 4.358 3.343 3.468
#> rep 2.885 4.090 3.860 3.010 2.675 2.945 3.754 3.404 3.729 2.474 3.014
#> 1 week 4.090 3.860 3.010 2.675 2.945 3.754 3.404 3.729 2.474 3.014 2.449
#> rep 3.021 3.976 4.286 4.186 2.866 3.663 2.988 3.793 4.358 3.343 3.468
#> prov
#> inoc 35 36 37 38 39 40 41 42 45 46 47
#> 7 weeks 3.848 3.688 3.625 3.772 3.132 3.972 3.545 3.705 3.439 2.599 4.119
#> rep 2.449 3.424 3.064 3.183 2.733 3.168 3.648 2.845 3.520 3.597 2.922
#> 1 week 3.424 3.064 3.183 2.733 3.168 3.648 2.845 3.520 3.597 2.922 3.347
#> rep 3.848 3.688 3.625 3.772 3.132 3.972 3.545 3.705 3.439 2.599 4.119
#> prov
#> inoc 48 50 51 52 53 54 55 56 57 58 59
#> 7 weeks 4.344 3.625 4.137 4.152 4.047 3.167 2.622 3.895 3.478 3.625 3.625
#> rep 3.347 2.867 3.183 3.233 2.318 3.463 4.043 2.858 3.430 3.339 3.183
#> 1 week 2.867 3.183 3.233 2.318 3.463 4.043 2.858 3.430 3.339 3.183 3.183
#> rep 4.344 3.625 4.137 4.152 4.047 3.167 2.622 3.895 3.478 3.625 3.625
#> prov
#> inoc 60 61 62 63
#> 7 weeks 3.685 3.630 3.560 3.235
#> rep 3.183 3.371 2.726 3.451
#> 1 week 3.371 2.726 3.451 3.601
#> rep 3.685 3.630 3.560 3.235
RESFit <-
data.frame(
fittedvalue = fitted.aovlist(fm8.4)
, residualvalue = proj(fm8.4)$Within[,"Residuals"]
)
ggplot(
data = RESFit
, mapping = aes(x = fittedvalue, y = residualvalue)
) +
geom_point(size = 2) +
labs(
x = "Residuals vs Fitted Values"
, y = ""
) +
theme_bw()
# Pg. 153
fm8.6 <-
aov(
formula = terms(
dbh ~ inoc + repl + col +
repl:row + repl:col +
prov + inoc:prov
, keep.order = TRUE
)
, data = DataExam8.1
)
summary(fm8.6)
#> Df Sum Sq Mean Sq F value Pr(>F)
#> inoc 1 11.54 11.542 48.054 0.00000000327 ***
#> repl 2 2.01 1.007 4.193 0.019746 *
#> col 9 65.24 7.249 30.182 < 2e-16 ***
#> repl:row 20 16.59 0.830 3.454 0.000105 ***
#> repl:col 27 16.41 0.608 2.530 0.001443 **
#> prov 58 53.89 0.929 3.869 0.00000026687 ***
#> inoc:prov 58 8.47 0.146 0.608 0.970544
#> Residuals 60 14.41 0.240
#> ---
#> Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1