まず， 「日本晴」品種の 2000年度収量データ を題材にして，施肥量と栽植密度が収量に与える影響を解析してみる．
　いま，収量の値を処理組み合わせでラベルをつけ，X_ijk と表すことにする． ここで，i は栽植密度の水準で，j は施肥量，k はブロックである． 試験区 2 の収量値 307.8 に対するラベルは，X₁₁₁ となる．つまり，この試験区 2 では， 密度が水準 1 (疎植)で，施肥が水準 1 (施肥無)，ブロック1という処理を行った結果， 収量が 307.8 であったことになる．

施肥量水準に対する収量

　興味が持たれるのは，栽植密度の違いや施肥量の違いにより，収量が異なる かどうかということである．これを調べるためには，各処理ごとの平均の表 を見ればよい． 栽植密度の違いによって 収量はあまり違わないようにみえるが，施肥量が多くなると収量も 増加しているようにみえる．また，ブロック2の方が多少「でき」がよかったように みえる．ブロック全体ではまったく同じ処理組み合わせを行っているのだから ブロックの間には差が無いのが望ましい．しかしながら，地力の差や 学生の栽培管理の巧拙，測定誤差などのもろもろの影響により差がでてくる．
　次に，密度と施肥の処理組み合わせによる収量の違いもみてみよう． 表のデータや右の図 をみると，密度水準にかかわらず施肥量を増加させると収量も増加している 様子がわかる． 処理平均を計算した結果類推されることは，「日本晴の収量は，栽植密度によっては違わないが， 施肥量によって異なり，その異なり方は栽植密度水準の違いには影響されない．」 ということであるが，この類推がどの程度正しいかはこのままではわからない． これに答えるのが分散分析として知られる統計手法である．

Rのコマンド

• 「日本晴」データの読み込み（ データダウンロード）
  rice.dfr <- read.csv("ricecul.csv")
• 収量データの取り出し
  yield0 <- rice.dfr$gy[1:12]        #gy 列の12行目まで(2000年データ)を取り出す
• 水準のラベル化

  dens0 <- factor(rice.dfr$density[1:12])	#密度水準12行目までラベル化
  fert0 <- factor(rice.dfr$fert[1:12])	#施肥料水準
  blk0 <- factor(rice.dfr$rep[1:12])	#ブロック水準

• 処理水準ごとの平均

  tapply(yield0, dens0, mean)	#密度水準ごとの平均
  tapply(yield0, fert0, mean)	#施肥料水準ごとの平均
  tapply(yield0, blk0, mean)	#ブロック水準ごとの平均

  yield0 の値を dens0 の水準ラベルごとにグルーピングして，mean で平均をとる．
• 処理組み合わせごとの平均
  x <- tapply(yield0, dens0:fert0, mean)          #密度水準と施肥量水準の組み合わせの平均
  x
• 施肥量水準に対する収量のグラフ
  plot(1:3, x[1:3], type="b", lwd=2, cex=1.5, xaxt="n", xlab="fertility", ylab="yield",
           ylim=c(300, 600), pch=0, cex.lab=1.0, cex.axis=1.0, col="blue")    #密度水準1のグラフを青色で表示
  axis(1, 1:3, labels=c("fert1","fert2","fert3"), cex.axis=0.8)
  box(lwd=3)
  par(new=T)                                   #グラフの重ね書き
  plot(1:3, x[4:6], type="b",lwd=2, cex=1.5,
           axes=F, xlab="", ylab="" ,ylim=c(300,600), pch=2, col="red")    #密度水準2のグラフを赤色で表示
  legend(locator(1), legend=c("density1","density2"), pch=c(0,2),
  col=c("blue","red"), cex=0.8)          #凡例をクリックで場所を指定して表示