【數(shù)學(xué)建模算法】(32)方差分析（下）

2.雙因素方差分析

如果要考慮兩個(gè)因素A,B對(duì)指標(biāo)的影響，A,B各劃分幾個(gè)水平，對(duì)每一個(gè)水平組合作若干次試驗(yàn)，對(duì)所得數(shù)據(jù)進(jìn)行方差分析，檢驗(yàn)兩因素是否分別對(duì)指標(biāo)有顯著影響，或者還要進(jìn)一步檢驗(yàn)兩因素是否對(duì)指標(biāo)有顯著的交互影響。

2.1.數(shù)學(xué)模型

設(shè)A取r個(gè)水平A_{1}, A_{2}, \cdots, A_{r}，B取s個(gè)水平B_{1}, B_{2}, \cdots, B_{s}，在水平組合\left(A_{i}, B_{j}\right)下總體\mathcal{X}_{i j}服從正態(tài)分布N\left(\mu_{i j}, \sigma^{2}\right), i=1, \cdots, r, j=1, \cdots, s。又設(shè)在水平組合\left(A_{i}, B_{j}\right)下作了t個(gè)試驗(yàn)，所得結(jié)果記作x_{i j k}, \quad x_{i j k}服從N\left(\mu_{i j}, \sigma^{2}\right)，i=1, \cdots, r, \quad j=1, \cdots, s，k=1, \cdots, t，且相互獨(dú)立。將這些數(shù)據(jù)列成下表的形式。

雙因素試驗(yàn)數(shù)據(jù)表

B_{1} B_{2} ... B_{s}
A_{1} x_{111} \cdots x_{11 t} x_{121} \cdots x_{12 t} ... x_{1s1} \cdots x_{1s t}
A_{2} x_{211} \cdots x_{21 t} x_{221} \cdots x_{22 t} ... x_{2s1} \cdots x_{2s t}
... ... ... ... ...
A_{r} x_{r11} \cdots x_{r1 t} x_{r21} \cdots x_{r2 t} ... x_{rs1} \cdots x_{rs t}

將x_{i j k}分解為：
x_{i j k}=\mu_{i j}+\varepsilon_{i j k}, \quad i=1, \cdots, r, \quad j=1, \cdots, s, \quad k=1, \cdots, t（14）

其中\varepsilon_{i j k} \sim N\left(0, \sigma^{2}\right)，且相互獨(dú)立。記
\mu=\frac{1}{r S} \sum_{i=1}^{r} \sum_{j=1}^{s} \mu_{i j}, \quad \mu_{i .}=\frac{1}{s} \sum_{j=1}^{s} \mu_{i j}, \quad \alpha_{i}=\mu_{i \bullet}-\mu
\mu_{ . j}=\frac{1}{r} \sum_{i=1}^{r} \mu_{i j}, \quad \beta_{j}=\mu_{ \bullet j}-\mu, \quad \gamma_{i j}=\mu_{i j}-\mu-\alpha_{i}-\beta_{j}（15）

\mu是總均值，\alpha_{i}是水平A_{i}對(duì)指標(biāo)的效應(yīng)，\beta_{j}是水平B_{j}對(duì)指標(biāo)的效應(yīng)，\gamma_{i j}是水平A_{i}與B_{j}對(duì)指標(biāo)的交互效應(yīng)。模型表為：
\left\{\begin{array}{l}{x_{i j k}=\mu+\alpha_{i}+\beta_{j}+\gamma_{i j}+\varepsilon_{i j k}} \\ {\sum_{i=1}^{r} \alpha_{i}=0, \sum_{j=1}^{s} \beta_{j}=0, \sum_{i=1}^{r} \gamma_{i j}=\sum_{j=1}^{s} \gamma_{i j}=0} \\ {\varepsilon_{i j k} \sim N\left(0, \sigma^{2}\right), i=1, \cdots, r, j=1, \cdots, s, k=1, \cdots, t}\end{array}\right.（16）

原假設(shè)為：
H_{01} : \alpha_{i}=0(i=1, \cdots, r)（17）
H_{02} : \beta_{j}=0(j=1, \cdots, s)（18）
H_{03} : \gamma_{i i}=0(i=1, \cdots, r ; j=1, \cdots, s)（19）

2.2.無(wú)交互影響的雙因素方差分析

如果根據(jù)經(jīng)驗(yàn)或某種分析能夠事先判定兩因素之間沒有交互影響，每組試驗(yàn)就不必重復(fù)，即可令t=1，過(guò)程大為簡(jiǎn)化。

假設(shè)\gamma_{i j}=0，于是：
\mu_{i j}=\mu+\alpha_{i}+\beta_{j}, \quad i=1, \cdots, r, \quad j=1, \cdots, s
此時(shí)，模型（16）可寫成：
\left\{\begin{array}{l}{x_{i j}=\mu+\alpha_{i}+\beta_{j}+\varepsilon_{i j}} \\ {\sum_{i=1}^{r} \alpha_{i}=0, \sum_{j=1}^{s} \beta_{j}=0} \\ {\varepsilon_{i j} \sim N\left(0, \sigma^{2}\right), i=1, \cdots, r, j=1, \cdots, s}\end{array}\right.（20）

對(duì)這個(gè)模型我們所要檢驗(yàn)的假設(shè)為式（17）和式（18）。下面采用與單因素方差分析模型類似的方法導(dǎo)出檢驗(yàn)統(tǒng)計(jì)量。
記：
\overline{x}=\frac{1}{r S} \sum_{i=1}^{r} \sum_{j=1}^{s} x_{i j}, \quad \overline{x}_{i .}=\frac{1}{S} \sum_{j=1}^{s} x_{i j}, \overline{x}_{\bullet j}=\frac{1}{r} \sum_{i=1}^{r} x_{i j}
S_{T}=\sum_{i=1}^{r} \sum_{j=1}^{s}\left(x_{i j}-\overline{x}\right)^{2}

其中S_{T}為全部試驗(yàn)數(shù)據(jù)的總變差，稱為總平方和，對(duì)其進(jìn)行分解：
\begin{aligned} S_{T} &=\sum_{i=1}^{r} \sum_{j=1}^{s}\left(x_{i j}-\overline{x}\right)^{2} \\ &=\sum_{i=1}^{r} \sum_{s=1}^{s}\left(x_{i j}-\overline{x}_{i \cdot}-\overline{x}_{\cdot j}+\overline{x}\right)^{2}+s \sum_{i=1}^{r}\left(\overline{x}_{i \cdot}-\overline{x}\right)^{2}+r \sum_{j=1}^{s}\left(\overline{x}_{\cdot j}-\overline{x}\right)^{2} \\ &=S_{E}+S_{A}+S_{B} \end{aligned}

可以驗(yàn)證，在上述平方和分解中交叉項(xiàng)均為 0。其中：
S_{E}=\sum_{i=1}^{r} \sum_{s=1}^{s}\left(x_{i j}-\overline{x}_{i .}-\overline{x}_{. j}+\overline{x}\right)^{2}
S_{A}=s \sum_{i=1}^{r}\left(\overline{x}_{i .}-\overline{x}\right)^{2}, \quad S_{B}=r \sum_{j=1}^{s}\left(\overline{x}_{ . j}-\overline{x}\right)^{2}

我們先來(lái)看看S_{A}的統(tǒng)計(jì)意義。因?yàn)?img class="math-inline" src="https://math.jianshu.com/math?formula=%5Coverline%7Bx%7D_%7Bi%20%5Cbullet%7D" alt="\overline{x}_{i \bullet}" mathimg="1">是水平A_{i}下所有觀測(cè)值的平均，所以\sum_{i=1}^{r}\left(\overline{x}_{i .}-\overline{x}\right)^{2}反映了\overline{x}_{1 .}, \overline{x}_{2 .}, \cdots, \overline{x}_{r .}差異的程度。這種差異是由于因素A的不同水平所引起的，因此S_{A}稱為因素A平方和。類似地，S_{B}稱為因素B的平方和。至于S_{E}的意義不甚明顯，我們可以這樣來(lái)理解：因?yàn)?br> S_{E}=S_{T}-S_{A}-S_{B}（21）

在我們所考慮的兩因素問(wèn)題中，除了因素A和B之外，剩余的再?zèng)]有其它系統(tǒng)性因素的影響，因此從總平方和中減去S_{A}和S_{B}之后，剩下的數(shù)據(jù)變差只能歸入隨機(jī)誤差，故S_{E}反映了試驗(yàn)的隨機(jī)誤差。
有了總平方和的分解式：
S_{T}=S_{E}+S_{A}+S_{B}
以及各個(gè)平方和的統(tǒng)計(jì)意義，我們就可以明白，假設(shè)（17）的檢驗(yàn)統(tǒng)計(jì)量應(yīng)取為S_{A}與S_{E}的比。

和一元方差分析相類似，可以證明，當(dāng)H_{01}成立時(shí)，
F_{A}=\frac{\frac{S_{A}}{r-1}}{\frac{S_{E}}{(r-1)(s-1)}} \sim F(r-1,(r-1)(s-1))（22）
當(dāng)H_{02}成立時(shí)：
F_{B}=\frac{\frac{S_{B}}{s-1}}{\frac{S_{E}}{(r-1)(s-1)}} \sim F(s-1,(r-1)(s-1))（23）
檢驗(yàn)規(guī)則為：
F_{A}<F_{1-\alpha}(r-1,(r-1)(s-1))時(shí)接受H_{01}，否則拒絕H_{01}；
F_{B}<F_{1-\alpha}(s-1,(r-1)(s-1))時(shí)接受H_{02}，否則拒絕H_{02}；
可以寫出方差分析表：

無(wú)交互效應(yīng)的兩因素方差分析表

方差來(lái)源 平方和 自由度 均方 F比
因素A S_{A} r-1 \overline{S}_{A}=\frac{S_{A}}{r-1} \frac{\overline{S}_{A}}{\overline{S}_{E}}
因素B S_{B} s-1 \overline{S}_{B}=\frac{S_{B}}{s-1} \frac{\overline{S}_{B}}{\overline{S}_{E}}
誤差 S_{E} (r-1)(s-1) \overline{S}_{E}=\frac{S_{E}}{(r-1)(s-1)}
總和 S_{T} rs-1

2.3.關(guān)于交互效應(yīng)的雙因素方差分析

與前面方法類似，記：
\overline{x}=\frac{1}{r s t} \sum_{i=1}^{r} \sum_{j=1}^{s} \sum_{k=1}^{t} x_{i j k}, \quad \overline{x}_{i j \bullet}=\frac{1}{t} \sum_{k=1}^{t} x_{i j k}
\overline{x}_{i \bullet \bullet}=\frac{1}{s t} \sum_{j=1}^{s} \sum_{k=1}^{t} x_{i j k}, \quad \overline{x}_{\bullet j \bullet}=\frac{1}{r t} \sum_{i=1}^{r} \sum_{k=1}^{t} x_{i j k}

將全體數(shù)據(jù)對(duì)\overline{x}的偏差平方和：
S_{T}=\sum_{i=1}^{r} \sum_{j=1}^{s} \sum_{k=1}^{t}\left(x_{i j k}-\overline{x}\right)^{2}（24）
進(jìn)行分解，可得：
S_{T}=S_{E}+S_{A}+S_{B}+S_{A B}（25）
其中：
S_{E}=\sum_{i=1}^{r} \sum_{j=1}^{s} \sum_{k=1}^{t}\left(x_{i j k}-\overline{x}_{i j .}\right)^{2}, \quad S_{A}=s t \sum_{i=1}^{r}\left(\overline{x}_{i . .}-\overline{x}\right)^{2}
S_{B}=r t \sum_{j=1}^{s}\left(\overline{x}_{. j .}-\overline{x}\right)^{2}, \quad S_{A B}=t \sum_{i=1}^{r} \sum_{j=1}^{s}\left(\overline{x}_{i j .}-\overline{x}_{i . .}-\overline{x}_{c_{j}, c}+\overline{x}\right)^{2}
稱S_{E}為誤差平方和，S_{A}為因素A的平方和（或行間平方和），S_{B}為因素B的平方和（或列間平方和），S_{A B}為交互作用的平方和（或格間平方和）。
可以證明，當(dāng)H_{03}成立時(shí)：
F_{A B}=\frac{\frac{S_{A B}}{(r-1)(s-1)}}{\frac{S_{E}}{r s(t-1)}} \sim F((r-1)(s-1), r s(t-1))（26）
據(jù)此統(tǒng)計(jì)量，可以檢驗(yàn)H_{03}。

檢驗(yàn)因子A和B各個(gè)水平的效應(yīng)是否有差異，與 2.2 中的檢驗(yàn)是一樣的。

雙因素方差分析表

方差來(lái)源 平方和 自由度 均方 F比
因素A S_{A} r-1 \overline{S}_{A}=\frac{S_{A}}{r-1} \frac{\overline{S}_{A}}{\overline{S}_{E}}
因素B S_{B} s-1 \overline{S}_{B}=\frac{S_{B}}{s-1} \frac{\overline{S}_{B}}{\overline{S}_{E}}
交互效應(yīng) S_{A B} (r-1)(s-1) \overline{S}_{A B}=\frac{S_{A B}}{(r-1)(s-1)} \frac{\overline{S}_{A B}}{\overline{S}_{E}}
誤差 S_{E} rs(t-1) \overline{S}=\frac{S_{E}}{r s(t-1)}
總和 S_{T} rst-1

2.4.Matlab實(shí)現(xiàn)

統(tǒng)計(jì)工具箱中用 anova2 作雙因素方差分析。命令為

p=anova2(x,reps)

其中 x 不同列的數(shù)據(jù)表示單一因素的變化情況，不同行中的數(shù)據(jù)表示另一因素的變化情況。如果每種行—列對(duì)（“單元”）有不止一個(gè)的觀測(cè)值，則用參數(shù) reps 來(lái)表明每個(gè)“單元”多個(gè)觀測(cè)值的不同標(biāo)號(hào)，即 reps 給出重復(fù)試驗(yàn)的次數(shù) t 。下面的矩陣中，列因素有 3 種水平，行因素有兩種水平，但每組水平有兩組樣本，相應(yīng)地用下標(biāo)來(lái)標(biāo)識(shí)：
\left[\begin{array}{lll}{x_{111}} & {x_{121}} & {x_{131}} \\ {x_{112}} & {x_{122}} & {x_{132}} \\ {x_{211}} & {x_{221}} & {x_{231}} \\ {x_{212}} & {x_{222}} & {x_{232}}\end{array}\right]

例3 一種火箭使用了四種燃料、三種推進(jìn)器，進(jìn)行射程試驗(yàn)，對(duì)于每種燃料與每種推進(jìn)器的組合作一次試驗(yàn)，得到試驗(yàn)數(shù)據(jù)如下表。問(wèn)各種燃料之間及各種推進(jìn)器之間有無(wú)顯著差異？

火箭試驗(yàn)數(shù)據(jù)

B_{1} B_{2} B_{3}
A_{1} 58.2 56.2 65.3
A_{2} 49.1 54.1 51.6
A_{3} 60.1 70.9 39.2
A_{4} 75.8 58.2 48.7

解：記燃料為因素A，它有4個(gè)水平，水平效應(yīng)為\alpha_{i}，i=1,2,3,4。推進(jìn)器為因素B ，它有 3 個(gè)水平，水平效應(yīng)為\beta_{j}, j=1,2,3。我們?cè)陲@著性水平\alpha=0.05下檢驗(yàn)：
H_{1} : \alpha_{1}=\alpha_{2}=\alpha_{3}=\alpha_{4}=0
H_{2} : \beta_{1}=\beta_{2}=\beta_{3}=0
編寫如下的Matlab程序：

x=[58.2 56.2 65.3
49.1 54.1 51.6
60.1 70.9 39.2
75.8 58.2 48.7];
[p,t,st]=anova2(x)

求得p=0.4491 0.7387，表明各種燃料和各種推進(jìn)器之間的差異對(duì)于火箭射程無(wú)顯著影響。

例4 一火箭使用了 4 種燃料，3 種推進(jìn)器作射程試驗(yàn)，每種燃料與每種推進(jìn)器的組合各發(fā)射火箭 2 次，得到如下表結(jié)果。

B_{1} B_{2} B_{3}
A_{1} 58.2,52.6 56.2,41.2 65.3,60.8
A_{2} 49.1,42.8 54.1,50.5 51.6,48.4
A_{3} 60.1,58.3 70.9,73.2 39.2,40.7
A_{4} 75.8,71.5 58.2,51.0 48.7,41.4

試在水平 0.05 下，檢驗(yàn)不同燃料（因素 A ）、不同推進(jìn)器（因素 B ）下的射程是否有顯著差異？交互作用是否顯著？

解 編寫程序如下：

clc,clear
x0=[58.2,52.6 56.2,41.2 65.3,60.8
49.1,42.8 54.1,50.5 51.6,48.4
60.1,58.3 70.9,73.2 39.2,40.7
75.8,71.5 58.2,51.0 48.7,41.4];
x1=x0(:,1:2:5);x2=x0(:,2:2:6);
for i=1:4
    x(2*i-1,:)=x1(i,:);
    x(2*i,:)=x2(i,:);
end
[p,t,st]=anova2(x,2)

求得\mathrm{p}=0.0035 \quad 0.0260 \quad 0.0001，表明各試驗(yàn)均值相等的概率都為小概率，故可拒絕均值相等假設(shè)。即認(rèn)為不同燃料（因素A）、不同推進(jìn)器（因素B）下的射程有顯著差異，交互作用也是顯著的。

數(shù)據(jù)非均衡的雙因素方差分析的 Matlab 命令要使用多因素方差分析的命令anovan，具體使用方法參見下面的例 5。

3.正交試驗(yàn)設(shè)計(jì)與方差分析

前面介紹了一個(gè)或兩個(gè)因素的試驗(yàn)，由于因素較少，我們可以對(duì)不同因素的所有可能的水平組合做試驗(yàn)，這叫做全面試驗(yàn)。當(dāng)因素較多時(shí)，雖然理論上仍可采用前面的方法進(jìn)行全面試驗(yàn)后再做相應(yīng)的方差分析，但是在實(shí)際中有時(shí)會(huì)遇到試驗(yàn)次數(shù)太多的問(wèn)題。如三因素四水平的問(wèn)題，所有不同水平的組合有4^{3}=64種，在每一種組合下只進(jìn)行一次試驗(yàn)，也需做 64 次。如果考慮更多的因素及水平，則全面試驗(yàn)的次數(shù)可能會(huì)大得驚人。因此在實(shí)際應(yīng)用中，對(duì)于多因素做全面試驗(yàn)是不現(xiàn)實(shí)的。于是我們考慮是否可以選擇其中一部分組合進(jìn)行試驗(yàn)，這就要用到試驗(yàn)設(shè)計(jì)方法選擇合理的試驗(yàn)方案，使得試驗(yàn)次數(shù)不多，但也能得到比較滿意的結(jié)果。

3.1.用正交表安排試驗(yàn)

正交表是一系列規(guī)格化的表格，每個(gè)表都有一個(gè)記號(hào)，如L_{9}\left(3^{4}\right)，見下表：

1 2 3 4
1 1 1 3 2
2 2 1 1 1
3 3 1 2 3
4 1 2 2 1
5 2 2 3 3
6 3 2 1 2
7 1 3 1 3
8 2 3 2 2
9 3 3 3 1

從上表可見，L_{9}\left(3^{4}\right)有9行，4列，表的主體中只有1,2,3三個(gè)數(shù)字組成。
正交表的組成：
（1）每列中數(shù)字出現(xiàn)的次數(shù)相同，如L_{9}\left(3^{4}\right)表每列中數(shù)字 1，2，3 均出現(xiàn)三次。
（2）任取兩列數(shù)字的搭配是均衡的，如L_{9}\left(3^{4}\right)表里每?jī)闪兄?img class="math-inline" src="https://math.jianshu.com/math?formula=(1%2C1)%2C%20%5Cquad(1%2C2)%2C%20%5Cquad%20%5Ccdots%EF%BC%8C(3%2C3)" alt="(1,1), \quad(1,2), \quad \cdots，(3,3)" mathimg="1">，九種組合各出現(xiàn)一次。

這種均衡性是一般正交表構(gòu)造的特點(diǎn)，它使得根據(jù)正交表安排的試驗(yàn)，其試驗(yàn)結(jié)果
具有很好的可比性，易于進(jìn)行統(tǒng)計(jì)分析。
用正交表安排試驗(yàn)時(shí)，根據(jù)因素和水平個(gè)數(shù)的多少以及試驗(yàn)工作量的大小來(lái)考慮選用哪張正交表，下面舉例說(shuō)明。

例5 為提高某種化學(xué)產(chǎn)品的轉(zhuǎn)化率（％），考慮三個(gè)有關(guān)因素：反應(yīng)溫度A（℃），反應(yīng)時(shí)間B（min）和使用催化劑的含量C（％）。各因素選取三個(gè)水平，如下表所示。

溫度A 時(shí)間B 催化劑含量C
1 80 90 5
2 85 120 6
3 90 150 7

如果做全面試驗(yàn)，則需3^{3}=27次，若用正交表L_{9}\left(3^{4}\right)，僅做 9 次試驗(yàn)。將三個(gè)因素A, B, C分別放在L_{9}\left(3^{4}\right)表的任意三列上，如將A, B分別放在L_{9}\left(3^{4}\right)的第 1，2 列上，C放在L_{9}\left(3^{4}\right)的第 4 列上。將表中A, B, C所在的三列上的數(shù)字 1，2，3 分別用相應(yīng)的因素水平去替代，得 9 次試驗(yàn)方案。以上工作稱為表頭設(shè)計(jì)。再將 9 次試驗(yàn)結(jié)果轉(zhuǎn)化率數(shù)據(jù)列于表上（見下表）。

反應(yīng)溫度A 反應(yīng)時(shí)間B 催化劑含量C 轉(zhuǎn)化率
1 80(1) 90(1) 6(2) 31
2 85(2) 90(1) 5(1) 54
3 90(3) 90(1) 7(3) 38
4 80(1) 120(2) 5(1) 53
5 85(2) 120(2) 7(3) 49
6 90(3) 120(2) 6(2) 42
7 80(1) 150(3) 7(3) 57
8 85(2) 150(3) 6(2) 62
9 90(3) 150(3) 5(1) 64

這里不做統(tǒng)計(jì)分析，直接利用 Matlab 多因素方差分析的函數(shù) anovan 進(jìn)行求解，程序如下：

y=[31 54 38 53 49 42 57 62 64];
g1=[1 2 3 1 2 3 1 2 3];
g2=[1 1 1 2 2 2 3 3 3];
g3=[2 1 3 1 3 2 3 2 1];
[p,t,st]=anovan(y,{g1,g2,g3})

求得概率\mathrm{p}=0.1364 \quad 0.0283 \quad 0.0714，可見因素B,C的各水平對(duì)指標(biāo)值的影響有顯著差異（顯著性水平取 0.1），而因素A 的各水平對(duì)指標(biāo)值的影響無(wú)顯著差異。

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請(qǐng)聯(lián)系作者
【社區(qū)內(nèi)容提示】社區(qū)部分內(nèi)容疑似由AI輔助生成，瀏覽時(shí)請(qǐng)結(jié)合常識(shí)與多方信息審慎甄別。
平臺(tái)聲明：文章內(nèi)容（如有圖片或視頻亦包括在內(nèi)）由作者上傳并發(fā)布，文章內(nèi)容僅代表作者本人觀點(diǎn)，簡(jiǎn)書系信息發(fā)布平臺(tái)，僅提供信息存儲(chǔ)服務(wù)。

相關(guān)閱讀更多精彩內(nèi)容

友情鏈接更多精彩內(nèi)容