R语言基础 | 方差分析（2）：多因素方差分析（上）
專注系列化、高質(zhì)量的R語言教程
推文索引 | 聯(lián)系小編 | 付費合集
當(dāng)存在兩個及以上的分組變量時，可以使用多因素方差分析（N-way ?ANOVA、Multifactor ANOVA）檢驗各組的樣本均值是否存在顯著差異。本篇主要以雙因素方差分析（Two-way ?ANOVA）為例介紹相關(guān)內(nèi)容。
本篇的目錄如下：
2 多因素方差分析
2.1 示例數(shù)據(jù)
2.2 平衡試驗設(shè)計
2.3 I型方差分析
2.4 交互效應(yīng)
2.5 多因素方差分析
未完待續(xù)
2 多因素方差分析
2.1 示例數(shù)據(jù)
本篇使用兩個示例數(shù)據(jù)。
第一個示例數(shù)據(jù)是來自基礎(chǔ)包datasets的npk：
npk ##????block?N?P?K?yield ##?1??????1?0?1?1??49.5 ##?2??????1?1?1?0??62.8 ##?3??????1?0?0?0??46.8 ##?4??????1?1?0?1??57.0 ##?5??????2?1?0?0??59.8 ##?6??????2?1?1?1??58.5 ##?7??????2?0?0?1??55.5 ##?...npk共有5個變量：yield表示產(chǎn)量，為觀測值；block為土地編號，取值為1-6，因子類型；N、P、K分別表示氮、磷、鉀的使用情況，均為二分變量。
第二個示例數(shù)據(jù)是由mtcars數(shù)據(jù)集的四個變量組成，記為mpg：
library(dplyr)? mpg?<-?mtcars?%>%select(mpg,?cyl,?gear,?carb)?%>%mutate(cyl?=?factor(cyl),gear?=?factor(gear),carb?=?factor(carb)) mpg ##??????????????????????mpg?cyl?gear?carb ##?Mazda?RX4???????????21.0???6????4????4 ##?Mazda?RX4?Wag???????21.0???6????4????4 ##?Datsun?710??????????22.8???4????4????1 ##?Hornet?4?Drive??????21.4???6????3????1 ##?Hornet?Sportabout???18.7???8????3????2 ##?Valiant?????????????18.1???6????3????1 ##?Duster?360??????????14.3???8????3????4 ##?...以mpg為觀測變量，cyl、gear、gear為分組變量，它們分別將樣本分為3、3、6組。
當(dāng)分組數(shù)量在三個及以上時，分組變量應(yīng)為因子類型。
2.2 平衡試驗設(shè)計
參照單因素方差分析的模型形式，可以很自然地推測雙因素方差分析的模型形式。
例2.1：觀察并比較如下兩組模型的結(jié)果
##?npk aov(yield?~?N?+?P,?npk) ##????????????????????????N????????P?Residuals ##?Sum?of?Squares??189.2817???8.4017??678.6817 ##?Deg.?of?Freedom????????1????????1????????21 aov(yield?~?P?+?N,?npk)? ##?Terms: ##????????????????????????P????????N?Residuals ##?Sum?of?Squares????8.4017?189.2817??678.6817 ##?Deg.?of?Freedom????????1????????1????????21##?mpg aov(mpg?~?cyl?+?gear,?mpg)? ##?Terms: ##??????????????????????cyl?????gear?Residuals ##?Sum?of?Squares??824.7846???8.2519??293.0107 ##?Deg.?of?Freedom????????2????????2????????27 aov(mpg?~?gear?+?cyl,?mpg)? ##?Terms: ##?????????????????????gear??????cyl?Residuals ##?Sum?of?Squares??483.2432?349.7933??293.0107 ##?Deg.?of?Freedom????????2????????2????????27上述兩組模型各包含兩個方差分析，模型形式的區(qū)別僅在于分組變量的順序不同：
對于第一組模型（npk）來說，兩個方差分析的結(jié)果完全相同；
對于第二組模型（mpg）來說，cyl和gear變量所對應(yīng)的離差平方和在兩個方差分析中是不同的，但二者之和仍然相同（824.7846 + 8.2519 = 483.2432 + 349.7933 = 833.0365），自由度和組內(nèi)離差平方和（Residuals）也相同。
以上結(jié)果意味著多因素方差分析可能與分組變量順序有關(guān)。實際上，這與樣本分組是否符合平衡試驗設(shè)計（balanced design）有關(guān)。平衡設(shè)計是指，由多個分組變量形成的交叉組所對應(yīng)的樣本數(shù)量完全一致。
對于平衡設(shè)計來說，多因素方差分析的結(jié)果與變量順序無關(guān)；而對于非平衡設(shè)計來說，存在三種類型的多因素方差分析：I型、II型和III型，而aov()函數(shù)使用的是I型方差分析，它的結(jié)果與變量順序有關(guān)。
在aov()函數(shù)的幫助文檔中作者寫道：aov()函數(shù)是為平衡設(shè)計準(zhǔn)備的，對于非平衡設(shè)計它的結(jié)果很難解釋。
統(tǒng)計兩個示例數(shù)據(jù)的交叉組樣本數(shù)量：
with(npk,?{tapply(yield,?list(N,?P),?length)}) ##???0?1 ##?0?6?6 ##?1?6?6with(mpg,?{tapply(mpg,?list(cyl,?gear),?length)}) ##????3??4?5 ##?4??1??8?2 ##?6??2??4?1 ##?8?12?NA?2可以看出，npk數(shù)據(jù)的4個交叉組的樣本均為6，屬于平衡設(shè)計；而mpg數(shù)據(jù)的9個交叉組的樣本并不完全一致，并且有一個交叉組沒有樣本（NA），屬于非平衡設(shè)計。
2.3 I型方差分析
對于平衡設(shè)計來說，三種類型的方差分析的結(jié)果一致，因此下文僅針對非平衡設(shè)計的方差分析，示例數(shù)據(jù)為mpg。
例2.2：aov()函數(shù)執(zhí)行的是I型方差分析，使用它運行如下模型
aov(mpg?~?cyl,?mpg) aov(mpg?~?gear,?mpg) aov(mpg?~?cyl?+?gear,?mpg) aov(mpg?~?gear?+?cyl,?mpg)為了方便比較，我們將模型的離差平方和匯總?cè)缦?#xff1a;
| cyl | 824.7846 | - | 301.2626 |
| gear | - | 483.2432 | 642.8040 |
| cyl + gear | 824.7846 | 8.2519 | 293.0107 |
| gear + cyl | 349.7933 | 483.2432 | 293.0107 |
SS表示離差平方和（Sum of Squares）。
通過比較可以得到如下結(jié)論：
所有模型的離差平方和之和（即總離差平方和）是相等的；
在雙因素方差分析中，分組變量的離差平方和與它的順序有關(guān)，其中第一個分組變量的離差平方和與它對應(yīng)的單因素方差分析的結(jié)果相同；
在雙因素方差分析中，cyl和gear的離差平方和與順序無關(guān)。
總離差平方和（）：
只與觀察變量有關(guān)，而與分組變量無關(guān)，其計算方法和單因素方差分析完全一致。使用表示樣本標(biāo)識，表示全體樣本的均值，有
組間離差平方和（和）：
這里組間離差平方和由兩部分組成，分別對應(yīng)兩個分組變量的離差平方和，記為、。使用、表示兩個分組變量的水平數(shù)，、表示分組標(biāo)識。
假設(shè)對應(yīng)的是第一個分組變量的離差平方和，我們已經(jīng)知道它和單因素方差分析的結(jié)果一樣：
其中，表示第組的樣本均值。
不易直接計算，但是我們知道它與之和為定值：
因此只要計算出即可反推出。
組內(nèi)離差平方和（）：
相當(dāng)于“剩余”部分，在輸出結(jié)果中也使用Residuals標(biāo)識。方差分析與線性回歸存在內(nèi)在聯(lián)系聯(lián)系，實際就等于線性回歸的殘差平方和。
fit.lm?<-?lm(mpg?~?cyl?+?gear,?mpg) sum(residuals(fit.lm)^2) ##?[1]?293.0107自由度和F統(tǒng)計量：
的自由度為（為樣本總量），、的自由度分別為和，的自由度為。
2.4 交互效應(yīng)
分組變量之間可能存在交互效應(yīng)。
例2.3：比較不含和含有交互項的雙因素方差分析
aov(mpg?~?cyl?+?gear,?mpg)? ##?Terms: ##??????????????????????cyl?????gear?Residuals ##?Sum?of?Squares??824.7846???8.2519??293.0107 ##?Deg.?of?Freedom????????2????????2????????27aov(mpg?~?cyl*gear,?mpg) ##?Terms: ##??????????????????????cyl?????gear?cyl:gear?Residuals ##?Sum?of?Squares??824.7846???8.2519??23.8907??269.1200 ##?Deg.?of?Freedom????????2????????2????????3????????24在有交互項的情況下，組間離差平方和由三部分組成：、、，其中表示交互項的離差平方和。
可以看出，交互項的出現(xiàn)不影響、、的值。那的值和自由度應(yīng)如何確定呢？
我們依然可以先使用線性回歸求出再反推：
fit2.lm?<-?lm(mpg?~?cyl*gear,?mpg) sum(residuals(fit2.lm)^2) ##?[1]?269.12不過也可以直接求出、、的和：
其中，表示交叉組的樣本均值。
將交叉組看作是由單個變量分組形成的，也可以使用類似單因素方差分析的方法直接計算：
假設(shè)交叉組的個數(shù)為，則、、的自由度之和為，因此的自由度為。
顯然，交叉組的最大數(shù)目為，在此情況下：
前文已經(jīng)提到，mpg數(shù)據(jù)中有一個交叉組不存在，因此，而對于平衡設(shè)計來說，總是成立的。
例2.4：手動計算和
data?<-?mpg?%>%mutate(mean?=?mean(mpg))?%>%group_by(cyl,?gear)?%>%mutate(mean2?=?mean(mpg))?%>%ungroup()?%>%mutate(SSAB?=?(mean2?-?mean)^2,SSE?=?(mpg?-?mean2)^2)sum(data$SSAB)?-?824.7846?-?8.2519? ##?[1]?23.89069 sum(data$SSE) ##?[1]?269.122.5 多因素方差分析
從前文可以看出，即使只有兩個分組變量，I型方差分析的都很難直接計算。那對于多因素方差分析來說該如何計算各個離差平方和呢？下面通過兩個例子進(jìn)行說明。
例2.5：運行形式最簡單的三因素方差分析
aov(mpg?~?cyl?+?gear?+?carb,?mpg) ##?Terms: ##??????????????????????cyl?????gear?????carb?Residuals ##?Sum?of?Squares??824.7846???8.2519??88.5199??204.4908 ##?Deg.?of?Freedom????????2????????2????????5????????22cyl的離差平方和：
fit1?<-?lm(mpg?~?cyl,?mpg)? sum((fitted(fit1)?-?mean(mpg$mpg))^2) ##?[1]?824.7846gear的離差平方和：
fit2?<-?lm(mpg?~?cyl?+?gear,?mpg)? sum((fitted(fit2)?-?mean(mpg$mpg))^2)?-?sum((fitted(fit1)?-?mean(mpg$mpg))^2) ##?[1]?8.251855carb的離差平方和：
fit3?<-?lm(mpg?~?cyl?+?gear?+?carb,?mpg)? sum((fitted(fit3)?-?mean(mpg$mpg))^2)?-?sum((fitted(fit2)?-?mean(mpg$mpg))^2) ##?[1]?88.5199Residuals的離差平方和：
sum(residuals(fit3)^2) ##?[1]?204.4908因此，I型方差分析可以使用逐步回歸的方法計算離差平方和。I型平方和又稱為順序平方和。
例2.6：運行含有交互項的三因素方差分析
aov(mpg?~?cyl*gear?+?carb,?mpg) ##?Terms: ##??????????????????????cyl?????gear?????carb?cyl:gear?Residuals ##?Sum?of?Squares??824.7846???8.2519??88.5199??25.3878??179.1030 ##?Deg.?of?Freedom????????2????????2????????5????????2????????20aov(terms(mpg?~?cyl*gear?+?carb,?keep.order?=?T),?mpg) ##?Terms: ##??????????????????????cyl?????gear?cyl:gear?????carb?Residuals ##?Sum?of?Squares??824.7846???8.2519??23.8907??90.0170??179.1030 ##?Deg.?of?Freedom????????2????????2????????3????????4????????20默認(rèn)情況下，交互項會排在所有主效應(yīng)變量之后；而使用terms()函數(shù)可以更改設(shè)置。位置不同會導(dǎo)致離差平方和發(fā)生變化。
兩種情況下cyl:gear的離差平方和：
##?默認(rèn)情況 fit11?<-?lm(mpg?~?cyl?+?gear?+?carb,?mpg)? fit22?<-?lm(mpg?~?cyl*gear?+?carb,?mpg)? sum((fitted(fit22)?-?mean(mpg$mpg))^2)?-?sum((fitted(fit11)?-?mean(mpg$mpg))^2)? ##?[1]?25.38784##?更改后的情況 fit33?<-?lm(mpg?~?cyl?+?gear,?mpg)? fit44?<-?lm(mpg?~?cyl*gear,?mpg)? sum((fitted(fit44)?-?mean(mpg$mpg))^2)?-?sum((fitted(fit33)?-?mean(mpg$mpg))^2)? ##?[1]?23.89074正如aov()函數(shù)的幫助文檔所言，非平衡設(shè)計的I型方差分析的結(jié)果很難解釋。學(xué)堂君注意到<醫(yī)學(xué)方>公眾號的一篇推文對其做了非常直觀的解釋，有興趣的讀者可以閱讀：
關(guān)于方差分析，您必須知道的四種類型平方和
總結(jié)
以上是生活随笔為你收集整理的R语言基础 | 方差分析（2）：多因素方差分析（上）的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： 硬件nat关闭还是开启_潜在的NAT加速
- 下一篇： zabbix添加微信报警