Python中怎么實現(xiàn)方差分析

這篇文章給大家介紹Python中怎么實現(xiàn)方差分析，內(nèi)容非常詳細(xì)，感興趣的小伙伴們可以參考借鑒，希望對大家能有所幫助。

首先，還是先簡介一下方差分析。

方差分析（Analysis of Variance，ANOVA）又稱“變異數(shù)分析”或“F檢驗”，是由羅納德·費舍爾（Ronald Aylmer Fisher）發(fā)明的，用于兩個及兩個以上樣本均數(shù)差別的顯著性檢驗，其原理是認(rèn)為不同處理組的均數(shù)間的差別基本來源有兩個：

(1) 實驗條件，即不同的處理造成的差異，稱為組間差異。用變量在各組的均值與總均值之偏差平方和的總和表示，記作SSa，組間自由度dfa。

(2) 隨機誤差，如測量誤差造成的差異或個體間的差異，稱為組內(nèi)差異，用變量在各組的均值與該組內(nèi)變量值之偏差平方和的總和表示， 記作SSe，組內(nèi)自由度dfe。

總偏差平方和 SSt = SSa + SSe。

組內(nèi)SSe、組間SSa除以各自的自由度（組內(nèi)dfe =n-m，組間dfa=m-1，其中n為樣本總數(shù)，m為組數(shù)），得到其均方MSe和MSa，一種情況是處理沒有作用，即各組樣本均來自同一總體，MSa/MSe≈1。另一種情況是處理確實有作用，組間均方是由于誤差與不同處理共同導(dǎo)致的結(jié)果，即各樣本來自不同總體。那么，MSa>>MSe（遠(yuǎn)大于）。

MSa/MSe比值構(gòu)成F分布。用F值與其臨界值比較，推斷各樣本是否來自相同的總體。

然后，我們再說明一下數(shù)據(jù)集。

數(shù)據(jù)集非常簡單，只有5組數(shù)值，每組數(shù)值有4個，共20個數(shù)字。分別命名為group1、group2、group3、group4和group5，數(shù)值都是隨意設(shè)置的，沒有什么要求，這里大家也可以根據(jù)自己的意愿設(shè)置數(shù)據(jù)。在這里，筆者專門將數(shù)據(jù)量設(shè)置得比較小，這樣方便觀察數(shù)據(jù)的之間的差異，我們的重點是方差分析的方法，而這里我們主要講的是單因素方差分析法。

group1 = [29.6, 24.3, 28.5, 32.0]

group2 = [27.3, 32.6, 30.8, 34.8]

group3 = [5.8, 6.2,11.0, 8.3]

group4 = [21.6, 17.4, 18.3, 19.0]

group5 = [29.2, 32.8, 25.0, 24.2]

設(shè)u1、u2、u3、u4和u5分別是這5個樣本所屬總體的均值，我們用單因素方差分析來檢驗下面的假設(shè)。

H0：u1=u2=u3=u4=u5

H1：u1、u2、u3、u4和u5不全相等

為了能更直觀了解這5組數(shù)據(jù)，我們首先手工計算一下這些數(shù)據(jù)的相關(guān)參數(shù)。這5組數(shù)據(jù)的總體情況如圖1所示。

圖1. 所用數(shù)據(jù)的基本情況

在圖1中，每列數(shù)據(jù)就是一個水平，這是一個統(tǒng)計學(xué)用語，水平和就是每組4個數(shù)值的總和，每組數(shù)據(jù)平均值分別是a1=28.6，a2=31.375，a3=7.825，a4=19.075，a5=27.8，全部20個數(shù)據(jù)的平均值為A=(a1+a2+a3+a4+a5)/5=114.675/5=22.935。所以總偏差平方和為ST=1616.65，此值為20個數(shù)據(jù)中每個數(shù)據(jù)與A的差的平方的總和，誤差平方和為SE=135.82，此值為每組數(shù)據(jù)中每個數(shù)據(jù)與這組數(shù)據(jù)的平均值的差的平方之和，效應(yīng)平方和為SA=1480.83，此值為每組數(shù)據(jù)的平均值與A的差的平方之和，也等于ST減去SE的差。由此我們可以得出本例的方差分析表，如圖2所示。

圖2. 方差分析表

圖2中的因素就是各組數(shù)據(jù)間的差異，這個可以是隨機的，也可以是人為的，而誤差就是每組數(shù)據(jù)的之間差異。我們可以看到本例中得到的F值為40.8848，遠(yuǎn)大于查表得到的F值F0.05(4,15)，其值為3.06，至于F0.05(4,15)的值我們同樣可以用python得出，后面會有講解。

以上就是這個例子的手工計算過程，下面我們用python來計算一下該例。

方法1：scipy

方法1用的庫是scipy，這是python中科學(xué)計算最常用的庫，其代碼如下，記得輸入前面的5組數(shù)據(jù)。

from scipy import stats  F, p = stats.f_oneway(group1, group2, group3, group4, group5)  F_test = stats.f.ppf((1-0.05), 4, 15) print('F值是%.2f，p值是%.9f' % (F,p))  print('F_test的值是%.2f' % (F_test))  if F>=F_test:  print('拒絕原假設(shè)，u1、u2、u3、u4、u5不全相等')  else:  print('接受原假設(shè)，u1=u2=u3=u4=u5')

結(jié)果如圖3所示。

圖3. 方法1的計算結(jié)果

scipy的單因素方差分析比較簡單，只要調(diào)用stats模塊的f_oneway方法即可，在f_oneway中輸入各組數(shù)據(jù)，然后會自動返回兩個數(shù)值F與p，第一個數(shù)值F就表示我們算出的F值，和圖2中的F值一樣，而第二個值p就是這個F值所對應(yīng)的概率，也就是假設(shè)檢驗問題中，由檢驗統(tǒng)計量的樣本觀察值得出的原假設(shè)可被拒絕的最小顯著性水平。在這里我們既可以通過F值來判斷，也可以通過p值來判斷，因為F大于F_test，落入了拒絕域，所以拒絕原假設(shè)，而p值也遠(yuǎn)小于α分位數(shù)（這里為0.05），所以也拒絕原假設(shè)。而這里的F_test就是圖2中的F0.05(4,15)，計算方法就是用stats.f.ppf((1-0.05), 4, 15)，這里ppf的意思是Percent point function，也就是百分點函數(shù)，它是Cumulative distribution function（累積分布函數(shù)）的逆運算，這里需要注意的是ppf的第一個參數(shù)要輸入1-0.05，0.05也就是我們設(shè)定的顯著性水平α，其值通常取0.05，而第二個和第三個參數(shù)是兩個自由度，這兩個自由度分別是4和15，其求法在前面原理部分已經(jīng)講過。

方法2：statsmodels

方法2用的是python的另一個統(tǒng)計學(xué)庫statsmodels，其代碼如下。

import statsmodels.api as sm  import pandas as pd  from statsmodels.formula.api import ols  num = sorted(['g1', 'g2', 'g3','g4', 'g5']*4)  data = group1 + group2 + group3 + group4 + group5  df = pd.DataFrame({'num':num, 'data': data})  mod = ols('data ~ num', data=df).fit()          ano_table = sm.stats.anova_lm(mod, typ=2)  print(ano_table)

結(jié)果如圖4所示。

圖4. 方法2的計算結(jié)果

從圖4中我們可以看到，得出的結(jié)果和前面手算以及scipy的結(jié)果一樣（部分小數(shù)精度問題可以忽略不計），圖中sum_sq列就表示平方和，df列就代表了自由度，這里還給出了p值就是PR(>F)列，信息比scipy要豐富一些。

從代碼上來看，statsmodels也同樣很簡單，只比scipy稍微復(fù)雜了一點，但卻提供了更多的信息。這里有幾點要注意的。一是我們生成了一個名為num的變量和一個名為data的變量，這兩個都是list類型，又用二者生成了名為df的pandas.DataFrame變量，這樣做的原因是statsmodels中普遍使用DataFrame數(shù)據(jù)格式，如果使用list類型會更麻煩一些。而data是把前面group1到group5中的數(shù)據(jù)放在了一個list中，num則是存放每個數(shù)據(jù)所對應(yīng)的數(shù)據(jù)組信息，g1就代表這個數(shù)值屬于group1，g2則是對應(yīng)group2，以此類推。這里還有一點要注意，就是num中數(shù)據(jù)格式最好是字符格式的，比如’a1’、‘num3’這樣的，不要是數(shù)字格式的，比如1、3、6.9這樣的，因為數(shù)字格式的數(shù)據(jù)很有可能會參與計算，最終的結(jié)果可能會出錯。第二點是mod = ols('data ~ num', data=df).fit()中的公式data ~ num，很多人對這一點很困惑，這種公式的使用方法來自于python的另一個庫patsy，其主要用于描述統(tǒng)計模型（尤其是線性模型），符號~前面的部分代表了y軸數(shù)據(jù)，后面的部分代表了x軸數(shù)據(jù)，根據(jù)這二者生成一個線性模型，ols中第二個參數(shù)data則是要輸入的數(shù)據(jù)源，一般是DataFrame格式，前面公式中符號~前后的名稱都要是data中的列名，這種方法確實有些奇怪，部分原因是patsy借鑒了R語言的一些用法。第三點是ano_table = sm.stats.anova_lm(mod, typ=2)中，typ=2的意思是DataFrame，typ共有3個值，分別是1、2和3，其中2代表了DataFrame格式。

關(guān)于Python中怎么實現(xiàn)方差分析就分享到這里了，希望以上內(nèi)容可以對大家有一定的幫助，可以學(xué)到更多知識。如果覺得文章不錯，可以把它分享出去讓更多的人看到。