两因素方差分析

1 两因素方差分析

1. 引子
   
   • 考虑如下两个变量的关系, 不同的种子和不同的肥料之间的关系
   • 设有三种不同的种子A,B,C, 肥料有四种 1,2,3,4
   • 将农田分为很多小块, 分别采用不同的种子和肥料, 考虑亩产量
   • 农作物设为小麦
   • 计算每单位农田产小麦的颗粒数还是重量？
2. 例- 火箭推进器和燃料问题
   

   一火箭采用4种燃料, 三种推进器做射程试验, 每种燃料与每种推进器的组合各发射两次, 得射程如下表所示

   燃料 / 推进器B1B2B3A158.256.265.3 52.641.260.8A249.154.151.6 42.850.548.4A360.170.939.2 58.373.240.7A475.858.248.7 71.551.041.4
3. 模型分析
   
   • 这里试验指标是射程, 推进器和燃料是因素, 分别有3, 4个水平
   • 这是一个双因素的试验
   • 试验的目的在于考察在各种因素的各个水平下射程有误显著的差别
4. 理论模型
   
   • 设有两个因素 A,B 作用于试验的指标,
   • 因素A有 r 个水平 A1,⋯,Ar ,
   • 因素 B 有 s 个水平, B1,⋯,Bs .
   • 现对因素 A,B 的每对组合 (Ai,Bj) 进 行 t 次试验,所得的结果记为Xijk,i=1,⋯,r,j=1⋯,s,k=1,⋯,t
   • 设 Xijk∼N(μij,σ2) 各 Xijk 相互独立
   • σ2 为未知方差
5. 理论分析
   
   • 设总体的全部平均值为 μ , 则 μ 的估计为
     X¯...=1rst∑i=1r∑j=1s∑k=1rXijk
   • 设因素 A 取第 i 个水平 Ai 时, 数据的均值为 μi⋅ , 则 μi⋅ 的估计为
     X¯i..=1rs∑j=1s∑k=1rXijk
   • 因素A 的第 i 个水平 μi⋅ 相对于总体均值 μ 的提升(可能为负值) 为 αi=μi⋅−μ 的估计为 Xi⋅⋅−X¯⋯
6. 理论分析续
   
   • 欲检验 因素 A 的各个水平之间是否有显著差异, 可以考虑
     SSA=∑i=1rst(Xi⋅⋅−X¯⋯)2

     的大小, 作为度量标准

   • 另设 因素B取第 j 个水平时总体的均值为 μ⋅j , 提升为 βj=μ⋅j−μ , 估计值为 X¯⋅j⋅=1rt∑i=1r∑k=1tXijk
   • 欲检验 因素 B 的各个水平之间是否有显著差异, 可以考虑
     SSB=∑j=1srt(X⋅j⋅−X¯⋯)2

     的大小, 作为度量标准

7. 理论分析之误差平方和
   
   • 记误差平方和
     SSE=∑i=1r∑j=1s∑k=1t(Xijk−X¯ij⋅)2
   • 其中 X¯ij⋅=1t∑k=1tXijk
   • 可以证明 SSE/(rs(t−1)) 为方差 σ2 的无偏估计
8. 检验统计量的构造
   
   • 利用单因素方差分析的思想, 可以考虑
     FA=SSA/(r−1)SSE/(rs(t−1))

     作为检验统计量

   • 可以证明当 αi=μi⋅−μ=0,∀i=1,⋯,r 时, 有
     FA∼F(r−1,rs(t−1))
   • 关于因素B的效应可以类似处理
9. 因素A和因素B之间有无交互作用
   
   • 设两者之间没有关系, 此时 μij 可以看做是在均值 μ 的基 础上有两次提升即 μ+(μi⋅−μ)+(μ⋅j−μ)
   • 两者的差异 μij−(μi⋅+μ⋅j−μ) 记为γij 称为因素 A的第 i 个水平和 因 素B 的 第j 个水平的交互效应
   • 在实际数据中, 往往两者的估计值并不安全相等, 就是说 X¯ij⋅ 和 X¯i⋅⋅+X¯⋅j⋅−X¯⋯的值有差异
10. 检验交互作用统计量的构造
    
    • 考虑其差值的平方和作为两者是否有交互作用的度量，记
      SSAB=∑i=1r∑j=1st(X¯ij⋅−X¯i⋅⋅+X¯⋅j⋅−X¯⋯)2
    • 在假设 γij≡0, i=1,⋯,r, j=1,⋯,s 成立的条件下, 有
      SSABσ2∼χ2((r−1)(s−1))
    • 其中 (r−1)(s−1)=rst−1−(r−1)−(s−1)−[rs(t−1)]
    • 将 σ2 用估计值 SSErs(t−1) 代替，并将分子除以相应的自由度，得
      SSAB/(r−1)(s−1)SSE/rs(t−1)∼F((r−1)(s−1),rs(t−1))

2 两因素方差分析的R实现

1. 例
   

   四种不同的机床, 操作员有三位, 轮换进行了三件元件的制作, 制作时间共36个, 数据如下, 请问机器和操作员不同对元件制作有无显著的影响，另外机器和操作员之间有无交互效应

   ## 数据集的生成Y<-c(15,15,17,19,19,16,16,18,21,17,17,17,15,15,15,19,22,22,15,17,16,18,17,16,18,18,18,18,20,22,15,16,17,17,17,17)A<-gl(4,9,36)B<-gl(3,3,36)jiqi<-data.frame(Y,A,B)write.table(jiqi,file="data/jiqi.csv")

   ## 机器操作员数据的结构，输出前6行数据  jiqi<-read.table("data/jiqi.csv")  head(jiqi)

      Y A B1 15 1 12 15 1 13 17 1 14 19 1 25 19 1 26 16 1 2

   data/jiqi.csv

2. 两因素方差分析的代码
   

   jiqi$A<-as.factor(jiqi$A)jiqi$B<-as.factor(jiqi$B)

   jiqi.aov<-aov(Y~A+B+A:B,data=jiqi)summary(jiqi.aov)

               Df Sum Sq Mean Sq F value   Pr(>F)A            3   2.75   0.917   0.532 0.664528B            2  27.17  13.583   7.887 0.002330 **A:B          6  73.50  12.250   7.113 0.000192 ***Residuals   24  41.33   1.722---Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

   • 说明 Y~A+B+A:B 为模型，表示考虑主效应 A ,B 和交互效应 A:B 对响应变量 Y 的影响