統計與數據科學方差分析簡介(以疫情為例)

「事實是每個人都相信的簡單陳述。這是無辜的除非被判有罪。假設是一種沒有人願意相信的新建議。這是有罪的，直到發現有效為止。「-愛德華·泰勒

我們正在努力應對一種前所未有的大流行。世界各地的研究人員都在瘋狂地試圖研製出新冠肺炎的疫苗或解藥，而醫生們只是想防止這種流行病席捲整個世界。

最近，我有了一個想法，把我的統計知識應用到最近正在生成的大量COVID數據中。

考慮一種情況，醫生有四種醫療方法來治療病人。一旦我們有測試結果，一種方法是假設治療時間最少的病人是其中最好的。

但是，如果這些病人中的一些已經部分治癒了，或者其他藥物已經在起作用了呢？

為了作出一個有信心和可靠的決定，我們需要有證據來支持我們的做法。這就是變異數概念發揮作用的地方。

在本文中，我將向您介紹ANOVA測試及其用於做出更好決策的不同類型。蛋糕上的糖霜？我將用Python演示每種類型的ANOVA測試，以可視化它們在新冠肺炎數據上的工作方式。我們走吧！

注意：要理解這個主題，您必須了解統計的基本知識。對.的認識T檢驗和假設檢驗將是一個額外的好處。

什麼是ANOVA測試？

方差分析(ANOVA)可以看作是2組以上t檢驗的推廣。採用獨立t檢驗方法比較兩組患者的病情均值.當我們想要比較兩個以上組之間的狀況時，就會使用Anova。

Anova檢驗模型中的平均值是否存在差異(測試是否有總體效果)，但它不能告訴我們差異在哪裡(如果有)。為了找出各組之間的不同之處，我們必須進行臨時測試。

要執行任何測試，我們首先需要定義NULL和備用假設：

• 零假設-各組間沒有顯著性差異。
• 交替假設-各組間存在顯著差異。

基本上，ANOVA是通過比較兩種類型的變異，即樣本均值之間的變化以及每個樣本內部的變化來實現的。下面的公式表示單向Anova測試統計數據.

方差分析公式F統計量(也稱為F比)的結果允許對多組數據進行分析，以確定樣本之間和樣本內部的差異。

單向方差分析的公式可以寫成這樣：

當我們繪製ANOVA表時，所有上述成分都可以在其中看到如下：

通常，如果與F相關的p值小於0.05，則空假設被拒絕，替代假設被支持。如果零假設被拒絕，我們可以得出結論，所有組的均值都是不相等的。

註：如果被測組之間沒有真正的差異，即所謂的零假設，方差分析的F-比率統計量的結果將接近1。

方差分析假設

在進行方差分析之前，我們需要做一些假設：

1. 這些觀測是從由因子水平定義的種群中獨立而隨機地得到的。
2. 每個因素級別的數據都是正態分布的。
3. 案例的獨立性：樣本案例應該是相互獨立的。
4. 方差的同質性：同質性是指各組間的方差近似相等。

方差齊性假設可以用Levene檢驗或Brown-Forsythe檢驗來檢驗。分數分布的正態性可以用直方圖、偏度和峰度值來檢驗，也可以用Shapiro-Wilk或Kolmogorov-Smirnov或Q-Q圖等測試來檢驗。獨立的假設可以從研究的設計中確定。

必須指出的是，方差分析對違反獨立假設的情況並不是很好。這就是說，即使你違反了同質性或正態性的假設，你也可以進行測試，基本上相信調查結果。

但是，如果違反了獨立性假設，則方差分析的結果是無效的。通常，如果您有相同大小的組，則如果違反同質性，則該分析被認為是可靠的。如果你有一個大樣本的話，如果你違反了正常情況，那麼繼續進行方差分析是可以的。

方差分析類型s

1. 單向方差分析：單因素方差分析只有一個自變量。例如，電暈病例的差異可以根據國家進行評估，一個國家可以有2、20種或更多的不同類別進行比較。
2. 雙向方差分析：雙向方差(也稱為階乘方差)是指使用兩個自變量的方差。擴大上面的例子，一個雙向的方差分析可以檢查不同年齡組(自變量1)和性別(自變量2)的日冕病例(因變量)的差異。雙向方差分析可以用來檢驗兩個自變量之間的相互作用.交互作用表明，不同類別的自變量之間的差異並不是一致的。例如，與年輕組相比，老年組的總體冠狀病例可能更高，但與歐洲國家相比，亞洲國家的這一差異可能更大(或更小)。
3. 正向方差研究人員也可以使用兩個以上的自變量，這是一個n向方差(n是自變量的數目)，也就是Manova檢驗。例如，可以同時按國家、性別、年齡組、種族等審查電暈病例的潛在差異。ANOVA會給你一個單一的(單變量)f值，而Manova會給出一個多變量的F值。

使用複製與不複製

您可能經常聽到與複製和不複製有關的方差分析測試。讓我們了解一下這些是什麼：

1. 雙向複製方差分析： 兩組這些團體的成員是做不止一件事例如，假設新冠肺炎還沒有研製出疫苗，醫生們正在嘗試兩種不同的治療方法來治療新冠肺炎感染的兩組病人。
2. 雙向變異數不複製：當你只有一組而你雙重測試同一組例如，假設已經為新冠肺炎研製了一種疫苗，研究人員正在測試一組志願者接種疫苗前後的情況，以確定疫苗是否有效。

後方差分析

當我們進行方差分析時，我們試圖確定各組間是否存在統計學上的顯著差異。如果我們發現存在差異，那麼我們就需要檢查組間的差異在哪裡。

因此，基本上，特設測試告訴研究人員哪些群體是不同的。

此時，您可以運行以下臨時測試：T檢驗檢查各組之間的平均差異。可以進行多個比較測試，以控制I類錯誤率，包括Bonferroni、Scheffe、Dunnet和Tukey測試。

現在，讓我們用一些真實的數據來理解每種類型的ANOVA測試，並使用Python來探索它。

Python中的單向方差分析

我從正在進行的Kaggle競賽中下載了這些數據。請隨意使用它。

在這個測試中，我們將嘗試分析一個區域或狀態的密度與日冕例數之間的關係。因此，我們將根據居住在其中的人口密度繪製每個州的地圖。

讓我們從導入所有必需的庫和數據開始：

import pandas as pdimport numpy as npimport scipy.stats as statsimport osimport randomimport statsmodels.api as smimport statsmodels.stats.multicompfrom statsmodels.formula.api import olsfrom statsmodels.stats.anova import anova_lmimport matplotlib.pyplot as pltfrom scipy import statsimport seaborn as sns

從目錄加載數據：

StatewiseTestingDetails=pd.read_csv(&39;)population_india_census2011=pd.read_csv(&39;)

StatewiseTestingDetails包含關於每個州一天中總陽性和陰性病例的信息。鑑於人口_印度_人口普查2011年載有關於每個邦的密度的信息和關於人口的其他相關信息。

population_india_census2011.head() take glimpse of dataStatewiseTestingDetails[&39;].sort_values().head() 39;State&39;Positive&39;State&39;Positive&39;Positive&39;Median&39;Positive&39;Positive&39;Median&39;State&39;State&Merge StatewiseTestingDetails & population_india_census2011 dataframesdata=pd.merge(StatewiseTestingDetails,population_india_census2011,on=&39;)

現在，我們可以編寫一個函數，根據每個狀態的密度創建一個密度組桶，其中dense 1<dense 2<dense 3<dense 4：

def densityCheck(data):    data[&39;]=0    for index,row in data.iterrows():        status=None        i=row[&39;].split(&39;)[0]        try:            if (&39; in i):                i=int(i.split(&39;)[0]+i.split(&39;)[1])            elif (&39; in i):                i=round(float(i))            else:                i=int(i)        except ValueError as err:            pass        try:            if (0<i<=300):                status=&39;            elif (300<i<=600):                status=&39;            elif (600<i<=900):                status=&39;            else:                status=&39;        except ValueError as err:            pass        data[&39;].iloc[index]=status    return data

現在，用它的密度組映射每個狀態。我們可以導出這些數據，這樣我們可以在以後的雙向方差分析中使用這些數據：

data=densityCheck(data)39;ll export this data so we can use it for Two - way ANOVA test.stateDensity=data[[&39;,&39;]].drop_duplicates().sort_values(by=&39;)

讓我們對可以用於ANOVA測試的數據集進行子集和重新排序：

我們的方差分析假設之一是，樣本應該是隨機選擇的，並且應該接近高斯分布。因此，讓我們從每個因素或級別中選擇10個隨機樣本：

np.random.seed(1234)dataNew=pd.DataFrame({&39;:random.sample(list(data[&39;][data[&39;]==&39;]), 10),&39;:random.sample(list(data[&39;][data[&39;]==&39;]), 10),&39;:random.sample(list(data[&39;][data[&39;]==&39;]), 10),&39;:random.sample(list(data[&39;][data[&39;]==&39;]), 10)})

讓我們繪製一個電暈病例數量的密度分布圖，以檢查它們在不同密度組中的分布：

我們可以清楚地看到，數據不遵循高斯分布。

有不同的數據轉換方法可以使數據接近高斯分布。我們將繼續進行Box Cox變換：

dataNew[&39;],fitted_lambda = stats.boxcox(dataNew[&39;])dataNew[&39;],fitted_lambda = stats.boxcox(dataNew[&39;])dataNew[&39;],fitted_lambda = stats.boxcox(dataNew[&39;])dataNew[&39;],fitted_lambda = stats.boxcox(dataNew[&39;])

現在，讓我們再一次繪製它們的分布圖，以檢查：

方法1：使用狀態模型模塊進行單向方差分析

Python中有幾種執行ANOVA測試的方法。一個是stats.f_oneway()方法：

F, p = stats.f_oneway(dataNew[&39;],dataNew[&39;],dataNew[&39;],dataNew[&39;])39;F-Statistic=%.3f, p=%.3f&39;Count ~ C(density_Group)& Seeing if the overall model is significantprint(f&34;)39;Count&39;density_Group& Normality Assumption checkw, pvalue = stats.shapiro(model.resid)print(w, pvalue)

從上面的代碼片段中，我們可以看到，對於所有的密度組，p值都>0.05。因此，我們可以得出這樣的結論：它們遵循高斯分布。

方法2：Q-Q圖測試：

我們可以使用正常的Q-Q圖來檢驗這個假設：

res = model.residfig = sm.qqplot(res, line=&39;)plt.show()

從上面的數字中，我們可以看到所有的數據點都靠近45度線，因此我們可以得出它服從正態分布的結論。

方差假設檢驗的同質性

對於分類變量的每個級別，應檢查方差假設的同質性。我們可以用Levene試驗來檢驗各組間的等方差。

w, pvalue = stats.bartlett(newDf[&39;][newDf[&39;]==&39;], newDf[&39;][newDf[&39;]==&39;]                           , newDf[&39;][newDf[&39;]==&39;], newDf[&39;][newDf[&39;]==&39;])print(w, pvalue)39;Dense1&39;Dense2&39;Dense3&39;Dense4&39;./individualDetails.csv&39;./stateDensity.csv&39;./individualDetails.csv&39;./stateDensity.csv&39;Percentage of missing values in age & gender columns respectively :&39;age&39;%&39;gender&39;%&39;age&39;age&39;age&39;State&39;age&39;age&39;age&Find the most frequent infected gender by COVID-19 for each stategenderModePerState=individualDetails.groupby([&39;])[&39;].agg(pd.Series.mode).to_frame().reset_index()39;State&39;Arunachal Pradesh&Impute missing values in age & gender columns nowfor index,row in individualDetails.iterrows():    if row[&39;]==&39;:            individualDetails.drop([index],inplace=True)        continue    if pd.isnull(row[&39;]):        individualDetails[&39;][index]=list(ageMedianPerState[&39;][ageMedianPerState[&39;]==&39;].values)[0]    if pd.isnull(row[&39;]):        if len(genderModePerState[&39;][genderModePerState[&39;]==row[&39;]].values)>0:            individualDetails[&39;][index]=(genderModePerState[&39;][genderModePerState[&39;]==row[&39;]].values[0])

現在，讓我們合併個性化Details&StateDensityDataaframe，為我們創建一個整體數據集：

data = pd.merge(individualDetails,stateDensity,on=&39;,how=&39;).reset_index(drop=True)

現在我們可以創建年齡組桶：

data.dropna(subset=[&39;],inplace=True)data.reset_index(drop=True,inplace=True)

將數據合併以獲得與年齡組和各自狀態密度組對應的每個人的數據集：

patient_Count=data.groupby([&39;,&39;])[[&39;]].count().\                        rename(columns={&39;:&39;}).reset_index()data=pd.merge(data,patient_Count,on=[&39;,&39;],how=&39;)newData=data.groupby([&39;,&39;])[&39;].apply(list).reset_index()newData.head()np.random.seed(1234)AnovaData=pd.DataFrame(columns=[&39;,&39;,&39;])for index,row in newData.iterrows():    count=17    tempDf=pd.DataFrame(index=range(0,count),columns=[&39;,&39;,&39;])       tempDf[&39;]=newData[&39;][index]    tempDf[&39;]=newData[&39;][index]    tempDf[&39;]=random.sample(list(newData[&39;][index]),count)    AnovaData=pd.concat([AnovaData,tempDf],axis=0)

檢查數據中計數列的分布情況，並使用方框繪圖方法檢查數據中是否存在異常值：

plt.hist(AnovaData[&39;])plt.show() sns.kdeplot(AnovaData[&39;],cumulative=False,bw=2)

從上面的代碼片段中，我們可以看到沒有感染嬰兒的記錄。接下來，檢查數據中缺少的值：

individualDetails.isna().sum()print(&39;, \      (individualDetails[&39;].isna().sum()/individualDetails.shape[0])*100,&39;,\      (individualDetails[&39;].isna().sum()/individualDetails.shape[0])*100,&39;)

我們發現在年齡和性別欄中分別有91%和80%的條目缺失。因此，我們需要設計一種方法來推演它們。

在這裡，我將根據每個州的男性和女性的模式計算出每個州的中位數和性別。因此，我將計算中值和模式：

ageMedianPerState=individualDetails[~individualDetails[&39;].isna()]ageMedianPerState[&39;]=ageMedianPerState[&39;].astype(str).astype(int)ageMedianPerState=ageMedianPerState.groupby(&39;)[[&39;]].median().reset_index()ageMedianPerState[&39;]=ageMedianPerState[&39;].apply(lambda x:math.ceil(x))39;State&39;gender&Drop Arunachal Pradesh since this has got no informatio about gender overallgenderModePerState=genderModePerState[genderModePerState[&39;]!=&39;]

39;State&39;Arunachal Pradesh&39;age&39;age&39;age&39;State&39;West Bengal&39;gender&39;gender&39;State&39;State&39;gender&39;gender&39;State&39;State&39;State&39;left&39;density_Group&39;diagnosed_date&39;density_Group&39;diagnosed_date&39;diagnosed_date&39;Count&39;diagnosed_date&39;density_Group&39;inner&39;density_Group&39;age_Group&39;Count&39;density_Group&39;age_Group&39;Count&39;density_Group&39;age_Group&39;Count&39;age_Group&39;age_Group&39;density_Group&39;density_Group&39;Count&39;Count&39;Count&39;Count&39;age_Group&39;Count&34;colorblind&39;Count&39;Count& Data TransformationAnovaData[&39;],fitted_lambda = stats.boxcox(AnovaData[&39;])import matplotlib.pyplot as pltsns.kdeplot(AnovaData[&39;],cumulative=False,bw=2)

現在讓我們使用OLS模型來驗證我們的假設：

39;newCount ~ C(age_Group)+ C(density_Group)&34;Overall model F({model2.df_model: },{model2.df_resid: }) = {model2.fvalue: }, p = {model2.f_pvalue: }& Create the ANOVA tableres2 = sm.stats.anova_lm(model2, typ=2)res239;s&Approach 2 - Check the interaction of groupsformula = &39;model = ols(formula, AnovaData).fit()model.summary()

aov_table = anova_lm(model, typ=2)print(aov_table.round(4))

對結果的解釋

經方差分析得到的日冕病例數、年齡組、密度組及相互作用的P值均有統計學意義(P<0.05)。我們的結論是，密度群的類型對電暈病例的結局有很大的影響。

年齡組對日冕病例結局有顯著影響，年齡組與密度組的交互作用對日冕病例結局也有顯著影響。

後特設測試

最後，讓我們確定哪些組在統計上是不同的。我們將使用Tuckey HSD方法：

mc = statsmodels.stats.multicomp.MultiComparison(AnovaData[&39;],AnovaData[&39;])mc_results = mc.tukeyhsd()print(mc_results)

mc = statsmodels.stats.multicomp.MultiComparison(AnovaData[&39;],AnovaData[&39;])mc_results = mc.tukeyhsd()print(mc_results)

從上述的Tuckey HSD測試結果可以清楚地看出，密度組的1-組3、1-組-4組與年齡組中的青年和青年-老年組之間存在顯著差異。

因此，TUKEY HSD的上述結果表明，除上述各組外，所有其他對冠病例數的兩兩比較都否定了零假設，表明沒有統計學上的顯著性差異。

尾注

在這些大流行時期，我試著用相關的案例研究來解釋方差分析。這是一個有趣的經驗，把這一切為我們的社區！

請隨時在下面評論你的問題