当前位置: 首页 > 科技观察

使用Python的两种方法进行方差分析

时间:2023-03-14 01:08:35 科技观察

在进行数据分析的时候,我们经常会遇到这样的情况,我们要分析某个变量的影响因素,而影响一件事的因素往往有很多。例如在化工生产中,有温度、压力、用量、反应时间等因素。每个因素的变化都可能影响产品的数量和质量。我们常常想找出那些对产品质量有重大影响的因素。方差分析是分析测试结果和识别各种相关因素对测试结果影响的有效方法。本文主要介绍如何使用python中的两种方法进行方差分析。首先,我们先介绍一下方差分析。方差分析(ANOVA),又称“方差分析”或“F检验”,由RonaldAylmerFisher发明,用于两个或多个样本。差异显着性检验的原理是,不同处理组均值差异的基本来源有两个:(1)实验条件,即不同处理引起的差异,称为组间差异。表示为各组均值与变量总均值的偏差平方和之和,记为SSa,组间自由度dfa。(2)随机误差,如因测量误差或个体间的差异而引起的差异,称为组内差异,用各组内变量均值与偏差的平方和表示群中的变量值,记为SSe,群内自由度dfe。总偏差平方和SSt=SSA+SSe。组内SSe和组间SSA除以各自的自由度(组内dfe=n-m,组间dfa=m-1,其中n为样本总数,m为样本数groups),求均方MSe和MSa,a的情况是处理没有效果,即每组样本来自同一个种群,MSa/MSe≈1。另一种情况是处理确实有影响,组间均方是误差和不同处理的结果,即每个样本来自不同的人群。然后,MSa>>MSe(远大于)。MSa/MSe比率构成F分布。将F值与其临界值进行比较,以推断每个样本是否来自同一总体。然后,让我们再次说明数据集。数据集很简单,只有5组值,每组值有4个,一共20个数。它们分别命名为group1、group2、group3、group4和group5。数值都是随意设置的,没有什么要求。在这里,你也可以根据自己的意愿设置数据。在这里,作者特意将数据量设置的比较小,方便观察数据之间的差异。我们的重点是方差分析方法,这里主要讲单因素方差分析方法。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.,25.0,24.2]假设u1,u2,u3,u4,u5为这五个样本所属总体的均值,我们使用单因素方差分析来检验以下假设。H0:u1=u2=u3=u4=u5H1:u1,u2,u3,u4,u5不都相等为了更直观的理解这5组数据,我们先手动计算这5组数据的相关参数。这五组数据的总体情况如图1所示。图1所用数据的基本情况图1中每一列数据都是一个水平,是一个统计项,水平之和就是每组4个值的总和,每组数据的平均值为a1=28.6,a2=31.375,a3=7.825,a4=19.075,a5=27.8,所有20个数据的平均值为A=(a1+a2+a3+a4+a5)/5=114.675/5=22.935。因此,总偏差的平方和为ST=1616.65,即20个数据中各数据与A之差的平方和,误差的平方和为SE=135.82,即是每组数据中每个数据的值与这组数据的均值差的平方和,效果的平方和为SA=1480.83,也就是差的平方和每组数据的均值与A之间,也等于ST减去SE的差值。由此我们可以得到本例的方差分析表,如图2所示。图2.方差分析表图2中的因素是各组数据之间的差异,可以是随机的也可以是人为的,误差是每组数据之间的差异。我们可以看到,本例得到的F值为40.8848,远大于查表得到的F值F0.05(4,15),其值为3.06。至于F0.05(4,15)的值,我们也可以用python获取,后面会说明。以上就是本例的人工计算过程。让我们使用python来计算这个例子。方法一:scipy方法一使用的库是scipy,是python中科学计算最常用的库。代码如下。记得输入前5组数据。fromscipyimportstatsF,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))ifF>=F_test:print('拒绝原假设,u1,u2,u3,u4,u5为notallequal')else:print('Acceptthenullhypothesis,u1=u2=u3=u4=u5')结果如图3所示。图3.方法1的计算结果.scipy的单向方差分析比较简单。只需要调用stats模块的f_oneway方法,在f_oneway中输入每组数据,然后自动返回两个值F和p,第一个值F表示我们计算的F值,与中的F值相同图2,而第二个值p就是这个F值对应的概率,它是从假设检验问题中检验统计量的样本观测值中得到的可以拒绝原假设的最小显着性水平。这里我们不仅可以通过F值来判断,还可以通过p值来判断。因为F大于F_test,落入拒绝域,所以拒绝原??假设,p值也远小于α分位数(这里是0.05),所以也拒绝原假设。这里的F_test是图2中的F0.05(4,15),计算方法是使用stats.f.ppf((1-0.05),4,15),其中ppf表示Percentpointfunction,即Percentagepointfunction,它是Cumulativedistributionfunction(累积分布函数)的逆运算,这里要注意ppf的第一个参数要输入1-0.05,0.05是我们设置的显着性水平α,它的值为通常取0.05,第二个和第三个参数是两个自由度,这两个自由度分别是4和15,计算方法在前面的原理部分已经讲过了。方法二:statsmodels方法二使用了python的另一个统计库statsmodels,其代码如下。importstatsmodels.apiassmimportpandasaspdfromstatsmodels.formula.apiimportolsnum=sorted(['g1','g2','g3','g4','g5']*4)data=group1+group2+group3+group4+group5df=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)结果如图4所示。图4.方法2的计算结果从图4可以看出,得到的结果与前面人工计算和scipy的结果是一样的(可以忽略一些小数精度问题),并且图中sum_sq列表示平方和,df列表示自由度,这里也给出了p值作为PR(>F)列,比scipy信息量更多。从代码上看,statsmodels也很简单,只是比scipy稍微复杂一点,但是提供的信息更多。这里有几点需要注意。一种是我们生成了一个名为num的变量和一个名为data的变量,都是list类型,用这两个生成了一个pandas.DataFrame变量,名为df。这样做的原因是statsmodels中常见的使用DataFrame数据格式,如果使用list类型会比较麻烦。而data就是把之前group1到group5的数据放在一个list中,num就是存放每条数据对应的数据分组信息,g1表示这个值属于group1,g2表示对应group2,以此类推.这里还有一点要注意,num中的数据格式最好是字符格式,比如'a1'和'num3',不要是数字格式,比如1,3,6.9,因为数字格式的数据是很可能会涉及到计算,最后的结果可能是错误的。第二点是mod=ols('data~num',data=df).fit()中的公式data~num。很多人对这一点感到困惑。这个公式的使用来自python的另一个库patsy,主要用来描述统计模型(尤其是线性模型),符号~之前的部分代表y轴数据,后面的部分代表x轴数据.根据两者生成线性模型,ols中的第二个参数data是要输入的数据源,一般是DataFrame格式。上式中符号~前后的名称必须是data中的列名。这个方法确实有点奇怪,部分原因是patsy借鉴了R语言的一些用法。第三点,在ano_table=sm.stats.anova_lm(mod,typ=2)中,typ=2表示DataFrame,typ有3个值,分别是1、2、3,2代表DataFrame格式。总结比较一下scipy和statsmodels这两种方法,可以说各有千秋。scipy是一个通用库,包含了科学计算的各个模块,统计分析只是其中的一部分,而statsmodels是专门用于统计分析的库。两者在功能上有一些区别,statsmodels在统计分析方面更专业。一些。scipy的语法比较符合python常用的语法,statsmodels的语法有些接近R语言,初学者可能比较陌生。所以大家可以根据自己的需要选择合适的方法。