方差分析python实现

单因素方差分析

在此之间我们先导入数据，

案例：我们探究施肥与否是否对植物的生长有影响，试验为： - 对照组：清水 - 实验组：某肥料四个浓度梯度，分别是A,B,C,D，施肥一段时间之后测量树高（要控制其他变量保持一致，比如施肥之前的树高要基本保持一致，生长势基本保持一致等等）

在方差分析满足，独立性，正态性，方差齐性，虽然没有满足方差齐性也可以进行非参数的检验，下面我们开始对数据进行处理

import numpy as np
import pandas as pd
import matplotlib.pyplot as pltdf  = {'baseline':list(np.random.normal(10, 5, 10)),'treat1':list(np.random.normal(15, 5, 10)),'treat2':list(np.random.normal(20, 5, 10)),'treat3':list(np.random.normal(30, 5, 10)),'treat4':list(np.random.normal(31, 5, 10))}data = pd.DataFrame(df)

我们看一下随机生成的数据长什么样

大概了解一下描述性统计数据

我们为了方便计算，将所有的数据合成一列，并画一下箱线图看一下

import seaborn as snsdata_melt = data.melt()
data_melt.columns = ['Treat', 'value']sns.boxplot(x = 'Treat', y = 'value', data = data_melt)

由于我们在生成随机数的时候为了让结果明显一点，所以从这个图我们也很容易看到差异性还是很明显的

方差分析

# 方差分析
from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm
from statsmodels.stats.multicomp import pairwise_tukeyhsdmodel = ols('value ~C(Treat)', data = data_melt).fit()
anova_table = anova_lm(model, type = 2)
pd.DataFrame(anova_table)# 进行事后比较分析
print(pairwise_tukeyhsd(data_melt['value'], data_melt['Treat']))

我们进行方差分析之后得到的结果如下：

我们可以看到因子（Treat）的不同水平，对于因变量的影响还是很显著的，p值显然小于0.05.

然后我们进一步进行多重比较，事后分析：

group1以及group2表示的是因子的不同水平，然后分析他们两个组是否有显著性差异，最后面的reject表示是否拒绝原假设，True表示的是拒绝原假设，说明两组均值有显著性差异

多因素方差分析

首先我们将数据导进来，看一下分布

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm
from statsmodels.stats.multicomp import pairwise_tukeyhsddata = pd.read_csv('anovaa.csv')
data = pd.DataFrame(data)
data.head()

便于理解，我们可以将people因子认为是不同的操作员，group因子可以理解成使用不同的方法操作，因变量就是value结果怎么样

按照惯例，我们看一下数据的描述性统计数据

我们看一下箱线图情况

sns.boxplot(x = 'group', y = 'value', data = data)

sns.boxplot(x = 'people', y = 'value', data = data)

model = ols('value ~C(group) + C(people) + C(group):C(people)', data = data).fit()
anova_table = anova_lm(model, type = 2)
pd.DataFrame(anova_table)

这里得到每一个因子不同水平的均值差异与否的检验
还有就是两个因子的交互效应

结果得到：group组是小于0.05的，people和两者的交互效应是大于0.05的，所以这样就明白了他们的交互效应

接下来进行进一步多重检验

print(pairwise_tukeyhsd(data['value'], data['people']))
print(pairwise_tukeyhsd(data['value'], data['group']))

这是people因子不同水平的比较结果

这是group不同水平的比较结果