4.3 格拉科-拉丁方设计
考虑一个 p × p p \times p p × p 拉丁方,并在其上叠加第二个 p × p p \times p p × p 拉丁方,后者中的处理用希腊字母表示。若叠加后两个方具有如下性质:每个希腊字母与每个拉丁字母都恰好一同出现一次,则称这两个拉丁方正交 (orthogonal),所得设计称为格拉科-拉丁方 (Graeco-Latin square)。表 4.18 给出了一个 4 × 4 4 \times 4 4 × 4 格拉科-拉丁方的例子。
格拉科-拉丁方设计可以系统地控制三个外部变异来源,即在三个方向上区组化。该设计允许在仅 p 2 p^{2} p 2 次运行中考察四个因子(行、列、拉丁字母和希腊字母),每个因子取 p p p 个水平。除 p = 6 p = 6 p = 6 以外,对所有 p ≥ 3 p \geq 3 p ≥ 3 都存在格拉科-拉丁方。
格拉科-拉丁方设计的统计模型为
y i j k l = μ + θ i + τ j + ω k + Ψ l + ε i j k l { i = 1 , 2 , … , p j = 1 , 2 , … , p k = 1 , 2 , … , p l = 1 , 2 , … , p (4.29) y_{ijkl} = \mu + \theta_{i} + \tau_{j} + \omega_{k} + \Psi_{l} + \varepsilon_{ijkl} \quad \left\{ \begin{array}{l} i = 1, 2, \ldots, p \\ j = 1, 2, \ldots, p \\ k = 1, 2, \ldots, p \\ l = 1, 2, \ldots, p \end{array} \right. \tag{4.29} y ijk l = μ + θ i + τ j + ω k + Ψ l + ε ijk l ⎩ ⎨ ⎧ i = 1 , 2 , … , p j = 1 , 2 , … , p k = 1 , 2 , … , p l = 1 , 2 , … , p ( 4.29 ) 其中 y i j k l y_{ijkl} y ijk l 是第 i i i 行、第 l l l 列中拉丁字母 j j j 、希腊字母 k k k 对应的观测值,θ i \theta_{i} θ i 是第 i i i 个行效应,τ j \tau_{j} τ j 是拉丁字母处理 j j j 的效应,ω k \omega_{k} ω k 是希腊字母处理 k k k 的效应,Ψ l \Psi_{l} Ψ l 是第 l l l 个列效应,ε i j k l \varepsilon_{ijkl} ε ijk l 是 N I D ( 0 , σ 2 ) \mathrm{NID}(0,\sigma^{2}) NID ( 0 , σ 2 ) 随机误差成分。完全确定一个观测值只需四个下标中的两个。
表 4.18 4 × 4 4 \times 4 4 × 4 格拉科-拉丁方设计
行 列 1 列 2 列 3 列 4 1 Aα Bβ Cγ Dδ 2 Bδ Aγ Dβ Cα 3 Cβ Dα Aδ Bγ 4 Dγ Cδ Bα Aβ
表 4.19 格拉科-拉丁方设计的方差分析
变异来源 平方和 自由度 拉丁字母处理 S S L = 1 p ∑ j = 1 p y . j . . 2 − y . . . . 2 N SS_{L} = \frac{1}{p} \sum_{j=1}^{p} y_{.j..}^{2} - \frac{y^{2}_{....}}{N} S S L = p 1 ∑ j = 1 p y . j .. 2 − N y .... 2 p − 1 p-1 p − 1 希腊字母处理 S S G = 1 p ∑ k = 1 p y . . k . 2 − y . . . . 2 N SS_{G} = \frac{1}{p} \sum_{k=1}^{p} y_{..k.}^{2} - \frac{y^{2}_{....}}{N} S S G = p 1 ∑ k = 1 p y .. k . 2 − N y .... 2 p − 1 p-1 p − 1 行 S S 行 = 1 p ∑ i = 1 p y i . . . 2 − y . . . . 2 N SS_{\text{行}} = \frac{1}{p} \sum_{i=1}^{p} y_{i...}^{2} - \frac{y^{2}_{....}}{N} S S 行 = p 1 ∑ i = 1 p y i ... 2 − N y .... 2 p − 1 p-1 p − 1 列 S S 列 = 1 p ∑ l = 1 p y . . . l 2 − y . . . . 2 N SS_{\text{列}} = \frac{1}{p} \sum_{l=1}^{p} y_{...l}^{2} - \frac{y^{2}_{....}}{N} S S 列 = p 1 ∑ l = 1 p y ... l 2 − N y .... 2 p − 1 p-1 p − 1 误差 S S E SS_{E} S S E (相减得到)( p − 3 ) ( p − 1 ) (p-3)(p-1) ( p − 3 ) ( p − 1 ) 总计 S S T = ∑ i ∑ j ∑ k ∑ l y i j k l 2 − y . . . . 2 N SS_{T} = \sum_{i} \sum_{j} \sum_{k} \sum_{l} y_{ijkl}^{2} - \frac{y^{2}_{....}}{N} S S T = ∑ i ∑ j ∑ k ∑ l y ijk l 2 − N y .... 2 p 2 − 1 p^{2}-1 p 2 − 1
其方差分析与拉丁方非常相似。由于希腊字母在每行每列中都恰好出现一次,并且与每个拉丁字母都恰好一同出现一次,希腊字母所代表的因子与行、列以及拉丁字母处理都正交。因此可以由希腊字母总和计算希腊字母因子对应的平方和,而试验误差又因这一数量而进一步减小。计算细节见表 4.19。行、列、拉丁字母处理和希腊字母处理均值相等的原假设,通过把相应的均方除以误差均方来检验。拒绝域是 F p − 1 , ( p − 3 ) ( p − 1 ) F_{p-1,(p-3)(p-1)} F p − 1 , ( p − 3 ) ( p − 1 ) 分布的上尾分位点。
例 4.3 ¶ 假设在例 4.2 的火箭推进剂试验中,还有一个因子——测试组件——可能很重要。设有五个测试组件,用希腊字母 α , β , γ , δ \alpha,\beta,\gamma,\delta α , β , γ , δ 和 ε \varepsilon ε 表示。所得的 5 × 5 5\times5 5 × 5 格拉科-拉丁方设计见表 4.20。
注意原材料批次(行)、操作人员(列)和配方(拉丁字母)的总和与例 4.2 中完全相同,因此有
S S 批次 = 68.00 , S S 操作人员 = 150.00 , 以及 S S 配方 = 330.00 \begin{array}{r l} SS_{\text{批次}} & = 68.00, \quad SS_{\text{操作人员}} = 150.00, \\ & \text{以及} \quad SS_{\text{配方}} = 330.00 \end{array} S S 批次 = 68.00 , S S 操作人员 = 150.00 , 以及 S S 配方 = 330.00 测试组件(希腊字母)的总和为
希腊字母 测试组件总和 α \alpha α y . . 1. = 10 y_{..1.} = 10 y ..1. = 10 β \beta β y . . 2. = − 6 y_{..2.} = -6 y ..2. = − 6 γ \gamma γ y . . 3. = − 3 y_{..3.} = -3 y ..3. = − 3 δ \delta δ y . . 4. = − 4 y_{..4.} = -4 y ..4. = − 4 ε \varepsilon ε y . . 5. = 13 y_{..5.} = 13 y ..5. = 13
于是测试组件对应的平方和为
S S 组件 = 1 p ∑ k = 1 p y . . k . 2 − y . . . . 2 N = 1 5 [ 1 0 2 + ( − 6 ) 2 + ( − 3 ) 2 + ( − 4 ) 2 + 1 3 2 ] − ( 10 ) 2 25 = 62.00 \begin{array}{r l} SS_{\text{组件}} & = \frac{1}{p} \sum_{k=1}^{p} y_{..k.}^{2} - \frac{y_{....}^{2}}{N} \\ & = \frac{1}{5} \left[ 10^{2} + (-6)^{2} + (-3)^{2} + (-4)^{2} + 13^{2} \right] - \frac{(10)^{2}}{25} = 62.00 \end{array} S S 组件 = p 1 ∑ k = 1 p y .. k . 2 − N y .... 2 = 5 1 [ 1 0 2 + ( − 6 ) 2 + ( − 3 ) 2 + ( − 4 ) 2 + 1 3 2 ] − 25 ( 10 ) 2 = 62.00 完整的方差分析汇总于表 4.21。各配方在 1% 水平上显著不同。比较表 4.21 与表 4.12 可以看出,剔除测试组件引起的变异后试验误差减小了。然而,在减小试验误差的同时,误差自由度也从 12(例 4.2 的拉丁方设计)减少到 8。因此我们对误差的估计自由度更少,检验的灵敏度可能降低。
表 4.20 火箭推进剂问题的格拉科-拉丁方设计
原材料批次 操作人员 1 操作人员 2 操作人员 3 操作人员 4 操作人员 5 y i . . . y_{i...} y i ... 1 A α = − 1 A\alpha = -1 A α = − 1 B γ = − 5 B\gamma = -5 B γ = − 5 C ε = − 6 C\varepsilon = -6 Cε = − 6 D β = − 1 D\beta = -1 D β = − 1 E δ = − 1 E\delta = -1 E δ = − 1 -14 2 B β = − 8 B\beta = -8 Bβ = − 8 C δ = − 1 C\delta = -1 C δ = − 1 D α = 5 D\alpha = 5 D α = 5 E γ = 2 E\gamma = 2 E γ = 2 A ε = 11 A\varepsilon = 11 A ε = 11 9 3 C γ = − 7 C\gamma = -7 C γ = − 7 D ε = 13 D\varepsilon = 13 D ε = 13 E β = 1 E\beta = 1 Eβ = 1 A δ = 2 A\delta = 2 A δ = 2 B α = − 4 B\alpha = -4 B α = − 4 5 4 D δ = 1 D\delta = 1 Dδ = 1 E α = 6 E\alpha = 6 E α = 6 A γ = 1 A\gamma = 1 A γ = 1 B ε = − 2 B\varepsilon = -2 Bε = − 2 C β = − 3 C\beta = -3 Cβ = − 3 3 5 E ε = − 3 E\varepsilon = -3 Eε = − 3 A β = 5 A\beta = 5 A β = 5 B δ = − 5 B\delta = -5 B δ = − 5 C α = 4 C\alpha = 4 C α = 4 D γ = 6 D\gamma = 6 D γ = 6 7 y . . . l y_{...l} y ... l -18 18 -4 5 9 10 = y . . . . 10 = y_{....} 10 = y ....
表 4.21 火箭推进剂问题的方差分析
变异来源 平方和 自由度 均方 F 0 F_{0} F 0 P P P 值配方 330.00 4 82.50 10.00 0.0033 原材料批次 68.00 4 17.00 操作人员 150.00 4 37.50 测试组件 62.00 4 15.50 误差 66.00 8 8.25 总计 676.00 24
由正交拉丁方对构成格拉科-拉丁方这一概念可以稍作推广。p × p p \times p p × p 超方 (hypersquare)是叠加三个或更多正交 p × p p \times p p × p 拉丁方所得的设计。一般地,如果可以获得一整套 p − 1 p-1 p − 1 个正交拉丁方,则最多可以研究 p + 1 p+1 p + 1 个因子。这样的设计会用尽全部 ( p + 1 ) ( p − 1 ) = p 2 − 1 (p+1)(p-1) = p^{2}-1 ( p + 1 ) ( p − 1 ) = p 2 − 1 个自由度,因此必须独立估计误差方差。当然,使用超方时因子之间必须没有交互作用。