SPSS：單因素重複測量方差分析

知識 08-17

一、問題與數據

研究者想知道鍛煉是如何降低心臟病發生風險的，於是他關注了一種與心臟病相關的蛋白——C反應蛋白（C-Reactive Protein, CRP）：CRP的濃度越高，發生心臟病的風險越高。研究者想知道鍛煉對CRP濃度的影響。

研究者招募了10名研究對象，研究對象進行了6個月的鍛煉干預。CRP濃度共測量了3次：干預前的CRP濃度——crp_pre；干預中（3個月）——crp_mid；干預後（6個月）——crp_post。這三個時間點代表了受試者內因素「時間」的三個水平，因變數是CRP的濃度，單位是mg/L。部分數據如下：

二、對問題的分析

使用One-way Repeated Measures Anova進行分析時，需要考慮6個假設。

對研究設計的假設：

假設1：因變數唯一，且為連續變數；

假設2：受試者內因素（Within-Subject Factor）有3個或以上的水平。

註：在重複測量的方差分析模型中，對同一個體相同變數的不同次觀測結果被視為一組，用於區分重複測量次數的變數被稱為受試者內因素，受試者內因素實際上是自變數。

對數據的假設：

假設3：受試者內因素的各個水平，因變數沒有極端異常值；

假設4：受試者內因素的各個水平，因變數需服從近似正態分布；

假設5：對於受試者內因素的各個水平組合而言，因變數的方差協方差矩陣相等，也稱為球形假設。

三、思維導圖

四、對假設的判斷

在分析時，如何考慮和處理這5個假設呢？

由於假設1-2都是對研究設計的假設，需要研究者根據研究設計進行判斷，所以我們主要對數據的假設3-5進行檢驗。

(一)檢驗假設3和假設4的SPSS操作

1.在主菜單點擊Analyze > Descriptive Statistics > Explore...，如下圖：

2.出現Explore對話框，將crp_pre、crp_mid和crp_post選入Dependent List，點擊Plots；

3.出現下圖Plots對話框；

4.在Boxplots下選擇Dependents together，去掉Descriptive下Stem-和-leaf，選擇Normality plots with tests，點擊Continue；

5.回到Explore主對話框，在Display下方選擇Plots，點擊OK。

(二)檢驗假設3：受試者內因素的各個水平，因變數沒有極端異常值

1.在輸出的結果中，如下圖所示，在箱式圖中，距離箱子邊緣超過1.5倍箱身長度的數據點定義為異常值，在本例中，未發現異常值。

2.為了方便進一步的理解，下面圖示是有異常值的箱式圖，上下邊緣超過1.5倍箱式長度為異常值，如下圖所示，用「圓圈」表示，右上標為異常值在數據表中所對應的行數；將距離箱子邊緣超過3倍箱身長度的數據點定義為極端值（極端異常值），用「*」表示，右上標代表異常值在數據表中所對應的行數。因此，在箱式圖中查看異常值時，可以直接看「圓圈」或「*」。

3. 異常值的處理

導致數據中存在異常值的原因有3種：

(1) 數據錄入錯誤：首先應該考慮異常值是否由於數據錄入錯誤所致。如果是，用正確值進行替換並重新進行檢驗；

(2) 測量誤差：如果不是由於數據錄入錯誤，接下來考慮是否因為測量誤差導致（如儀器故障或超過量程），測量誤差往往不能修正，需要把測量錯誤的數據刪除；

(3) 真實存在的異常值：如果以上兩種原因都不是，那最有可能是一種真實的異常數據。這種異常值不好處理，但也沒有理由將其當作無效值看待。目前它的處理方法比較有爭議，尚沒有一種特別推薦的方法。

需要注意的是，如果存在多個異常值，應先把最極端的異常值去掉後，重新檢查異常值情況。這是因為有時最極端異常值去掉後，其他異常值可能會回歸正常。

異常值的處理方法分為2種：

(1) 保留異常值：

採用非參數Friedman test檢驗；

用非最近端的值代替極端異常值（如用第二大的值代替極端異常值）；

因變數轉換成其他形式；

將異常值納入分析，並堅信其對結果不會產生實質影響。

(2) 剔除異常值：

直接刪除異常值很簡單，但卻是沒有辦法的辦法。當我們需要刪掉異常值時，應報告異常值大小及其對結果的影響，最好分別報告刪除異常值前後的結果。而且，應該考慮有異常值的個體是否符合研究的納入標準。如果其不屬於合格的研究對象，應將其剔除，否則會影響結果的推論。

（更多閱讀：怎麼判別我的數據中存在特異值？教你幾招！）

(三) 檢驗假設4：受試者內因素的各個水平，因變數需服從近似正態分布

1.對於樣本量較小（

2.如果數據不服從正態分布，可以有如下4種方法進行處理：

(1) 數據轉換：對轉換後呈正態分布的數據進行單因素方差分析。當各組因變數的分布形狀相同時，正態轉換才有可能成功。對於一些常見的分布，有特定的轉換形式，但是對於轉換後數據的結果解釋可能比較複雜。

(2) 使用非參數檢驗：可以使用Friedman test等非參數檢驗方法，但是要注意Friedman test和單因素重複測量方差分析的無效假設和備擇假設不太一致。

(3) 直接進行分析：由於單因素重複測量方差分析對於偏離正態分布比較穩健，尤其是在各組樣本量相等或近似相等的情況下，而且非正態分布實質上並不影響犯I型錯誤的概率。因此可以直接進行檢驗，但是結果中仍需報告對正態分布的偏離。

(4) 檢驗結果的比較：將轉換後和未轉換的原始數據分別進行單因素重複測量方差分析，如果二者結論相同，則再對未轉換的原始數據進行分析。

（更多閱讀：SPSS教程：判斷數據正態分布的超多方法！；SPSS詳細操作：正態轉換的多種方法）

五、SPSS操作

(一) 單因素重複測量方差分析的操作

1.在主菜單下點擊Analyze > General Linear Model > Repeated measures...，如下圖所示：

2.出現Repeated Measures Define Factor(s)對話框，如下圖所示：

3.在Within-Subject Factor Name:中將「factor1」更改為time，因為共測量了3次CRP的水平，所以在Number of Levels：中填入3；

4.點擊Add，出現下圖：

5.在Measure Name:中賦予因變數一個合理的名字，本例中因變數為CRP的水平，所以填入CRP，點擊下方的Add，出現下圖：

6.點擊Define，出現下圖Repeated Measures對話框：

7.如下圖所示，Within-Subjects Variables後面的括弧內是受試者內因素的名字，框內三條分別代表「time」對應的各個水平；

8.將crp_pre、crp_mid和crp_post一起選入右側的框中，如下圖所示：

9.點擊Plots，出現Repeated Measures: Profile Plots 對話框，如下圖所示：

10.將time選入Horizontal Axis：框中；

11.點擊Add，出現下圖，點擊Continue；

12.點擊Options，出現Repeated Measures: Options對話框；

13.將time選入Display Means for：中；下方Compare main effects為勾選狀態；在Confidence interval adjustment：下選擇Bonferroni；在Display下方勾選Descriptive statistics 和Estimates of effect size；點擊Continue，點擊OK。

(二)單因素重複測量方差分析自定義組間比較（custom contrasts）

如果只關心特定組別間的差異，你需要知道如何進行自定義比較（custom contrasts），以及如何對多重比較結果進行調整，這就要用到SPSS軟體中的Syntax Editor窗口編寫相應程序語句。當滿足方差齊性條件時，推薦採用單因素重複測量方差分析進行自定義組間比較。

1.操作同上述1-13步的操作，這裡不再贅述。接下來：

(1) 點擊Paste，出現IBM SPSS Statistics Syntax Editor窗口：

(2) 在 /PRINT 和 /CRITERIA兩行中間，輸入/MMATRIX = "0 vs 6 Months" ALL 1 0 -1;

/MMATRIX=旨在告訴SPSS我們要做一個自定義假設； 1 0 -1表示要進行比較的自變數組別，組別的順序和SPSS里輸入的組別順序有關：這裡從左到右（1 0 -1）分別對應著crp_pre、crp_mid和crp_post，表示將crp_pre與crp_post進行比較。「ALL」代表受試者內因素的所有水平。

註：自定義比較包括了簡單比較（simple contrasts）和複合比較（complex contrasts）。簡單比較為只比較自變數某兩個組別間的差異，需要建立線性比較函數（linear contrast，φ）。它包含一系列組別和每個組別對應的均數，組別取值只能為1，-1，0。我們把要比較的兩組的組別分別賦值為1和-1，其他不比較的組別組別賦值為0。

本例中crp_pre組別為1，crp_post組別為-1，其他組別均為0，則是要比較干預前CRP水平和干預後CRP水平的差異，看二者的平均CRP水平的差值是否為0（即組別為1的組別減去組別為-1的組別，以組別為-1的組別為參照組，組別賦值的正負與研究設計和研究假設有關）。

(3) 用/MMATRIX指令增加另外比較：

/MMATRIX = "0 vs 3&6 Months" ALL 1 -1/2 -1/2

註：複合比較為比較自變數超過2個組別的組合間的差異。同樣採用線性比較函數的方法，某組合的組別賦值為1或-1除以組合內的組數，但是要保證要比較的組間組合與另一組（組合）的所有組別加起來為0，組別賦值的正負與研究設計和研究假設有關。

本例中， /MMATRIX = "0 vs 3&6 Months" ALL 1 -1/2 -1/2表示crp_pre與crp_mid、crp_post差異的比較。

(4) 多重比較的校正

接下來，我們需要校正顯著性水平（α），通常也可以校正每次比較的P值和可信區間，得到調整後P值和聯合可信區間（simultaneous confidence intervals）。我們首先採用Bonferroni方法對顯著性水平α進行校正，公式如下：

調整後α=調整前α ÷ 比較的次數

本例中我們需要進行2次比較，則調整後α=0.05÷2=0.025。

(5) 箭頭標註處為SPSS軟體默認的顯著性水平α=0.05：

/CRITERIA=ALPHA(.05)

我們將其改為調整後的顯著性水平α=0.025：

/CRITERIA=ALPHA(.025)

(6) 在菜單欄點擊Run > All。

六、結果解釋

1. 基本描述

(1) SPSS首先會給出Within-Subjects Factors表，該表提示了受試者內因素crp_pre、crp_mid和crp_post對應的標籤為1、2和3，在後面的Estimates表和Pairwise Comparisons表中會用到。

(2) Descriptive Statistics表給出了crp_pre、crp_mid和crp_post的均值、標準差和例數。受試者干預前、干預中和干預後的CRP濃度分別為4.33 ± 0.64 mg/L、3.94 ± 0.57 mg/L和3.65 ± 0.43mg/L。

(3) Estimates表中沒有再出現crp_pre、crp_mid和crp_post的變數名，而是給出了對應的3個時間點的標籤。該表中給出了crp_pre、crp_mid和crp_post的均值、標準誤和95%的置信區間。

(4) Estimated Marginal Means of CRP給出了三個時間點CRP的均值的折線圖，可以看到從干預前到干預後呈下降趨勢，但1到2的時間點下降程度比2到3的下降程度大。

2. 球形假設的檢驗結果

(1) 在Mauchly s Test of Sphericity表中，給出了球形假設的檢驗結果。如果P0.05，則滿足球形假設。本例中，χ2=6.270，P=0.043，所以不滿足球形假設。

(2) 當違背了球形假設條件時，需要進行epsilon (ε)校正。如下圖突出顯示，SPSS共用了三種方法進行校正，分別為：Greenhouse-Geisser、Huynh-Feldt 和Lower-bound。

在實際應用中，只用Greenhouse-Geisser和Huynh-Feldt兩種方法，這兩種方法計算的epsilon (ε)的值越低，說明違反球形假設的程度越大，當epsilon (ε)=1時，說明完美的服從了球形假設。

Maxwell & Delaney (2004)建議當epsilon (ε)0.75時，使用Huynh-Feldt方法校正。

(3) 滿足球形假設的結果

在Tests of Within-Subjects Effects表中，該表給出了crp_pre、crp_mid和crp_post在不同時間點的均值是否存在差異；偏η2（Partial Eta Squared）表示控制了其他自變數後，因變數被某一自變數解釋的百分比（單因素方差分析時，即自變數對因變數的解釋程度），代表樣本的效應量。

如果P>0.05，則表示各組間均數差異無統計學意義（本例中，P值顯示為0.000，不代表P值實際為0，而是表示P

以上具體結果見下表：

目前許多雜誌要求除了列出上述結果，還要報告效應量。計算單因素方差分析的效應量有很多方法，比較推薦的是計算偏ω2，公式如下：

公式里k =受試者內因素的個數，F為F統計量，n=研究對象的數量，計算結果如下：

可見偏ω2=0.63

(4) 不滿足球形假設的結果

當不滿足球形假設時，可以採用Greenhouse & Geisser方法進行校正，如下表中突出顯示的內容。可見，自由度（df）由原來符合球形假設時的2變成了1.296，均方（Mean Square）由原來的1.164變成1.797，方差不變。

Greenhouse & Geisser的校正結果顯示不同的時間點的CRP濃度的差異具有統計學意義，F（1.296, 11.663）=26.938，P

以上具體結果見下表：

(5) Bonferroni post hoc test檢驗結果

如果事前沒有對特定組間差異進行假設，或者關心所有組間的兩兩比較，則應該進行所有組間的兩兩比較（post hoc test）。推薦使用Bonferroni post hoc test方法進行組間兩兩比較。該檢驗不僅提供了每兩個組間比較的P值，也給出了均數差值的可信區間。

在前面我們已經知道球形假設不容易滿足，會影響單因素重複測量方差分析的結果。Bonferroni post hoc test 的實質是在進行Bonferroni對顯著性水平校正的基礎上，配對t檢驗的多重比較。結果見Pairwise Comparisons表格：

本例中受試者內因素有3個，因此會有3種不同的組間組合。當我們要比較crp_pre和crp_mid時，可以見如下標黃部分。兩行均為crp_pre和crp_mid比較的結果。

我們回顧一下之前的結果，可見上表中第一行標黃部分Mean Difference（I-J）為crp_pre與crp_mid之差：4.3300 - 3.9400 = 0.390。而第二行的標黃部分則為crp_mid與crp_pre之差： 3.9400 - 4.3300 = -0.390。

另外，見下圖標黃部分，可見標準誤和P值是相同的。

我們再對比一下95%的置信區間，見下表中標黃部分，置信區間的上限和下限符號不同，絕對值相同。

通過對上述數值的理解，我們可以報告以上結果：干預前CRP濃度為4.33 ± 0.64mg/L，干預中三個月時CRP濃度為3.94 ± 0.57mg/L，比干預前顯著降低了0.390（95%置信區間：0.242-0.538）mg/L（P

1) 自定義組間比較的結果

當我們在syntax編輯器輸入了自定義比較的組別後，我們會得到如下Contrast Results (K Matrix) 表格。此表為干預前和干預後CRP濃度比較的結果，與我們當時輸入的比較順序對應一致。

首先，我們來看Contrast Estimate、Lower Bound和Upper Bound這三行的結果。Contrast Estimate一行結果為0.680，這是所比較的2組間的均數差值；而我們看到置信區間是97.5%置信區間而非95%置信區間，這是因為我們進行了Bonferroni校正，我們還是報告成95%的置信區間。

對於P值，我們雙擊P值，出現下圖，由於進行了Bonferroni校正，顯著性水平為0.025，0.000229

另外，還需要報告效應量，在Univariate Test Results中，偏η2=0.795。

2) 複合比較的結果

如下圖標黃部分所示，是干預前和干預中、干預後兩者的平均值進行比較，干預前與干預中、干預後二者均值的差值為0.535（95%置信區間：0.338-0.732）mg/L，P

另外，還需要報告效應量，在Univariate Test Results中，偏η2=0.855。

七、撰寫結論

1. 假如：滿足球形假設，單因素重複測量方差分析顯示組間差異無統計學意義：

採用單因素重複測量方差分析方法，判斷6個月的鍛煉干預對受試者CRP濃度的影響。經箱線圖判斷，數據無異常值；經Shapiro-Wilk檢驗，各組數據服從正態分布（P>0.05）；經Mauchly s球形假設檢驗，因變數的方差協方差矩陣相等，χ2(2)=3.343,P=0.188。

數據以均數±標準差的形式表示。受試者干預前、干預中和干預後的CRP濃度分別為4.33 ± 0.64 mg/L、3.94 ± 0.57 mg/L和3.65 ± 0.43mg/L。干預前、干預中和干預後的CRP濃度差異不具有統計學意義，F(2, 18)=1.256,P=0.309，偏η2=0.02。

2. 假如：不滿足球形假設，單因素重複測量方差分析顯示組間差異無統計學意義：

採用單因素重複測量方差分析方法，判斷6個月的鍛煉干預對受試者CRP濃度的影響。經箱線圖判斷，數據無異常值；經Shapiro-Wilk檢驗，各組數據服從正態分布（P>0.05）；經Mauchly s球形假設檢驗，因變數的方差協方差矩陣不相等，χ2(2)=6.270,P=0.043，通過Greenhouse & Geisser方法校正ε=0.648。

數據以均數±標準差的形式表示。受試者干預前、干預中和干預後的CRP濃度分別為4.33 ± 0.64 mg/L、3.94 ± 0.57 mg/L和3.65 ± 0.43mg/L。干預前、干預中和干預後的CRP濃度差異不具有統計學意義，校正後F(1.296, 11.663)=1.256,P=0.300，偏η2=0.02。

3.不滿足球形假設，單因素重複測量方差分析顯示組間差異有統計學意義，並進行了兩兩比較（本文例子的結果）：

數據以均數±標準差的形式表示。受試者干預前、干預中和干預後的CRP濃度分別為4.33 ± 0.64 mg/L、3.94 ± 0.57 mg/L和3.65 ± 0.43mg/L。干預前、干預中和干預後的CRP濃度差異具有統計學意義，校正後F(1.296, 11.663)=26.938,P

干預中三個月時CRP濃度比干預前顯著降低了0.390（95%置信區間：0.242-0.538）mg/L（P

4.從無效假設和備擇假設的角度出發，當單因素重複測量方差分析顯示組間差異有統計學意義時：各組間均數差異有統計學意義（P

5.從無效假設和備擇假設的角度出發，當單因素重複測量方差分析顯示組間差異無統計學意義時：各組間均數差異無統計學意義（P>0.05）。因此，不能拒絕無效假設，不能接受備擇假設。

八、繪製圖表

最後，我們來學習如何在SPSS軟體中繪製簡單線圖，從而更好地展示展示單因素重複測量方差分析的結果，使其更適合於學術發表。

1. 繪製簡單線圖的操作

(1) 在菜單欄中，點擊Graphs > Chart Builder...：

出現如下Chart Builder對話框：

(2) 在Chart Builder對話框的左下角，Choose from:模塊中選擇「Line」：

(3) 在Chart Builder對話框的中下部，出現2個不同的線圖選項，把左側的第一個（Simple Line）拖進上面的主要圖表預覽窗口，並點擊Element Properties：

(4) 出現下圖，圖表預覽窗口的線圖橫縱軸分別顯示「X-Axis?」和「Y-Axis?」：

(5) 從Variables:模塊中把自變數time拖進「X-Axis?」，把因變數crp拖進「Y-Axis?」：

(6) 在Element Properties對話框中勾選Display error bars，激活–Error Bars Represent– 模塊，勾選Confidence intervals，Level (%):設定為95，當然也可以根據需要勾選Standard error或者Standard deviation，點擊Apply：