寫在前言

關於方差分析（Analysis of Variance，ANOVA），知乎的大神們都已經科普過其概念，簡單來說就是檢驗多組樣本之間均值的差異。而重複測量方差分析，顧名思義即多次測量的數據來自於同一個體，這些數據之間可能存在相關性，與方差分析的基本假設——數據獨立性相違背，從而導致常規的方差分析會大大增加犯I類錯誤的概率。

隨著R的使用越來越普遍，世界各地的研究者們根據不同的統計模型相繼開發出各種各樣的載入包，使得R在數據處理和分析的功能上越來越強大。然而，初學R的朋友可能不熟悉各種函數的定義，對於該用哪個函數去分析自己的數據就成了一個非常頭痛的問題。

僅應用於方差分析，最基本的就有aov( )，anova( )，Anova( ){car包}，這三種函數使用情況和公式都不一樣，各有利弊。關於它們的使用和對比也有很多人介紹過，今天要介紹的是另一個我比較常用的包——『ez』。

這個包在國外也越來越多人使用，正如其名，使用起來非常的easy，相比anova( ) 以及Anova( )需要先用lm( )定義線性模型，或者aov( )需要更改平方和的計算類型（Type I, II, III），ezANOVA( )的運算簡化到基本一個語句就搞定所有，並且輸出結果更為直觀。

詳見：ez package

示例使用的數據是來自某研究關於文章的主題熟悉性與生字密度對學生閱讀理解的影響。

  • reading.xlsx 提取碼：854g

與SPSS要求的數據格式wide-data不同，『ez』包對數據格式的要求是long-data，對於用慣SPSS的朋友們來說，相當於數據處理只做到Aggregate這一步就可以了。

> library(readxl) # 讀取數據
> data <- read_excel(file.choose())

> data$Subject <- as.factor(data$Subject)
> data$Familiarity <- as.factor(data$Familiarity)
> data$Density <- as.factor(data$Density)

> attach(data) # 將數據框的變數（列）直接添加到R的搜索路徑，調用時就不用加上$

首先對數據進行一些簡單的預覽，執行str( )函數我們可以了解數據包含的變數及其類型：

> str(data)
Classes 『tbl_df』, 『tbl』 and data.frame: 60 obs. of 4 variables:
$ Subject : Factor w/ 10 levels "1","2","3","4",..: 1 1 1 1 1 1 2 2 2 2 ...
$ Familiarity: Factor w/ 2 levels "familiar","unfamiliar": 2 2 2 1 1 1 2 2 2 1 ...
$ Density : Factor w/ 3 levels "high","low","mid": 1 3 2 1 3 2 1 3 2 1 ...
$ Score : num 3 4 5 4 8 12 6 6 7 5 ...

其中兩個被試內變數分別為主題熟悉性（Familiarity: unfamiliar vs. familiar）和生字密度（Density: high vs. mid vs. low），因變數為閱讀分數（Score）。

*執行head( )函數後我們可以看到其格式如下：

> head(data)
# A tibble: 6 x 4
Subject Familiarity Density Score
1 1 unfamiliar high 3
2 1 unfamiliar mid 4
3 1 unfamiliar low 5
4 1 familiar high 4
5 1 familiar mid 8
6 1 familiar low 12

了解完數據以後就可以開始我們的分析了。

> library(ez) # 載入包

分析使用的函數是ezANOVA( data, dv, wid, within = NULL, within_full = NULL, within_covariates = NULL, between = NULL, between_covariates = NULL, observed = NULL, diff = NULL, reverse_diff = FALSE, type = 2, white.adjust = FALSE, detailed = FALSE, return_aov = FALSE )，在此例子中需要聲明的變數如下：

data, 需要分析的數據集，如前文所述格式要求為long-data

dv, 研究中的因變數（dependent variable）wid, 數據集中的被試編號（case）

within, 數據集中的被試內自變數，如果只有一個時可以直接聲明，有多變數的情況

下其格式為 within = .(FactorA, FactorB) *注意括弧前面有個英文的句號type, 默認的計算類型的II型，只考慮自變數的主效應，不考慮交互作用。如果數 據是平衡的，且不存在交互作用，則Type I，II和III得出的結果是相同的。

detailed, 默認是FALSE，如果選擇TRUE，則返回自變數的error及Sum squares（用

於後期計算）

> model <- ezANOVA(data, dv = Score, wid = Subject, within = .(Familiarity, Density), type = 3, detailed = T)
> model
$`ANOVA`
Effect DFn DFd SSn SSd F p p<.05 ges
1 (Intercept) 1 9 2136.0667 37.60000 511.29255 3.067983e-09 * 0.9605192
2 Familiarity 1 9 166.6667 25.66667 58.44156 3.176033e-05 * 0.6549646
3 Density 2 18 176.5333 9.80000 162.12245 3.079071e-12 * 0.6678436
4 Familiarity:Density 2 18 102.9333 14.73333 62.87783 7.565333e-09 * 0.5396714

$`Mauchlys Test for Sphericity`
Effect W p p<.05
3 Density 0.6462585 0.1744316
4 Familiarity:Density 0.3596364 0.0167284 *

$`Sphericity Corrections`
Effect GGe p[GG] p[GG]<.05 HFe p[HF] p[HF]<.05
3 Density 0.7386935 1.537778e-09 * 0.8490314 1.111375e-10 *
4 Familiarity:Density 0.6096209 4.014620e-06 * 0.6549759 1.928903e-06

輸出結果包含三個表，第一個表ANOVA即方差分析的結果，其中ges（Generalized Eta-Squared）表示的eta^{2}，是SPSS中所沒有輸出的。

後面兩個表依次為球形檢驗的結果及其校正，由於生字密度和交互作用的水平數都大於2，需要進行球形檢驗。當球形檢驗的p 值大於0.05時，自變數各水平之間都符合球形假設，即相互獨立。反之則表明不符合球形假設，需要參考第三個表Sphericity Corrections中的GGe（Greenhouse-Geisser epsilon）或者HFe（Huynh-Feldt epsilon）對df 進行校正。關於採取哪種校正方案， Girden (1992)認為當epsilon > 0.75時，建議用Huynh-Feldt校正，當epsilon < 0.75時，建議用Greenhouse-Geisser校正。

Mauchlys Test表的結果顯示，只有交互作用不滿足球形假設，為此需要根據Sphericity Corrections 中的GGe對交互作用的df 進行校正：

> Int_DFn <- 2 * model$`Sphericity Corrections`[2,2]
> Int_DFd <- 18 * model$`Sphericity Corrections`[2,2]
> table(int_DFn,int_DFd)
int_DFn int_DFd
1.2192418 10.9731762

前文提到ezANOVA默認輸出的是 eta^{2} ，但有些期刊要求報告 eta_{p}^{2} ，在多變數模型中表示排除其他變數後自變數對因變數單獨的影響。根據eta_{p}^{2}的計算公式，我們將ANOVA列表中的SSn / SSn + SSd 即可得到eta_{p}^{2}的值（對應前文公式中提到的detailed）。

> model$ANOVA[4] / (model$ANOVA[4] + model$ANOVA[5])
SSn
1 0.9827020
2 0.8665511
3 0.9474061
4 0.8747875

到這一步我們的重複測量方差分析就做完了，從輸出的結果中我們可以看到，主題熟悉性（F (1, 9) = 58.44, p < 0.001, eta_{p}^{2}= 0.867）和生字密度（F (2, 18) = 162.12, p < 0.001, eta_{p}^{2} = 0.947）的主效應以及交互作用（F (1.22, 10.97) = 62.88, p < 0.001, eta_{p}^{2} = 0.875）均顯著，接下來就主要對交互作用進行分析了。具體可以參考網上的其它教程，這裡就不再多說。

結果顯著了可以發文章了噢耶！！！！

總的來說，ezANOVA( )的使用方法較為直觀，其輸出結果也與SPSS相似，輸出內容除了包括常規的F值，p值，還可以根據分析的結果調整自變數的df，並求出 eta_{p}^{2}，基本涵蓋報告數據結果時所需要的所有指標。相比anova( )，Anova( ){car}以及aov( )，個人更傾向使用ezANOVA做重複測量方差分析。

第一次分享學習的經驗，歡迎對R有興趣的小夥伴一起交流探討，如上述內容存在遺漏或者錯誤，還請多多批評指教。小小做個預告，下次還是用這個例子分享一下如何用『BayesFactor』這個包計算貝葉斯因子。

作者聲明：

本文原創性歸作者所有，轉載請註明來源和作者。內容只用於交流與學習，禁止用於任何商業用途。

推薦閱讀：

相关文章