KKP - KKCS

Introduction

This document explores the Oats dataset — a classic agricultural experiment on oat yield. We will create tables, figures, and a simple model to demonstrate how Quarto integrates code, text, and output in a single file.

For more examples of what Quarto can produce, visit the Quarto Gallery.

Load Data

data("Oats", package = "nlme")
Oats <- as.data.frame(Oats)
str(Oats)

'data.frame':   72 obs. of  4 variables:
 $ Block  : Ord.factor w/ 6 levels "VI"<"V"<"III"<..: 6 6 6 6 6 6 6 6 6 6 ...
 $ Variety: Factor w/ 3 levels "Golden Rain",..: 3 3 3 3 1 1 1 1 2 2 ...
 $ nitro  : num  0 0.2 0.4 0.6 0 0.2 0.4 0.6 0 0.2 ...
 $ yield  : num  111 130 157 174 117 114 161 141 105 140 ...
 - attr(*, "formula")=Class 'formula'  language yield ~ nitro | Block
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 
 - attr(*, "labels")=List of 2
  ..$ y: chr "Yield"
  ..$ x: chr "Nitrogen concentration"
 - attr(*, "units")=List of 2
  ..$ y: chr "(bushels/acre)"
  ..$ x: chr "(cwt/acre)"
 - attr(*, "inner")=Class 'formula'  language ~Variety
  .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv> 

head(Oats)

  Block     Variety nitro yield
1     I     Victory   0.0   111
2     I     Victory   0.2   130
3     I     Victory   0.4   157
4     I     Victory   0.6   174
5     I Golden Rain   0.0   117
6     I Golden Rain   0.2   114

This dataset comes from a split-plot experiment on oat yield (Yates, 1935). It contains 3 varieties (Golden Rain, Marvellous, Victory), 4 nitrogen levels (0, 0.2, 0.4, 0.6 cwt/acre), and 6 blocks — 72 observations in total.

Inline Code

The dataset contains 72 observations. The overall mean yield was 104 bushels/acre, ranging from 53 to 174.

Visualization

Base R

boxplot(yield ~ nitro, data = Oats,
        xlab = "Nitrogen (cwt/acre)", ylab = "Yield (bu/acre)",
        main = "Oat Yield by Nitrogen Level",
        col = c("#E8B4B8", "#B8D4E3", "#93C5B8", "#F5D6A8"))

Figure 1: Oat yield by nitrogen application rate (base R)

ggplot2

ggplot2::ggplot(Oats, ggplot2::aes(x = factor(nitro), y = yield,
                                    fill = factor(nitro))) +
  ggplot2::geom_boxplot(alpha = 0.7, show.legend = FALSE) +
  ggplot2::geom_jitter(width = 0.1, alpha = 0.3) +
  ggplot2::labs(x = "Nitrogen (cwt/acre)", y = "Yield (bu/acre)") +
  ggplot2::theme_minimal(base_size = 14)

Figure 2: Oat yield by nitrogen application rate

As shown in Figure 2, yield increases with nitrogen application.

ggplot2::ggplot(Oats, ggplot2::aes(x = nitro, y = yield, color = Variety)) +
  ggplot2::geom_point(size = 2, alpha = 0.6) +
  ggplot2::geom_smooth(method = "lm", se = F) +
  ggplot2::labs(x = "Nitrogen (cwt/acre)", y = "Yield (bu/acre)") +
  ggplot2::theme_minimal(base_size = 14)

`geom_smooth()` using formula = 'y ~ x'

Figure 3: Nitrogen response by oat variety

Figure 3 shows that all three varieties respond positively to nitrogen, with similar slopes.

Simple Model

model <- lm(yield ~ nitro + Variety, data = Oats)
summary(model)

Call:
lm(formula = yield ~ nitro + Variety, data = Oats)

Residuals:
    Min      1Q  Median      3Q     Max 
-33.725 -15.133  -1.392  12.333  54.275 

Coefficients:
                  Estimate Std. Error t value Pr(>|t|)    
(Intercept)         82.400      5.483  15.029  < 2e-16 ***
nitro               73.667     11.192   6.582 7.96e-09 ***
VarietyMarvellous    5.292      6.130   0.863    0.391    
VarietyVictory      -6.875      6.130  -1.122    0.266    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 21.24 on 68 degrees of freedom
Multiple R-squared:  0.4102,    Adjusted R-squared:  0.3841 
F-statistic: 15.76 on 3 and 68 DF,  p-value: 6.965e-08

summary(model)$coefficients |>
  knitr::kable(digits = 3)

Table 1: Linear model coefficients

	Estimate	Std. Error	t value	Pr(>|t|)
(Intercept)	82.400	5.483	15.029	0.000
nitro	73.667	11.192	6.582	0.000
VarietyMarvellous	5.292	6.130	0.863	0.391
VarietyVictory	-6.875	6.130	-1.122	0.266

The model (Table 1) suggests that each additional unit of nitrogen increases yield by approximately 73.7 bushels/acre.

Try Other Formats!

Try changing format: html to format: docx in the YAML header and render again. Same content, different output!