一文展示R語言中的方差分析常用模型 _#2021.9.11

方差分析是一個全新的思路，它采用的是變異分解的思路，將組內(nèi)組件分開，查看顯著性。

變異分解，和數(shù)量遺傳學(xué)的創(chuàng)立也密不可分，比如

下面看我的總結(jié)：

1. 方差分析的假定

上面這個思維導(dǎo)圖，也可以看出，方差分析有三大假定：正態(tài)，獨(dú)立和齊次，如果不滿足，可以使用廣義線性模型或者混合線性模型，或者廣義線性混合模型去分析。

「本次我們的主題有：」

2. 數(shù)據(jù)來源

這里，我們使用的數(shù)據(jù)來源于R包agridat，它是講農(nóng)業(yè)相關(guān)的論文，書籍中相關(guān)的數(shù)據(jù)收集在了一起，更加符合我們的背景。

包的下載地址：https://cran./web/packages/agridat/index.html

「包的介紹」

?

「agridat: Agricultural Datasets」Datasets from books, papers, and websites related to agriculture. Example graphics and analyses are included. Data come from small-plot trials, multi-environment trials, uniformity trials, yield monitors, and more.

?

「包的安裝方式：」

3. 單因素方差分析

?

比如一個處理有A，B，C三個水平，觀測值y，想看一下這個處理是否達(dá)到顯著性水平，這就可以用到方差分析了。

?

「數(shù)據(jù)描述：」

?

Corn yield and nitrogen fertilizer treatment with field characteristics for the Las Rosas farm, Rio Cuarto, Cordoba, Argentina.

?

「數(shù)據(jù)結(jié)構(gòu)：」

這里數(shù)據(jù)有很多列，但是我們要演示單因素方差分析，這里的因素為nf，自變量（Y變量）是yield，想要看一下nf的不同水平是否達(dá)到顯著性差異。

「建模：」

Y變量：yield

因子：nf

「R中的建模代碼：」

• m1為模型保存的名稱
• aov為R中的方差分析代碼
• yield為數(shù)據(jù)中的Y變量，這里是yield
• ~，波浪號前面為Y變量，后面為X變量
• nf為分析的因子變量
• dat為數(shù)據(jù)

「結(jié)果：」

結(jié)果可以看出，nf之間達(dá)到極顯著水平，可以進(jìn)行多重比較。

「方差分析的顯著性和多重比較有何關(guān)系？？？」

?

方差分析，一般會得到顯著性（即P值小于0.05，說明顯著，小于0.01，說明極顯著，大于0.05，說明不顯著），顯著的意思是因素之間至少有一對水平達(dá)到顯著性差異，具體是那一對呢？有幾對呢？這就需要用到多重比較。所以，多重比較是在方差分析達(dá)到顯著性之后進(jìn)行的，只有顯著了（P值小于0.05）才有能進(jìn)行多重比較。多重比較的方法有很多，比如LSD，Tukey，Bonferroni，Holm等等，我們后面系統(tǒng)的介紹。

?

4. 單因素隨機(jī)區(qū)組

這里區(qū)組的設(shè)置，可以控制一些環(huán)境誤差。

「數(shù)據(jù)描述：」

?

Switchback trial in dairy with three treatments

?

「數(shù)據(jù)結(jié)構(gòu)：」

「建模：」

Y變量：yield

因子：trt

區(qū)組：block

「R中的建模代碼：」

「結(jié)果：」

結(jié)果可以看出，trt之間達(dá)到?jīng)]有極顯著水平。

5. 二因素方差分析（無交互）

這里區(qū)組的設(shè)置，可以控制一些環(huán)境誤差。

無交互的意思是，二因素之間不考慮互作。

「數(shù)據(jù)描述：」

?

Switchback trial in dairy with three treatments

?

「數(shù)據(jù)結(jié)構(gòu)：」

「建模：」

• Y變量：yield
• 因子1：trt
• 因子2：period
• 區(qū)組：block

「R中的建模代碼：」

「結(jié)果：」

結(jié)果可以看出，trt之間達(dá)到?jīng)]有極顯著水平，period之間沒有達(dá)到顯著性水平。

6. 二因素方差分析（有交互）

這里區(qū)組的設(shè)置，可以控制一些環(huán)境誤差。

交互的意思是，二因素之間考慮互作。

「數(shù)據(jù)描述：」

?

Switchback trial in dairy with three treatments

?

「數(shù)據(jù)結(jié)構(gòu)：」

「建模：」

• Y變量：yield
• 因子1：trt
• 因子2：period
• 區(qū)組：block

「R中的建模代碼：」

注意，這里的trt*period是R中公式的簡寫，表示trt + period + trt:period，其中trt：period表示互作效應(yīng)。

「結(jié)果：」

結(jié)果可以看出，trt之間達(dá)到?jīng)]有極顯著水平，period之間沒有達(dá)到顯著性水平，trt和period交互沒有達(dá)到顯著水平。

7. 多因素方差分析

這里區(qū)組的設(shè)置，可以控制一些環(huán)境誤差。

多于兩個因素的方差分析

「數(shù)據(jù)描述：」

?

Switchback trial in dairy with three treatments

?

「數(shù)據(jù)結(jié)構(gòu)：」

「建模：」

• Y變量：yield
• 因子1：trt
• 因子2：period
• 因子3：cow
• 區(qū)組：block

「R中的建模代碼：」

注意，這里的trt*period是R中公式的簡寫，表示trt + period + trt:period，其中trt：period表示互作效應(yīng)。

「結(jié)果：」

結(jié)果可以看出，trt之間達(dá)到極顯著水平，period之間達(dá)到顯著性水平，trt和period交互沒有達(dá)到顯著水平，cow達(dá)到極顯著水平。

8. 裂區(qū)試驗方差分析

裂區(qū)試驗，包括主區(qū)（wplot）和裂區(qū)（splot），其中裂區(qū)是鑲嵌在主區(qū)中的，主區(qū)和裂區(qū)的殘差不一樣，在模型中需要特殊定義。

「數(shù)據(jù)描述：」

?

Strip-split plot of barley with fertilizer, calcium, and soil factors.

?

「數(shù)據(jù)結(jié)構(gòu)：」

「建模：」

• Y變量：yield
• 主區(qū)：fert
• 裂區(qū)：soil
• 區(qū)組：brep

「R中的建模代碼：」

注意，這里的Error(rep/fert)是R中公式的定義殘差，主要用于不同因素使用不同殘差的情況，這里fert是主區(qū)。

「結(jié)果：」

結(jié)果可以看出，soil達(dá)到極顯著，fert不顯著，soil和fert的互作不顯著。

9. 正態(tài)性檢驗

方差分析中，結(jié)果是否可信，在于數(shù)據(jù)是否滿足假定條件。方差分析的假定包括數(shù)據(jù)正態(tài)性，數(shù)據(jù)的方差齊性，數(shù)據(jù)的獨(dú)立性，其中可以檢驗的假定有：

• 數(shù)據(jù)的正態(tài)性
• 數(shù)據(jù)的齊性

這里，我們介紹如何對數(shù)據(jù)的正態(tài)性進(jìn)行檢驗。

可以使用球性檢驗（Shapiro-Wilk）檢驗數(shù)據(jù)的正態(tài)性，也可以用qqplot查看殘差的圖，判斷數(shù)據(jù)的正態(tài)性，也可以對數(shù)據(jù)做直方圖，查看數(shù)據(jù)的正態(tài)性。

「數(shù)據(jù)描述：」

?

A classical N, P, K (nitrogen, phosphate, potassium) factorial experiment on the growth of peas conducted on 6 blocks. Each half of a fractional factorial design confounding the NPK interaction was used on 3 of the plots.

?

「數(shù)據(jù)結(jié)構(gòu)：」

一般分析時，我們僅對Y變量進(jìn)行正態(tài)性檢驗，如果是單因素或者多因素，也可以根據(jù)因素分組進(jìn)行正態(tài)性檢驗，數(shù)據(jù)量大時，對于稍微偏態(tài)的數(shù)據(jù)，即使不太符合正態(tài)分布，也不影響結(jié)論。

這里，我們對yield進(jìn)行正態(tài)性檢驗

「作直方圖」

可以看到，yield大體符合正態(tài)分布。

「做qq圖」

使用car軟件包中的qqPlot函數(shù)。

可以看到，數(shù)據(jù)基本落在置信區(qū)間里面，數(shù)據(jù)符合正態(tài)分布。

「使用Shapiro-Wilk」檢驗數(shù)據(jù)正態(tài)分布

可以看到，P值為0.8735，數(shù)據(jù)符合正態(tài)分布，與上圖顯示的結(jié)果一致。

10. 齊性檢驗

方差分析中，我們對結(jié)果是否自信，在于數(shù)據(jù)是否滿足假定條件，方差分析的假定條件包括數(shù)據(jù)正態(tài)性，數(shù)據(jù)的方差齊性，數(shù)據(jù)的獨(dú)立性，其中可以檢驗的假定有：

• 數(shù)據(jù)的正態(tài)性
• 數(shù)據(jù)的齊性

這里，我們介紹如何對數(shù)據(jù)的齊性進(jìn)行檢驗。齊性檢驗，是檢驗不同樣本的總體方差是否相同，是根據(jù)特定的模型，需要考慮不同的因子放到模型中，不能單獨(dú)對某一列變量進(jìn)行齊性檢驗。

「齊性檢驗：」

• Bartlet檢驗
• Levene檢驗

「數(shù)據(jù)描述：」

?

A classical N, P, K (nitrogen, phosphate, potassium) factorial experiment on the growth of peas conducted on 6 blocks. Each half of a fractional factorial design confounding the NPK interaction was used on 3 of the plots.

?

「數(shù)據(jù)結(jié)構(gòu)：」

比如上面數(shù)據(jù)中，相對N進(jìn)行單因素方差分析，查看不同N的水平是否滿足方差齊性，可以這樣操作：

「Bartlett檢驗」

Bartlett檢驗可以比較多個總體的方差

結(jié)果：

結(jié)果可以看出，不同的N之間，方差滿足齊性要求。

「Levene檢驗」

Bartlett檢驗對數(shù)據(jù)的正態(tài)性非常敏感，而Levene檢驗是一種非參數(shù)檢驗方法，使用范圍更廣。

結(jié)果：

結(jié)果可以看出，P值為0.9421，大于0.05，說明數(shù)據(jù)符合方差齊性。

11. 多重比較

這里，我們介紹一下方差分析中的多重比較方法。

「agricolae包」

?

This package contains functionality for the Statistical Analysis of experimental designs applied specially for field experiments in agriculture and plant breeding.

?

「多重比較方法：」

• LSD
• scheffe
• HSD（Tukey）
• SNK
• Duncan

「數(shù)據(jù)描述：」

?

A classical N, P, K (nitrogen, phosphate, potassium) factorial experiment on the growth of peas conducted on 6 blocks. Each half of a fractional factorial design confounding the NPK interaction was used on 3 of the plots.

?

「數(shù)據(jù)結(jié)構(gòu)：」

「方差分析」這里對N進(jìn)行單因素方差分析，block為區(qū)組：

結(jié)果可以看到，N因素達(dá)到極顯著水平

11.1 LSD多重比較

?

Multiple comparisons of treatments by means of LSD and a grouping of treatments. The level by alpha default is 0.05. Returns p-values adjusted using one of several methods

?

結(jié)果：

11.2 scheffe多重比較

?

Scheffe 1959, method is very general in that all possible contrasts can be tested for significance and confidence intervals can be constructed for the corresponding linear. The test is conservative.

?

代碼：

結(jié)果

11.3 Tukey多重比較

?

It makes multiple comparisons of treatments by means of Tukey. The level by alpha default is 0.05.

?

代碼：

結(jié)果

11.4 SNK多重比較

?

SNK is derived from Tukey, but it is less conservative (finds more differences). Tukey controls the error for all comparisons, where SNK only controls for comparisons under consideration. The level by alpha default is 0.05.

?

代碼：

結(jié)果

11.5 Duncan多重比較

?

This test is adapted from the Newman-Keuls method. Duncan's test does not control family wise error rate at the specified alpha level. It has more power than the other post tests, but only because it doesn't control the error rate properly. The Experimentwise Error Rate at: 1-(1-alpha)^(a-1); where "a" is the number of means and is the Per-Comparison Error Rate. Duncan's procedure is only very slightly more conservative than LSD. The level by alpha default is 0.05.

?

代碼：

結(jié)果

12. 獲得所有代碼及注釋腳本