第 5 章补充材料
S5.1 两因子因子设计的期望均方 ¶ 考虑两因子固定效应模型
y i j = μ + τ i + β j + ( τ β ) i j + ε i j k { i = 1 , 2 , … , a j = 1 , 2 , … , b k = 1 , 2 , … , n y_{ij} = \mu + \tau_{i} + \beta_{j} + (\tau\beta)_{ij} + \varepsilon_{ijk} \left\{ \begin{array}{l} i = 1, 2, \dots, a \\ j = 1, 2, \dots, b \\ k = 1, 2, \dots, n \end{array} \right. y ij = μ + τ i + β j + ( τ β ) ij + ε ijk ⎩ ⎨ ⎧ i = 1 , 2 , … , a j = 1 , 2 , … , b k = 1 , 2 , … , n 即教材中的式 (5.1)。我们列出了该模型的期望均方 (expected mean square),但未给出推导。直接应用期望算子来推导期望均方是相对容易的。
考虑求
E ( M S A ) = E ( S S A a − 1 ) = 1 a − 1 E ( S S A ) E\left(MS_{A}\right) = E\left(\frac{SS_{A}}{a-1}\right) = \frac{1}{a-1}E\left(SS_{A}\right) E ( M S A ) = E ( a − 1 S S A ) = a − 1 1 E ( S S A ) 其中 S S A SS_{A} S S A 是行因子的平方和。由于
S S A = 1 b n ∑ i = 1 a y i . . 2 − y . . . 2 a b n SS_{A} = \frac{1}{bn}\sum_{i=1}^{a} y_{i..}^{2} - \frac{y_{...}^{2}}{abn} S S A = bn 1 i = 1 ∑ a y i .. 2 − abn y ... 2 E ( S S A ) = 1 b n E ∑ i = 1 a y i . . 2 − E ( y … 2 a b n ) E\left(SS_{A}\right) = \frac{1}{bn}E\sum_{i=1}^{a} y_{i..}^{2} - E\left(\frac{y_{\dots}^{2}}{abn}\right) E ( S S A ) = bn 1 E i = 1 ∑ a y i .. 2 − E ( abn y … 2 ) 回忆一下,τ . = 0 \tau_{.} = 0 τ . = 0 、β . = 0 \beta_{.} = 0 β . = 0 、( τ β ) . j = 0 (\tau\beta)_{.j} = 0 ( τ β ) . j = 0 、( τ β ) i . = 0 (\tau\beta)_{i.} = 0 ( τ β ) i . = 0 且 ( τ β ) . . = 0 (\tau\beta)_{..} = 0 ( τ β ) .. = 0 ,其中“点”下标表示对该下标求和。现在
y i . . = ∑ j = 1 b ∑ k = 1 n y i j k = b n μ + b n τ i + n β . + n ( τ β ) i . + ε i . . = b n μ + b n τ i + ε i . . \begin{array}{r l} y_{i..} & = \sum_{j=1}^{b}\sum_{k=1}^{n} y_{ijk} = bn\mu + bn\tau_{i} + n\beta_{.} + n(\tau\beta)_{i.} + \varepsilon_{i..} \\ & = bn\mu + bn\tau_{i} + \varepsilon_{i..} \end{array} y i .. = ∑ j = 1 b ∑ k = 1 n y ijk = bn μ + bn τ i + n β . + n ( τ β ) i . + ε i .. = bn μ + bn τ i + ε i .. 以及
1 b n E ∑ i = 1 a y i . . 2 = 1 b n E ∑ i = 1 a [ ( b n μ ) 2 + ( b n ) 2 τ i 2 + ε i . . 2 + 2 ( b n ) 2 μ τ i + 2 b n μ ε i . . + 2 b n τ i ε i . . ] = 1 b n [ a ( b n μ ) 2 + ( b n ) 2 ∑ i = 1 a τ i 2 + a b n σ 2 ] = a b n μ 2 + b n ∑ i = 1 a τ i 2 + a σ 2 \begin{array}{r l} \frac{1}{bn}E\sum_{i=1}^{a} y_{i..}^{2} & = \frac{1}{bn}E\sum_{i=1}^{a}\left[(bn\mu)^{2} + (bn)^{2}\tau_{i}^{2} + \varepsilon_{i..}^{2} + 2(bn)^{2}\mu\tau_{i} + 2bn\mu\varepsilon_{i..} + 2bn\tau_{i}\varepsilon_{i..}\right] \\ & = \frac{1}{bn}\left[a(bn\mu)^{2} + (bn)^{2}\sum_{i=1}^{a}\tau_{i}^{2} + abn\sigma^{2}\right] \\ & = abn\mu^{2} + bn\sum_{i=1}^{a}\tau_{i}^{2} + a\sigma^{2} \end{array} bn 1 E ∑ i = 1 a y i .. 2 = bn 1 E ∑ i = 1 a [ ( bn μ ) 2 + ( bn ) 2 τ i 2 + ε i .. 2 + 2 ( bn ) 2 μ τ i + 2 bn μ ε i .. + 2 bn τ i ε i .. ] = bn 1 [ a ( bn μ ) 2 + ( bn ) 2 ∑ i = 1 a τ i 2 + abn σ 2 ] = abn μ 2 + bn ∑ i = 1 a τ i 2 + a σ 2 此外,我们很容易证明
y … = a b n μ + ε … y_{\dots} = abn\mu + \varepsilon_{\dots} y … = abn μ + ε … 1 a b n E ( y … 2 ) = 1 a b n E ( a b n μ + ε … ) 2 = 1 a b n E [ ( a b n μ ) 2 + ε … 2 + 2 a b n μ ε … ) ] = 1 a b n [ ( a b n μ ) 2 + a b n σ 2 ] = a b n μ 2 + σ 2 \begin{array}{r l} \frac{1}{abn}E(y_{\dots}^{2}) & = \frac{1}{abn}E(abn\mu + \varepsilon_{\dots})^{2} \\ & = \frac{1}{abn}E\left[(abn\mu)^{2} + \varepsilon_{\dots}^{2} + 2abn\mu\varepsilon_{\dots})\right] \\ & = \frac{1}{abn}\left[(abn\mu)^{2} + abn\sigma^{2}\right] \\ & = abn\mu^{2} + \sigma^{2} \end{array} abn 1 E ( y … 2 ) = abn 1 E ( abn μ + ε … ) 2 = abn 1 E [ ( abn μ ) 2 + ε … 2 + 2 abn μ ε … ) ] = abn 1 [ ( abn μ ) 2 + abn σ 2 ] = abn μ 2 + σ 2 因此
E ( M S A ) = E ( S S A a − 1 ) = 1 a − 1 E ( S S A ) 1 a − 1 [ a b n μ 2 + b n ∑ i = 1 a τ i 2 + a σ 2 − ( a b n μ 2 + σ 2 ) ] = 1 a − 1 [ σ 2 ( a − 1 ) + b n ∑ i = 1 a τ i 2 ] = σ 2 + b n ∑ i = 1 a τ i 2 a − 1 \begin{array}{r l} E(MS_{A}) & = E\left(\frac{SS_{A}}{a-1}\right) \\ & = \frac{1}{a-1}E(SS_{A}) \\ & \frac{1}{a-1}\left[abn\mu^{2} + bn\sum_{i=1}^{a}\tau_{i}^{2} + a\sigma^{2} - (abn\mu^{2} + \sigma^{2})\right] \\ & = \frac{1}{a-1}\left[\sigma^{2}(a-1) + bn\sum_{i=1}^{a}\tau_{i}^{2}\right] \\ & = \sigma^{2} + \frac{bn\sum_{i=1}^{a}\tau_{i}^{2}}{a-1} \end{array} E ( M S A ) = E ( a − 1 S S A ) = a − 1 1 E ( S S A ) a − 1 1 [ abn μ 2 + bn ∑ i = 1 a τ i 2 + a σ 2 − ( abn μ 2 + σ 2 ) ] = a − 1 1 [ σ 2 ( a − 1 ) + bn ∑ i = 1 a τ i 2 ] = σ 2 + a − 1 bn ∑ i = 1 a τ i 2 这正是教材中给出的结果。其他期望均方可用类似方法推导。
S5.2 交互作用的定义 ¶ 在 5.1 节中,我们为两因子因子试验同时介绍了效应模型 (effects model)和均值模型 (means model)。若两因子模型中不存在交互作用 (interaction),则
μ i j = μ + τ i + β j \mu_{ij} = \mu + \tau_{i} + \beta_{j} μ ij = μ + τ i + β j 定义行均值和列均值为
μ i . = ∑ j = 1 b μ i j b \mu_{i.} = \frac{\sum_{j=1}^{b}\mu_{ij}}{b} μ i . = b ∑ j = 1 b μ ij μ . j = ∑ i = 1 a μ i j a \mu_{.j} = \frac{\sum_{i=1}^{a}\mu_{ij}}{a} μ . j = a ∑ i = 1 a μ ij 那么,若不存在交互作用,则
μ i j = μ i . + μ . j − μ \mu_{ij} = \mu_{i.} + \mu_{.j} - \mu μ ij = μ i . + μ . j − μ 其中 μ = ∑ i μ i . / a = ∑ j μ . j / b \mu = \sum_{i}\mu_{i.}/a = \sum_{j}\mu_{.j}/b μ = ∑ i μ i . / a = ∑ j μ . j / b 。还可以证明,若不存在交互作用,则每个单元格 (cell)均值都可用另外三个单元格均值表示:
μ i j = μ i j ′ + μ i ′ j − μ i ′ j ′ \mu_{ij} = \mu_{ij^{\prime}} + \mu_{i^{\prime}j} - \mu_{i^{\prime}j^{\prime}} μ ij = μ i j ′ + μ i ′ j − μ i ′ j ′ 这说明了为什么无交互作用的模型有时称为可加模型 (additive model),也说明了为什么我们说处理效应是可加的。
当存在交互作用时,上述关系不再成立。于是交互作用项 ( τ β ) i j (\tau\beta)_{ij} ( τ β ) ij 可以定义为
( τ β ) i j = μ i j − ( μ + τ i + β j ) (\tau\beta)_{ij} = \mu_{ij} - (\mu + \tau_{i} + \beta_{j}) ( τ β ) ij = μ ij − ( μ + τ i + β j ) 或等价地,
( τ β ) i j = μ i j − ( μ i j ′ + μ i ′ j − μ i ′ j ′ ) = μ i j − μ i j ′ − μ i ′ j + μ i ′ j ′ \begin{array}{r} (\tau\beta)_{ij} = \mu_{ij} - (\mu_{ij^{\prime}} + \mu_{i^{\prime}j} - \mu_{i^{\prime}j^{\prime}}) \\ = \mu_{ij} - \mu_{ij^{\prime}} - \mu_{i^{\prime}j} + \mu_{i^{\prime}j^{\prime}} \end{array} ( τ β ) ij = μ ij − ( μ i j ′ + μ i ′ j − μ i ′ j ′ ) = μ ij − μ i j ′ − μ i ′ j + μ i ′ j ′ 因此,只要判断是否所有单元格均值都能表示为 μ i j = μ + τ i + β j \mu_{ij} = \mu + \tau_{i} + \beta_{j} μ ij = μ + τ i + β j ,就可以判断是否存在交互作用。
有时交互作用是响应变量的测量尺度造成的。例如,假定因子效应是乘积形式的,
μ i j = μ τ i β j \mu_{ij} = \mu\tau_{i}\beta_{j} μ ij = μ τ i β j 如果我们假定各因子以相加的方式起作用,那么很快就会发现存在交互作用。这种交互作用可以通过取对数变换来消除,因为
log μ i j = log μ + log τ i + log β j \log\mu_{ij} = \log\mu + \log\tau_{i} + \log\beta_{j} log μ ij = log μ + log τ i + log β j 这说明,如果我们希望得到易于解释的结果(即不存在交互作用),那么响应变量原来的测量尺度并不是最合适的选择。响应变量取对数尺度会更合适。
最后我们注意到,两个因子完全可能发生交互作用,但其中一个(甚至两个)因子的主效应却很小、接近于零。为说明这一点,考虑教材图 5.1 中含交互作用的两因子因子试验。我们已经注意到交互作用很大,AB = -29。然而因子 A 的主效应为 A = 1。因此 A 的主效应小到可以忽略不计。这种情况并不经常出现,通常我们发现交互作用效应不会大于主效应。然而,大的两因子交互作用可能掩盖一个或两个主效应。审慎的试验者需要对这种可能性保持警觉。
S5.3 两因子因子模型中的可估函数 ¶ 两因子因子模型的最小二乘正规方程在教材中由式 (5.14) 给出:
a b n μ ^ + b n ∑ i = 1 a τ ^ i + a n ∑ j = 1 b β j + ∑ i = 1 a ∑ j = 1 b ( τ β ^ ) i j = y … abn\hat{\mu} + bn\sum_{i=1}^{a}\hat{\tau}_{i} + an\sum_{j=1}^{b}\beta_{j} + \sum_{i=1}^{a}\sum_{j=1}^{b}(\tau\hat{\beta})_{ij} = y_{\dots} abn μ ^ + bn i = 1 ∑ a τ ^ i + an j = 1 ∑ b β j + i = 1 ∑ a j = 1 ∑ b ( τ β ^ ) ij = y … b n μ ^ + b n τ ^ i + n ∑ j = 1 b β j + n ∑ j = 1 b ( τ β ^ ) i j = y i . . , i = 1 , 2 , … , a bn\hat{\mu} + bn\hat{\tau}_{i} + n\sum_{j=1}^{b}\beta_{j} + n\sum_{j=1}^{b}(\tau\hat{\beta})_{ij} = y_{i..}, i = 1, 2, \dots, a bn μ ^ + bn τ ^ i + n j = 1 ∑ b β j + n j = 1 ∑ b ( τ β ^ ) ij = y i .. , i = 1 , 2 , … , a a n μ ^ + n ∑ i = 1 a τ ^ i + a n β ^ j + n ∑ i = 1 a ( τ β ^ ) i j = y . j . , j = 1 , 2 , … , b an\hat{\mu} + n\sum_{i=1}^{a}\hat{\tau}_{i} + an\hat{\beta}_{j} + n\sum_{i=1}^{a}(\tau\hat{\beta})_{ij} = y_{.j.}, j = 1, 2, \dots, b an μ ^ + n i = 1 ∑ a τ ^ i + an β ^ j + n i = 1 ∑ a ( τ β ^ ) ij = y . j . , j = 1 , 2 , … , b n μ ^ + n τ ^ i + n β ^ j + n ( τ β ^ ) i j = y i j . , { i = 1 , 2 , … , a j = 1 , 2 , … , b n\hat{\mu} + n\hat{\tau}_{i} + n\hat{\beta}_{j} + n(\tau\hat{\beta})_{ij} = y_{ij.}, \left\{ \begin{array}{l} i = 1, 2, \dots, a \\ j = 1, 2, \dots, b \end{array} \right. n μ ^ + n τ ^ i + n β ^ j + n ( τ β ^ ) ij = y ij . , { i = 1 , 2 , … , a j = 1 , 2 , … , b 回忆一下,一般地,可估函数 (estimable function)必须是正规方程左端的线性组合。考虑比较第 i i i 行与第 i ′ i^{\prime} i ′ 行处理效应的一个对比。该对比为
τ i − τ i ′ + ( τ β ‾ ) i . − ( τ β ‾ ) i ′ . \tau_{i} - \tau_{i^{\prime}} + (\tau\overline{\beta})_{i.} - (\tau\overline{\beta})_{i^{\prime}}. τ i − τ i ′ + ( τ β ) i . − ( τ β ) i ′ . 由于它恰好是两个正规方程之差,所以它是一个可估函数。注意,行因子任意两个水平之差还包含这两行中平均交互作用效应之差。类似地可以证明,任意一对列处理之差也包含这两列中平均交互作用效应之差。一个涉及交互作用的可估函数为
( τ β ) i j − ( τ β ‾ ) i . − ( τ β ‾ ) . j + ( τ β ‾ ) . . (\tau\beta)_{ij} - (\tau\overline{\beta})_{i.} - (\tau\overline{\beta})_{.j} + (\tau\overline{\beta})_{..} ( τ β ) ij − ( τ β ) i . − ( τ β ) . j + ( τ β ) .. 事实上,在效应模型中能够检验的假设必须涉及可估函数。因此,当我们检验无交互作用的假设时,实际上是在检验原假设
H 0 : ( τ β ) i j − ( τ β ‾ ) i . − ( τ β ‾ ) . j + ( τ β ‾ ) . . = 0 对一切 i , j H_{0}: (\tau\beta)_{ij} - (\tau\overline{\beta})_{i.} - (\tau\overline{\beta})_{.j} + (\tau\overline{\beta})_{..} = 0 \text{ 对一切 } i, j H 0 : ( τ β ) ij − ( τ β ) i . − ( τ β ) . j + ( τ β ) .. = 0 对一切 i , j 当我们检验主效应 A 和 B 的假设时,实际上是在检验原假设
H 0 : τ 1 + ( τ β ‾ ) 1. = τ 2 + ( τ β ‾ ) 2. = ⋯ = τ a + ( τ β ‾ ) a . H_{0}: \tau_{1} + (\tau\overline{\beta})_{1.} = \tau_{2} + (\tau\overline{\beta})_{2.} = \dots = \tau_{a} + (\tau\overline{\beta})_{a.} H 0 : τ 1 + ( τ β ) 1. = τ 2 + ( τ β ) 2. = ⋯ = τ a + ( τ β ) a . 以及
H 0 : β 1 + ( τ β ‾ ) . 1 = β 2 + ( τ β ‾ ) . 2 = ⋯ = β b + ( τ β ‾ ) . b H_{0}\colon \beta_{1} + (\tau\overline{\beta})_{.1} = \beta_{2} + (\tau\overline{\beta})_{.2} = \dots = \beta_{b} + (\tau\overline{\beta})_{.b} H 0 : β 1 + ( τ β ) .1 = β 2 + ( τ β ) .2 = ⋯ = β b + ( τ β ) . b 也就是说,我们实际检验的并不是只涉及处理效应相等的假设,而是把处理效应与相应行或列中平均交互作用效应加在一起进行比较的假设。显然,当交互作用很大时,这些假设可能没有多大意义,也没有多大实用价值。这就是教材(5.1 节)指出当交互作用很大时主效应可能没有多大实用价值的原因。此外,当交互作用很大时,关于主效应的统计检验可能并不能真正告诉我们多少关于单个处理效应的信息。有些统计学家在无交互作用的原假设被拒绝时,甚至不再进行主效应检验。
可以证明 [参见 Myers and Milton (1991)],原来的效应模型
y i j k = μ + τ i + β j + ( τ β ) i j + ε i j k y_{ijk} = \mu + \tau_{i} + \beta_{j} + (\tau\beta)_{ij} + \varepsilon_{ijk} y ijk = μ + τ i + β j + ( τ β ) ij + ε ijk 可以重新表示为
y i j k = [ μ + τ ‾ + β ‾ + ( τ β ‾ ) ] + [ τ i − τ ‾ + ( τ β ‾ ) i . − ( τ β ‾ ) ] + [ β j − β ‾ + ( τ β ‾ ) . j − ( τ β ‾ ) ] + [ ( τ β ) i j − ( τ β ‾ ) i . − ( τ β ‾ ) . j + ( τ β ‾ ) ] + ε i j k \begin{array}{c} y_{ijk} = [\mu + \overline{\tau} + \overline{\beta} + (\tau\overline{\beta})] + [\tau_{i} - \overline{\tau} + (\tau\overline{\beta})_{i.} - (\tau\overline{\beta})] + \\ [\beta_{j} - \overline{\beta} + (\tau\overline{\beta})_{.j} - (\tau\overline{\beta})] + [(\tau\beta)_{ij} - (\tau\overline{\beta})_{i.} - (\tau\overline{\beta})_{.j} + (\tau\overline{\beta})] + \varepsilon_{ijk} \end{array} y ijk = [ μ + τ + β + ( τ β )] + [ τ i − τ + ( τ β ) i . − ( τ β )] + [ β j − β + ( τ β ) . j − ( τ β )] + [( τ β ) ij − ( τ β ) i . − ( τ β ) . j + ( τ β )] + ε ijk 或
y i j k = μ ∗ + τ i ∗ + β j ∗ + ( τ β ) i j ∗ + ε i j k y_{ijk} = \mu^{*} + \tau_{i}^{*} + \beta_{j}^{*} + (\tau\beta)_{ij}^{*} + \varepsilon_{ijk} y ijk = μ ∗ + τ i ∗ + β j ∗ + ( τ β ) ij ∗ + ε ijk 可以证明,新参数 μ ∗ \mu^{*} μ ∗ 、τ i ∗ \tau_{i}^{*} τ i ∗ 、β j ∗ \beta_{j}^{*} β j ∗ 和 ( τ β ) i j ∗ (\tau\beta)_{ij}^{*} ( τ β ) ij ∗ 都是可估的。
因此,有理由期望感兴趣的假设可以简单地用这些重新定义的参数来表达。特别地,可以证明不存在交互作用的充要条件是 ( τ β ) i j ∗ = 0 (\tau\beta)_{ij}^{*} = 0 ( τ β ) ij ∗ = 0 。在教材中,我们把无交互作用的原假设表述为:对一切 i i i 和 j j j 有 H 0 : ( τ β ) i j = 0 H_{0}\colon(\tau\beta)_{ij} = 0 H 0 : ( τ β ) ij = 0 。只要理解为我们使用的是用重新定义的(即“带星号”的)参数表示的模型,这种表述并没有错。然而,重要的是要理解,一般地,交互作用并不是一个只与第 ( i j ) (ij) ( ij ) 个单元格有关的参数,它含有来自该单元格、第 i i i 行、第 j j j 列以及总体平均响应的信息。
最后一点是,由于定义了新的“带星号”参数,我们对它们施加了某些约束。特别地,我们有
τ . ∗ = 0 , β . ∗ = 0 , ( τ β ) i . ∗ = 0 , ( τ β ) . j ∗ = 0 且 ( τ β ) . . ∗ = 0 \tau_{.}^{*} = 0, \beta_{.}^{*} = 0, (\tau\beta)_{i.}^{*} = 0, (\tau\beta)_{.j}^{*} = 0 \text{ 且 } (\tau\beta)_{..}^{*} = 0 τ . ∗ = 0 , β . ∗ = 0 , ( τ β ) i . ∗ = 0 , ( τ β ) . j ∗ = 0 且 ( τ β ) .. ∗ = 0 这些就是施加在正规方程上的“通常约束”。此外,主效应的检验变为
H 0 : τ 1 ∗ = τ 2 ∗ = ⋯ = τ a ∗ = 0 H_{0}\colon \tau_{1}^{*} = \tau_{2}^{*} = \dots = \tau_{a}^{*} = 0 H 0 : τ 1 ∗ = τ 2 ∗ = ⋯ = τ a ∗ = 0 以及
H 0 : β 1 ∗ = β 2 ∗ = ⋯ = β b ∗ = 0 H_{0}\colon \beta_{1}^{*} = \beta_{2}^{*} = \dots = \beta_{b}^{*} = 0 H 0 : β 1 ∗ = β 2 ∗ = ⋯ = β b ∗ = 0 教材中就是这样陈述这些假设的,当然,没有带“星号”。
S5.4 两因子因子设计的回归模型形式 ¶ 我们在第 3 章注意到方差分析与回归之间有密切关系,并在第 3 章补充材料中说明了如何把单因子方差分析模型表述为回归模型。现在我们来说明如何把两因子模型表述为回归模型,并使用标准多元回归计算机程序来进行通常的方差分析。
我们用例 5.1 的电池寿命试验来说明这一做法。回忆一下,有三个感兴趣的材料类型(因子 A)和三个温度(因子 B),感兴趣的响应变量是电池寿命。方差分析模型的回归模型表述使用指示变量 (indicator variable)。我们对设计因子材料类型和温度定义指示变量如下:
材料类型 X 1 X_{1} X 1 X 2 X_{2} X 2 1 0 0 2 1 0 3 0 1
温度 X 3 X_{3} X 3 X 4 X_{4} X 4 15 0 0 70 1 0 125 0 1
回归模型为
y i j k = β 0 + β 1 x i j k 1 + β 2 x i j k 2 + β 3 x i j k 3 + β 4 x i j k 4 + β 5 x i j k 1 x i j k 3 + β 6 x i j k 1 x i j k 4 + β 7 x i j k 2 x i j k 3 + β 8 x i j k 2 x i j k 4 + ε i j k (1) \begin{array}{c} y_{ijk} = \beta_{0} + \beta_{1}x_{ijk1} + \beta_{2}x_{ijk2} + \beta_{3}x_{ijk3} + \beta_{4}x_{ijk4} \\ + \beta_{5}x_{ijk1}x_{ijk3} + \beta_{6}x_{ijk1}x_{ijk4} + \beta_{7}x_{ijk2}x_{ijk3} + \beta_{8}x_{ijk2}x_{ijk4} + \varepsilon_{ijk} \end{array}\tag{1} y ijk = β 0 + β 1 x ijk 1 + β 2 x ijk 2 + β 3 x ijk 3 + β 4 x ijk 4 + β 5 x ijk 1 x ijk 3 + β 6 x ijk 1 x ijk 4 + β 7 x ijk 2 x ijk 3 + β 8 x ijk 2 x ijk 4 + ε ijk ( 1 ) 其中 i , j = 1 , 2 , 3 i, j = 1, 2, 3 i , j = 1 , 2 , 3 ,重复次数 k = 1 , 2 , 3 , 4 k = 1, 2, 3, 4 k = 1 , 2 , 3 , 4 。在这个模型中,项 β 1 x i j k 1 + β 2 x i j k 2 \beta_{1}x_{ijk1} + \beta_{2}x_{ijk2} β 1 x ijk 1 + β 2 x ijk 2 表示因子 A(材料类型)的主效应,项 β 3 x i j k 3 + β 4 x i j k 4 \beta_{3}x_{ijk3} + \beta_{4}x_{ijk4} β 3 x ijk 3 + β 4 x ijk 4 表示温度的主效应。这两组项各含两个回归系数,给出两个自由度。式 (1) 中的项 β 5 x i j k 1 x i j k 3 + β 6 x i j k 1 x i j k 4 + β 7 x i j k 2 x i j k 3 + β 8 x i j k 2 x i j k 4 \beta_{5}x_{ijk1}x_{ijk3} + \beta_{6}x_{ijk1}x_{ijk4} + \beta_{7}x_{ijk2}x_{ijk3} + \beta_{8}x_{ijk2}x_{ijk4} β 5 x ijk 1 x ijk 3 + β 6 x ijk 1 x ijk 4 + β 7 x ijk 2 x ijk 3 + β 8 x ijk 2 x ijk 4 表示自由度为四的 AB 交互作用。注意这一项中有四个回归系数。
表 1 给出了该试验的数据,这些数据最初在教材的表 5-1 中给出。在表 1 中,我们给出了该试验 36 次试验中每一次的指示变量。该表中的记号是:对于上述回归模型中的主效应,X i = x i , i = 1 , 2 , 3 , 4 X_{i} = x_{i},\ i = 1, 2, 3, 4 X i = x i , i = 1 , 2 , 3 , 4 ;对于模型中的交互作用项,X 5 = x 1 x 3 X_{5} = x_{1}x_{3} X 5 = x 1 x 3 、X 6 = x 1 x 4 X_{6} = x_{1}x_{4} X 6 = x 1 x 4 、X 7 = x 2 x 3 X_{7} = x_{2}x_{3} X 7 = x 2 x 3 和 X 8 = x 2 x 4 X_{8} = x_{2}x_{4} X 8 = x 2 x 4 。
表 1 例 5.1 数据的回归模型形式
Y Y Y X 1 X_{1} X 1 X 2 X_{2} X 2 X 3 X_{3} X 3 X 4 X_{4} X 4 X 5 X_{5} X 5 X 6 X_{6} X 6 X 7 X_{7} X 7 X 8 X_{8} X 8 130 0 0 0 0 0 0 0 0 34 0 0 1 0 0 0 0 0 20 0 0 0 1 0 0 0 0 150 1 0 0 0 0 0 0 0 136 1 0 1 0 1 0 0 0 25 1 0 0 1 0 1 0 0 138 0 1 0 0 0 0 0 0 174 0 1 1 0 0 0 1 0 96 0 1 0 1 0 0 0 1 155 0 0 0 0 0 0 0 0 40 0 0 1 0 0 0 0 0 70 0 0 0 1 0 0 0 0 188 1 0 0 0 0 0 0 0 122 1 0 1 0 1 0 0 0 70 1 0 0 1 0 1 0 0 110 0 1 0 0 0 0 0 0 120 0 1 1 0 0 0 1 0 104 0 1 0 1 0 0 0 1 74 0 0 0 0 0 0 0 0 80 0 0 1 0 0 0 0 0 82 0 0 0 1 0 0 0 0 159 1 0 0 0 0 0 0 0 106 1 0 1 0 1 0 0 0 58 1 0 0 1 0 1 0 0 168 0 1 0 0 0 0 0 0 150 0 1 1 0 0 0 1 0 82 0 1 0 1 0 0 0 1 180 0 0 0 0 0 0 0 0 75 0 0 1 0 0 0 0 0 58 0 0 0 1 0 0 0 0 126 1 0 0 0 0 0 0 0 115 1 0 1 0 1 0 0 0 45 1 0 0 1 0 1 0 0 160 0 1 0 0 0 0 0 0 139 0 1 1 0 0 0 1 0 60 0 1 0 1 0 0 0 1
该表被用作 Minitab 回归程序的输入,拟合式 (1) 得到如下结果:
Regression Analysis
The regression equation is
y = 135 + 21.0 x1 + 9.2 x2 - 77.5 x3 - 77.2 x4 + 41.5 x5 - 29.0 x6 + 79.2 x7 + 18.7 x8minitab Output (Continued)
Predictor Coef StDev T P Constant 134.75 12.99 10.37 0.000 x1 21.00 18.37 1.14 0.263 x2 9.25 18.37 0.50 0.619 x3 -77.50 18.37 -4.22 0.000 x4 -77.25 18.37 -4.20 0.000 x5 41.50 25.98 1.60 0.122 x6 -29.00 25.98 -1.12 0.274 x7 79.25 25.98 3.05 0.005 x8 18.75 25.98 0.72 0.477 S = 25.98 R-Sq = 76.5% R-Sq(adj) = 69.6% Analysis of Variance Source DF SS MS F P Regression 8 59416.2 7427.0 11.00 0.000 Residual Error 27 18230.7 675.2 Total 35 77647.0 Source DF Seq SS x1 1 141.7 x2 1 10542.0 x3 1 76.1 x4 1 39042.7 x5 1 788.7 x6 1 1963.5 x7 1 6510.0 x8 1 351.6
首先查看上述显示中的方差分析信息。注意自由度为 8 的回归平方和等于教材表 5.5 中材料类型和温度主效应的平方和与交互作用平方和之和。此外,回归的自由度 (8) 等于表 5.5 中主效应与交互作用的自由度之和 (2 + 2 + 4)。上述方差分析显示中的 F 检验可以看作检验所有模型系数都为零这一原假设,即不存在显著的主效应或交互作用效应,备择假设为至少有一个模型参数非零。显然这一假设被拒绝。有些处理产生了显著效应。
现在考虑上述显示底部的“序贯平方和”。回忆一下,X 1 X_{1} X 1 与 X 2 X_{2} X 2 表示材料类型的主效应。序贯平方和是按“效应按顺序加入”的方法计算的,其中“按顺序”指变量在模型中的列出顺序。现在
S S 材料类型 = S S ( X 1 ) + S S ( X 2 ∣ X 1 ) = 141.7 + 10542.0 = 10683.7 SS_{\text{材料类型}} = SS(X_{1}) + SS(X_{2} \mid X_{1}) = 141.7 + 10542.0 = 10683.7 S S 材料类型 = SS ( X 1 ) + SS ( X 2 ∣ X 1 ) = 141.7 + 10542.0 = 10683.7 这正是表 5.5 中材料类型的平方和。记号 S S ( X 2 ∣ X 1 ) SS(X_{2} \mid X_{1}) SS ( X 2 ∣ X 1 ) 表示这是一个“序贯”平方和;也就是说,它是在变量 X 1 X_{1} X 1 已经在回归模型中的条件下变量 X 2 X_{2} X 2 的平方和。
类似地,
S S 温度 = S S ( X 3 ∣ X 1 , X 2 ) + S S ( X 4 ∣ X 1 , X 2 , X 3 ) = 76.1 + 39042.7 = 39118.8 SS_{\text{温度}} = SS(X_{3} \mid X_{1}, X_{2}) + SS(X_{4} \mid X_{1}, X_{2}, X_{3}) = 76.1 + 39042.7 = 39118.8 S S 温度 = SS ( X 3 ∣ X 1 , X 2 ) + SS ( X 4 ∣ X 1 , X 2 , X 3 ) = 76.1 + 39042.7 = 39118.8 这与表 5.5 中温度的平方和十分吻合。最后注意,表 5.5 中的交互作用平方和为
S S 交互作用 = S S ( X 5 ∣ X 1 , X 1 , X 3 , X 4 ) + S S ( X 6 ∣ X 1 , X 1 , X 3 , X 4 , X 5 ) + S S ( X 7 ∣ X 1 , X 2 , X 3 , X 4 , X 5 , X 6 ) + S S ( X 8 ∣ X 1 , X 2 , X 3 , X 4 , X 5 , X 6 , X 7 ) = 788.7 + 1963.5 + 6510.0 + 351.6 = 9613.8 \begin{array}{r l} SS_{\text{交互作用}} & = SS(X_{5} | X_{1}, X_{1}, X_{3}, X_{4}) + SS(X_{6} | X_{1}, X_{1}, X_{3}, X_{4}, X_{5}) \\ & \quad + SS(X_{7} | X_{1}, X_{2}, X_{3}, X_{4}, X_{5}, X_{6}) + SS(X_{8} | X_{1}, X_{2}, X_{3}, X_{4}, X_{5}, X_{6}, X_{7}) \\ & = 788.7 + 1963.5 + 6510.0 + 351.6 = 9613.8 \end{array} S S 交互作用 = SS ( X 5 ∣ X 1 , X 1 , X 3 , X 4 ) + SS ( X 6 ∣ X 1 , X 1 , X 3 , X 4 , X 5 ) + SS ( X 7 ∣ X 1 , X 2 , X 3 , X 4 , X 5 , X 6 ) + SS ( X 8 ∣ X 1 , X 2 , X 3 , X 4 , X 5 , X 6 , X 7 ) = 788.7 + 1963.5 + 6510.0 + 351.6 = 9613.8 当设计是平衡的,即每个单元格中的观测个数相等时,可以证明这种使用序贯平方和的模型回归方法给出的结果与“通常的”方差分析完全相同。此外,由于设计的平衡性,变量 A 与 B 的顺序无关紧要。
对总模型平方和按“效应按顺序加入”进行分解有时称为 Type 1 分析。这一术语在 SAS 统计软件包中很常见,但其他作者和软件系统也使用它。另一种分解方式是,把每个效应看作最后加入到一个已包含其他所有效应的模型中。这种“效应最后加入”的方法通常称为 Type 3 分析。
还有一种方法可以利用两因子因子设计的回归模型表述来生成主效应和交互作用的标准 F 检验。考虑拟合式 (1),并令上述 Minitab 输出中该模型的回归平方和为全模型的模型平方和。于是
S S 模型 ( F M ) = 59416.2 \mathrm{SS}_{\text{模型}}(\mathrm{FM}) = 59416.2 SS 模型 ( FM ) = 59416.2 自由度为 8。假设我们要检验不存在交互作用这一假设。用模型 (1) 表示,无交互作用的假设为
H 0 : β 5 = β 6 = β 7 = β 8 = 0 H 0 : 至少有一个 β j ≠ 0 , j = 5 , 6 , 7 , 8 (2) \begin{array}{l} H_{0}: \beta_{5} = \beta_{6} = \beta_{7} = \beta_{8} = 0 \\ H_{0}: \text{ 至少有一个 } \beta_{j} \neq 0, j = 5, 6, 7, 8 \end{array}\tag{2} H 0 : β 5 = β 6 = β 7 = β 8 = 0 H 0 : 至少有一个 β j = 0 , j = 5 , 6 , 7 , 8 ( 2 ) 当原假设为真时,一个约简模型为
y i j k = β 0 + β 1 x i j k 1 + β 2 x i j k 2 + β 3 x i j k 3 + β 4 x i j k 4 + ε i j k (3) y_{ijk} = \beta_{0} + \beta_{1}x_{ijk1} + \beta_{2}x_{ijk2} + \beta_{3}x_{ijk3} + \beta_{4}x_{ijk4} + \varepsilon_{ijk}\tag{3} y ijk = β 0 + β 1 x ijk 1 + β 2 x ijk 2 + β 3 x ijk 3 + β 4 x ijk 4 + ε ijk ( 3 ) 用 Minitab 拟合式 (2) 得到如下结果:
The regression equation is
y = 122 + 25.2 x1 + 41.9 x2 - 37.3 x3 - 80.7 x4Predictor Coef StDev T P
Constant 122.47 11.17 10.97 0.000
x1 25.17 12.24 2.06 0.048
x2 41.92 12.24 3.43 0.002
x3 -37.25 12.24 -3.04 0.005
x4 -80.67 12.24 -6.59 0.000
S = 29.97 R-Sq = 64.1% R-Sq(adj) = 59.5%
Analysis of Variance
Source DF SS MS F P
Regression 4 49802 12451 13.86 0.000
Residual Error 31 27845 898
Total 35 77647该约简模型的模型平方和为
S S 模型 ( R M ) = 49802.0 \mathrm{SS}_{\text{模型}}(\mathrm{RM}) = 49802.0 SS 模型 ( RM ) = 49802.0 自由度为 4。无交互作用假设 (2) 的检验用“额外”平方和来进行
S S 模型 ( 交互作用 ) = S S 模型 ( F M ) − S S 模型 ( R M ) = 59 , 416.2 − 49 , 812.0 = 9 , 604.2 \begin{array}{r l} \mathrm{SS}_{\text{模型}}(\text{交互作用}) & = \mathrm{SS}_{\text{模型}}(\mathrm{FM}) - \mathrm{SS}_{\text{模型}}(\mathrm{RM}) \\ & = 59,416.2 - 49,812.0 \\ & = 9,604.2 \end{array} SS 模型 ( 交互作用 ) = SS 模型 ( FM ) − SS 模型 ( RM ) = 59 , 416.2 − 49 , 812.0 = 9 , 604.2 自由度为 8 − 4 = 4 8 - 4 = 4 8 − 4 = 4 。除了 Minitab 报告结果时的舍入误差之外,这个量就是教材表 5.5 中原方差分析的交互作用平方和。它度量的是在拟合主效应之后再拟合交互作用的贡献。
现在考虑检验材料类型没有主效应。用式 (1) 表示,假设为
H 0 : β 1 = β 2 = 0 H 0 : 至少有一个 β j ≠ 0 , j = 1 , 2 (4) \begin{array}{l} H_{0}: \beta_{1} = \beta_{2} = 0 \\ H_{0}: \text{ 至少有一个 } \beta_{j} \neq 0, j = 1, 2 \end{array}\tag{4} H 0 : β 1 = β 2 = 0 H 0 : 至少有一个 β j = 0 , j = 1 , 2 ( 4 ) 由于我们使用的是平衡设计,检验这一假设只需拟合模型
y i j k = β 0 + β 1 x i j k 1 + β 2 x i j k 2 + ε i j k (5) y_{ijk} = \beta_{0} + \beta_{1}x_{ijk1} + \beta_{2}x_{ijk2} + \varepsilon_{ijk}\tag{5} y ijk = β 0 + β 1 x ijk 1 + β 2 x ijk 2 + ε ijk ( 5 ) 用 Minitab 拟合该模型得到
Regression Analysis
The regression equation is
y = 83.2 + 25.2 x1 + 41.9 x2
Predictor Coef StDev T P
Constant 83.17 13.00 6.40 0.000
x1 25.17 18.39 1.37 0.180
x2 41.92 18.39 2.28 0.029Analysis of Variance
Source DF SS MS F P
Regression 2 10684 5342 2.63 0.087
Residual Error 33 66963 2029
Total 35 77647注意该模型 [式 (5)] 的回归平方和与教材表 5-5 中材料类型的平方和基本相同。类似地,检验不存在温度效应等价于检验
H 0 : β 3 = β 4 = 0 H 0 : 至少有一个 β j ≠ 0 , j = 3 , 4 (6) \begin{array}{l} H_{0}: \beta_{3} = \beta_{4} = 0 \\ H_{0}: \text{ 至少有一个 } \beta_{j} \neq 0, j = 3, 4 \end{array}\tag{6} H 0 : β 3 = β 4 = 0 H 0 : 至少有一个 β j = 0 , j = 3 , 4 ( 6 ) 为检验 (6) 中的假设,我们只需拟合模型
y i j k = β 0 + β 3 x i j k 3 + β 4 x i j k 4 + ε i j k (7) y_{ijk} = \beta_{0} + \beta_{3}x_{ijk3} + \beta_{4}x_{ijk4} + \varepsilon_{ijk}\tag{7} y ijk = β 0 + β 3 x ijk 3 + β 4 x ijk 4 + ε ijk ( 7 ) Minitab 回归输出为
Regression Analysis
The regression equation is
y = 145 - 37.3 x3 - 80.7 x4
Predictor Coef StDev T P
Constant 144.833 9.864 14.68 0.000
x3 -37.25 13.95 -2.67 0.012
x4 -80.67 13.95 -5.78 0.000
S = 34.17 R-Sq = 50.4% R-Sq(adj) = 47.4%
Analysis of Variance
Source DF SS MS F P
Regression 2 39119 19559 16.75 0.000
Residual Error 33 38528 1168
Total 35 77647注意该模型即式 (7) 的回归平方和与表 5.5 中温度主效应的平方和基本相同。
S5.5 模型层次 ¶ 在例 5.4 中,我们用电池寿命试验(例 5.1)的数据演示了当两因子因子试验中的一个因子是定量因子、另一个是定性因子时如何拟合响应曲线。此时因子是温度 (A) 和材料类型 (B)。使用 Design-Expert 软件包,我们拟合了一个包含材料类型主效应、温度的线性与二次效应、材料类型与温度线性效应的交互作用,以及材料类型与温度二次效应的交互作用的模型。参见教材中的表 5.15。通过考察该表,我们观察到温度的二次效应以及材料类型与温度线性效应的交互作用不显著;也就是说,它们的 P 值相当大。我们把这两个不显著的项留在模型中,以保持层次。
层次原则 (hierarchy principle)指出,如果模型中含有较高阶的项,那么它也应含有组成该项的所有较低阶项。因此,如果模型中含有一个二阶项(例如交互作用),那么该交互作用涉及的所有主效应以及涉及这些因子的所有较低阶交互作用也应包含在模型中。
有时层次性是有道理的。一般地,如果模型将用于解释目的,那么层次模型相当合乎逻辑。另一方面,也可能存在非层次模型更有道理的场合。为说明这一点,考虑例 5.4 在表 2 中给出的另一个分析,它由 Design-Expert 得到。我们选用了一个非层次模型,其中不含温度的二次效应(它很可能是最弱的效应),但温度与材料类型交互作用的两个自由度为二的成分都在模型中。从表 2 可以看出,非层次模型的残差均方更小(653.81,而表 5.15 为 675.21)。这一点很重要,因为残差均方可以看作模型未能解释的残余变异的方差。也就是说,非层次模型实际上对试验数据拟合得更好。
还要注意,非层次模型中模型参数的标准误更小。这表明,去掉不显著的项可以使参数估计的精度更高,尽管这样得到的模型不满足层次原则。此外,注意层次模型中模型参数的 95% 置信区间总是比非层次模型中相应的置信区间更长。在这个例子中,非层次模型确实比层次模型给出了更好的因子效应估计。
表 2 非层次模型的 Design-Expert 输出,例 5.4
ANOVA for Response Surface Reduced Cubic Model
Analysis of variance table [Partial sum of squares]
Source Sum of Squares DF Mean Square F Value Prob > F Model 59340.17 7 8477.17 12.97 < 0.0001 A 10683.72 2 5341.86 8.17 0.0016 B 39042.67 1 39042.67 59.72 < 0.0001 AB 2315.08 2 1157.54 1.77 0.1888 A B 2 AB^{2} A B 2 7298.69 2 3649.35 5.58 0.0091 Residual 18306.81 28 653.81 Lack of Fit 76.06 1 76.06 0.11 0.7398 Pure Error 18230.75 27 675.21 Cor Total 77646.97 35
Std. Dev. 25.57 R-Squared 0.7642 Mean 105.53 Adj R-Squared 0.7053 C.V. 24.23 Pred R-Squared 0.6042 PRESS 30729.09 Adeq Precision 8.815
Term Estimate DF Standard Error 95% CI Low 95% CI High Intercept 105.53 1 4.26 96.80 114.26 A[1] -50.33 1 10.44 -71.72 -28.95 A[2] 12.17 1 10.44 -9.22 33.55 B-Temp -40.33 1 5.22 -51.02 -29.64 A[1]B 1.71 1 7.38 -13.41 16.83 A[2]B -12.79 1 7.38 -27.91 2.33 A [ 1 ] B 2 A[1]B^{2} A [ 1 ] B 2 41.96 1 12.78 15.77 68.15 A [ 2 ] B 2 A[2]B^{2} A [ 2 ] B 2 -14.04 1 12.78 -40.23 12.15
补充参考文献 ¶ Myers, R. H. and Milton, J. S. (1991), A First Course in the Theory of the Linear Model, PWS-Kent, Boston, MA.