---
title: "randomized block"
author: "Nancy"
date: "08/01/2020"
output:
  html_document: default
  pdf_document: default
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = FALSE)
```
```{r libraries}
library(faraway)
library(ggplot2)
library(dplyr)
library(kableExtra)
``` 
Blocking is very widely used in agricultural experiments (even today). In these experients blocks are physical areas of land in the field. In the FLM example in Ch 17 (2nd ed) or Ch 16 (1st ed) there are eight varieties of oats that have been planted in five blocks. The response variable is the yield in grams per row. (See `?oatvar`.) Every variety is planted in each block. So when we compare the average yield from  variety 1, say to the average yield from variety 2, it won't matter if, for example, Block V gives low yields, because that will impact all 8 varieties equally.

The data:
```{r data, echo = TRUE}
data(oatvar);head(oatvar);dim(oatvar)
is.factor(oatvar$block); is.factor(oatvar$variety)
``` 

There are 8 varieties and 5 blocks for a total of 40 observations. In the data frame there is one line per variety/block combination, but the data is more easily viewed in a table: 

```{r table, echo = TRUE}
xtabs(yield ~ variety + block, data = oatvar)
kable(xtabs(yield ~ variety + block, data = oatvar))
```

(I don't know much about kable, but it creates latex code, so only helps if you knit to pdf. It has many available adjustments.)

```{r plots}
plot(yield ~ variety, data = oatvar)
plot(yield ~ block, data = oatvar)
# interaction plot doesn't have a data argument
interaction.plot(oatvar$variety, oatvar$block,oatvar$yield)
```

See FLM-2 Figure 17.1 for a `ggplot` version.

```{r anova}
twoway <- lm(yield ~ variety + block, data = oatvar)
anova(twoway)
```

I will let you check that if you fit `yield ~ block + variety` you get the same anova table. (Why?)

What about the estimated coefficients? 
```{r summary}
summary(twoway)
```

These estimates are not very informative. For one thing, they are all relative to the baseline value of 0 for variety 1, block I. Also, the block effects aren't really of interest; they are there to ensure that the variety comparisons are fair. It is of more interest to construct a table of means:

```{r means}
tapply(oatvar$yield, oatvar$variety, mean)
```

(There must be more elegant ways to see the table of means, along with estimated standard errors.)

In this model checking constant variance is can be done any of 3 ways: 

```{r residuals}
plot(twoway$fitted, twoway$residuals)
qqnorm(twoway$residuals)
```

You can also plot the residuals against `block` and `variety`; the default plot is a side-by-side boxplot, but I think points look nicer. 
