一、教学内容
 
多因素方差分析，用于研究一个因变量是否受到多个自变量(也称为因素)的影响，它检验多个因素取值水平的不同组合之间，因变量的均值之间是否存在显著的差异。多因素方差分析既可以分析单个因素的作用(主效应)，也可以分析因素之间的交互作用(交互效应)，还可以进行协方差分析，以及各个因素变量与协变量的交互作用。
 根据观测变量(即因变量)的数目，可以把多因素方差分析分为：单变量多因素方差分析(也叫一元多因素方差分析)与多变量多因素方差分析(即多元多因素方差分析)。本文将重点讲述一元多因素方差分析，下篇文章将详细讲述多元多因素方差分析。
 一元多因素方差分析：只有一个因变量，考察多个自变量对该因变量的影响。例如：分析不同品种、不同施肥量对农作物产量的影响时，可将农作物产量作为观测变量，品种和施肥量作为控制变量。利用多因素方差分析方法，研究不同品种、不同施肥量是如何影响农作物产量的，并进一步研究哪种品种与哪种水平的施肥量是提高农作物产量的最优组合。
 01
 分析原理
 通过计算F统计量，进行F检验。F统计量是平均组间平方和与平均组内平方和的比。
 
这里，把总的影响平方和记为SST，它分为两个部分，一部分是由控制变量引起的离差，记为SSA(组间离差平方和)，另一部分是由随机变量引起的SSE(组内离差平方和)。即SST=SSA+SSE。
 组间离差平方和SSA是各水平均值和总体均值离差的平方和，反映了控制变量的影响。组内离差平方和是每个数据与本水平组平均值离差的平方和，反映了数据抽样误差的大小程度。
 通过F值看出，如果控制变量的不同水平对观测变量有显著影响，那观测变量的组间离差平方和就大，F值也大；相反，如果控制变量的不同水平没有对观测变量造成显著影响，那组内离差平方和就比较大，F值就比较小。
 同时，SPSS还会依据F分布表给出相应的相伴概率值sig。如果sig小于显著性水平(一般显著性水平设为0.05、0.01、或者0.001)，则认为控制变量不同水平下各总体均值有显著差异，反之，则不然。一般地，F值越大，则sig值越小。
 02
 SPSS分析案例
 现在有一个公司员工的工资表，想看一下员工性别“gender”与接受教育年限“edu”这两个控制变量对员工“当前工资”的影响。采用多因素方差分析法，则要分别考虑“gender”、“edu”对“当前工资”的影响，称为主效应，还要考虑“gender*edu”对“当前工资”的影响，称为交互效应。
 (1)分析步骤：将数据导入SPSS后，选择：分析—— 一般线性模型——单变量
 
(2)将“当前工资”选入因变量(也就是观测变量)，将性别“gender”与受教育年限“edu”选入固定因子(也就是控制变量)。
 
(3)选择“单变量”的“模型”，打开对话框后选择“全因子”，表示方差分析的模型包括所有因素的主效应，也包括因素之间的交互效应。然后“继续”。
 
(4)打开“单变量”的“绘制”对话框，选择“gender”为横轴变量，选择“edu”为分线变量，单击“添加”，即显示这两个因素变量的交互作用，即 “gender*edu”这个交互作用变量。
 由于此例中“gender”只有两个水平，即男、女；而“edu”有多种水平。因此，如果主效应显著，则表明因素两种或多种水平之间存在显著性差异。事后可以继续对同一因素多个水平之间的均值差异进行比较，该过程称为多重比较。
 但实际上如果主效应和交互效应都达到显著，我们更关心在多因素交互作用下，因变量有什么影响。
 因此，如果交互效应显著的话，通常需要进行简单效应检验。所谓简单效应检验，是指一个因素的水平在另一个因素的某个因素的某个水平上的变异。例如我们本例中的，如果gender与edu之间存在显著的交互作用，我们可以检验当gender为“女”时，edu的各个水平之间的差异，称为edu在“女”性水平上的简单效应；以及在“男”性水平上edu各水平之间的差异，称为edu在“男”性水平上的简单效应。
 简单效应检验，实际上是把其中一个自变量固定在某一个特定的水平上，考察另一个自变量对因变量的影响。简单效应检验在SPSS里是用一个“MANOVA”命令来实现的。
 同理，当我们检验三个自变量时，若这些自变量之间的交互作用显著，需要进行简单简单效应检验，即一个因素的水平在另外两个因素的水平结合上的效应。
 也就是把两个因素固定在各自的某一个水平上，考察第三个因素对因变量的影响。也是用“MANOVA”命令来实现的。我们观察简单效应显著与否，是通过F值与sig值来看的，一般用sig值与我们设定的一个数值(0.05、0.01、或者0.001)来比较，若sig值大于该数值，说明简单效应不显著；反之，若sig值小于该数值，说明简单效应显著。
 
(5)打开“选项”对话框，将左边三个控制变量均移入右边，“显示均值”，同时选中“描述统计”，选中“比较主效应”。
 
(6)点击“确定”以后，就会在SPSS查看器里显示出结果。其中，最上面的那部分代码是我们所做的操作在SPSS里具体实现的步骤的代码。下面的表格是我们想要的结果，从表格里得出结论。
 
(7)从下面的“主体间效应的检验”表格里，我们比较性别gender、受教育程度edu、及genderedu交互作用的F值及sig值，看到edu的F值最大，sig值最小，且sig<0.05。而gender与genderedu的sig值都大于0.05，得出结论：“gender”的主效应未达到显著，而“edu”的主效应达到显著，gender与edu的交互效应未达到显著(当交互效应达到显著时，进而可以进行简单效应检验结果)，就不需要进行简单效应检验。则该公司员工“受教育程度”对员工“当前工资”的影响显著，而“性别”对“当前工资”的影响不明显。
 
(8)下图为均值分布图，即为两因素edu与gender作用下，因变量员工工资的均值分布情况。通常，若交互效应不显著时，图中的因素分布线均为平行线；若交互效应显著，图中的因素分线不平行。
 此图中，将性别“gender”作为横轴变量，观察接受教育年限“edu”对因变量“当前工资”的影响。
 
图中得出结论：当受教育年限为20年,一般为研究生水平的时候，男女工资差别不大；受教育年限为14年，一般为专科生水平，男女工资差别不明显。但当受教育年限为8年、10年、12年、17年的时候，男女工资差别较大，尤其为8年、17年的时候，男女工资差别尤其明显。
 
二、备注
 
相关资料已上传我的资源，下载链接https://blog.csdn.net/TIQCmatlab?spm=1011.2124.3001.5343

1. 前言
 
背景：表格为随机挑选的不同性别与受教育程度的对象的幸福指数数据，
 
目的：现要求分析幸福指数是否受不同的性别和受教育程度影响。
 
分析方法：两个自变量是分类变量，因变量是连续变量，所以选择多因素方差分析。
 
方差分析需要满足的条件：
 
1.各样本须是相互独立的随机样本；
2.各样本来自正态分布总体；
3.各样本方差齐性。
 
显著性水平：选取为0.05
 
工具：Jupyter Notebook（Python 3.8）
 
2. Python
 
数据查看
 
data = [['Male', '高中及以下', 63.0],
       ['Male', '高中及以下', 64.0],
       ['Male', '高中及以下', 60.0],
       ['Male', '高中及以下', 63.0],
       ['Male', '高中及以下', 66.0],
       ['Male', '高中及以下', 65.0],
       ['Male', '高中及以下', 61.0],
       ['Male', '高中及以下', 62.0],
       ['Male', '高中及以下', 68.0],
       ['Male', '大学本科', 66.5],
       ['Male', '大学本科', 67.0],
       ['Male', '大学本科', 69.5],
       ['Male', '大学本科', 70.0],
       ['Male', '大学本科', 69.0],
       ['Male', '大学本科', 71.5],
       ['Male', '大学本科', 67.0],
       ['Male', '大学本科', 68.5],
       ['Male', '大学本科', 63.5],
       ['Male', '研究生及以上', 88.0],
       ['Male', '研究生及以上', 89.0],
       ['Male', '研究生及以上', 86.0],
       ['Male', '研究生及以上', 92.0],
       ['Male', '研究生及以上', 90.0],
       ['Male', '研究生及以上', 84.0],
       ['Male', '研究生及以上', 91.0],
       ['Male', '研究生及以上', 87.0],
       ['Male', '研究生及以上', 88.0],
       ['Male', '研究生及以上', 85.0],
       ['Female', '高中及以下', 65.0],
       ['Female', '高中及以下', 66.0],
       ['Female', '高中及以下', 61.0],
       ['Female', '高中及以下', 64.0],
       ['Female', '高中及以下', 69.0],
       ['Female', '高中及以下', 70.0],
       ['Female', '高中及以下', 67.0],
       ['Female', '高中及以下', 63.0],
       ['Female', '高中及以下', 63.0],
       ['Female', '高中及以下', 66.0],
       ['Female', '大学本科', 70.0],
       ['Female', '大学本科', 71.0],
       ['Female', '大学本科', 66.0],
       ['Female', '大学本科', 69.0],
       ['Female', '大学本科', 74.0],
       ['Female', '大学本科', 73.0],
       ['Female', '大学本科', 72.0],
       ['Female', '大学本科', 68.0],
       ['Female', '大学本科', 65.0],
       ['Female', '大学本科', 64.0],
       ['Female', '研究生及以上', 82.0],
       ['Female', '研究生及以上', 83.0],
       ['Female', '研究生及以上', 88.0],
       ['Female', '研究生及以上', 91.0],
       ['Female', '研究生及以上', 90.0],
       ['Female', '研究生及以上', 86.0],
       ['Female', '研究生及以上', 84.0],
       ['Female', '研究生及以上', 80.0],
       ['Female', '研究生及以上', 85.0],
       ['Female', '研究生及以上', 76.0]]

df = pd.DataFrame(data, columns = ['gender', 'education', 'Index'])
df.head()
 
 
转换
 
df1 = pd.DataFrame()
data_list = []
for i in df.gender.unique():
    for j in df.education.unique():
        data = df[(df.gender == i)&(df.education == j)]['Index'].values
        data_list.append(data)
        df1 = df1.append(pd.DataFrame(data, columns = pd.MultiIndex.from_arrays([[i],[j]])).T)
df1 = df1.T
df1
 
 
转换成更直观的形式，方便后续的操作
 
各组数量统计
 
# 查看各组数量分布
df1.count().to_frame()
 
 
各组数量存在差异
 
箱线图查看异常值
 
plt.figure(figsize = (12,8))
sns.boxplot(x = 'gender', y = 'Index', data = df, hue = 'education')
plt.show()
 
 
各组数据均没有异常值
 
正态性检验
 
# Shapiro-Wilk 检验
sw_test_res = pd.DataFrame()
for i in df1.columns:
    statistic,pvalue = stats.shapiro(df1[i].dropna())
    sw_test_res[i]  = [statistic, pvalue]
sw_test_res.index = ['statistic', 'p_value']
sw_test_res.T.round(3)
 
 
以上各组p值均大于0.05，满足正态分布
 
方差齐性检验
 
# levene test
print('基于中位数的levene test P值：', stats.levene(*data_list, center='median').pvalue)
 
P值为0.286，大于0.05，即任一分类都具有等方差性
 
判断交互作用
 
计算平均值
 
df_mean = df1.mean().to_frame().unstack().round(1)
df_mean.columns = ['大学本科', '研究生及以上', '高中及以下']
df_mean = df_mean[['高中及以下', '大学本科', '研究生及以上']]
df_mean
 

 绘制交互图
 
# 定义一个绘图函数
def draw_pics(data, feature):
    fig, ax = plt.subplots(figsize=(8, 6)) 
    for i in data.index:
        ax.plot(data.columns, data.loc[i,], label = i, marker='o')
        ax.legend()
    ax.set_title("幸福指数平均值")
    ax.set_xlabel(feature, fontdict={'fontsize': 14})
    ax.set_ylabel("平均值", fontdict={'fontsize': 14})
    plt.show()
# 绘制不同的性别在不同的教育程度下的均值变化
draw_pics(data_mean, 'education')
 
 
# 绘制不同的教育程度在不同的性别下的均值变化
draw_pics(data_mean.T, 'gender')
 
 
可以看到受教育程度与性别在对幸福指数的影响上可能存在交互作用。
男性和女性的幸福指数都随受教育程度的增加而增加。但两者的增加趋势不同。男性的受教育程度在高中及以下、大学本科时幸福指数比女性低；但当男性的受教育程度达到研究生及以上时，其幸福指数比女性高。可见，在提高受教育程度增加幸福指数的过程中，男性比女性获益更大。
 
多因素方差分析
 将交互项放入方差分析
 
anova = smf.ols('Index ~ C(gender) + C(education) + C(gender)*C(education)',data = df).fit()
sm.stats.anova_lm(anova, typ=1)
 

 结果显示：
 
性别P = 0.404大于0.05，表明性别对幸福指数没有影响
教育P < 0.001，表明教育对于幸福指数有显著影响
交互项具F（2，52）=4.148，P= 0.021，有统计学意义，表明性别和受教育程度在对幸福指数的影响上存在交互作用
 
事后比较
 不同教育程度的事后比较
 
# 事后多重比较 
sm.stats.multicomp.pairwise_tukeyhsd(groups = df.education, endog=df.Index).summary()
 

 由事后比较可得：
 
大学本科学历 的幸福指数比 高中及以下学历 的人高4.13（95%CI：1.55-6.72），P=0.001；
研究生及以上学历 的幸福指数比 高中及以下学历 的人高21.72（95%CI：19.17-24.28），P<0.001；
研究生及以上学历 的幸福指数比 大学本科学历 的人高17.59（95%CI：15.04-20.15），P<0.001。
 
单独效应的解释
 性别的单独效应
 
# 性别的单独效应
gender_pc_df = pd.DataFrame()
for i in df.gender.unique():
    pc = sm.stats.multicomp.pairwise_tukeyhsd(groups = df.query("gender == @i").education, 
                                              endog=df.query("gender == @i").Index).summary()
    pc_df = pd.DataFrame(pc, index = [i] * (df.education.nunique() + 1), )[1:]
    gender_pc_df = gender_pc_df.append(pc_df)

gender_pc_df.columns = ['group1', 'group2', 'meandiff', 'p-adj', 'lower', 'upper', 'reject']
gender_pc_df
 

 从结果中可得：
 
 
男性中：
 
  
高中及以下学历 的幸福指数评分比 大学本科学历 的低4.50（95%CI：1.58 - 7.42），P = 0.0021；
高中及以下学历 的幸福指数评分比 研究生及以上学历 的低24.44（95%CI：21.60 - 27.29），P = 0.001。
大学本科学历 的幸福指数评分比 研究生及以上学历 的低19.94（95%CI：17.10 - 22.79），P = 0.001。
 
 
女性中：
 
  
高中及以下学历 的幸福指数评分比 大学本科学历 没有显著性差异，P = 0.071；
高中及以下学历 比 研究生及以上学历 的低19.10（95%CI：15.02 - 23.17），P = 0.001。
大学本科学历 的幸福指数评分比 研究生及以上学历 的低15.30（95%CI：11.23 - 19.37），P = 0.001。
 
 
结论
 
1.采用因素方差分析性别和受教育程度对幸福指数的影响，显著性水平选取为P = 0.05。满足方差分析前提条件：用分箱图检验没有异常值，用Shapiro-Wilk检验均满足正态性，用Levene方差齐性检验均满足等方差性（P = 0.286）。
2.本研究中，性别对幸福指数没有影响（P = 0.404），教育对于幸福指数有显著影响（P < 0.001），性别和受教育程度在对幸福指数的影响上存在交互作用（P = 0.021）。
3.不同性别中不同教育水平对于幸福指数的影响：
 男性中：
 - 高中及以下学历 的幸福指数评分比 大学本科学历 的低4.50（95%CI：1.58 - 7.42），P = 0.0021；
 - 高中及以下学历 的幸福指数评分比 研究生及以上学历 的低24.44（95%CI：21.60 - 27.29），P = 0.001。
 - 大学本科学历 的幸福指数评分比 研究生及以上学历 的低19.94（95%CI：17.10 - 22.79），P = 0.001。
 女性中：
 - 高中及以下学历 的幸福指数评分比 大学本科学历 没有显著性差异，P = 0.071；
 - 高中及以下学历 比 研究生及以上学历 的低19.10（95%CI：15.02 - 23.17），P = 0.001。
 - 大学本科学历 的幸福指数评分比 研究生及以上学历 的低15.30（95%CI：11.23 - 19.37），P = 0.001。
 
参考链接：两因素方差分析-SPSS教程点击链接

 
在社会科学研究中，主要的多变量分析方法包括多变量方差分析（Multivariate analysis of variance，MANOVA）、主成分分析（Principal component analysis）、因子分析（Factor analysis）、典型相关（Canonical correlation analysis）、聚类分析（Cluster analysis）、判别分析（Discriminant analysis）、多维量表分析（Multidimensional scaling），以及近来颇受瞩目的验证性因子分析(Confirmatory factor analysis )或线性结构模型（LISREL）与逻辑斯蒂回归分析等，以下简单说明这些方法的观念和适用时机。
 
 
 
1 多变量方差分析
 
MANOVA适用于同时探讨一个或多个自变量与两个以上因变量间因果关系的统计方法，依照研究者所操作自变量的个数，可以分为单因素（一个自变量）或多因素（两个以上自变量）MANOVA。进行多变量方差分析时，自变量必须是离散的定类或定序变量，而因变量则必须是定距以上层次的变量。
 
 
 
2 主成分分析
 
主成分分析的主要功能在分析多个变量间的相关，以建构变量间的总体性指标（overall indicators）。当研究者测量一群彼此间具有高度相关的变量，则在进行显著性检验钱，为避免变量数过多，造成解释上的复杂与困扰，常会先进行主成分分析，在尽量不丧失原有信息的前提下，抽取少数几个主成分，作为代表原来变量的总体性指标，达到资料缩减（data reduction）的功能。进行主成分分析时，并无自变量和因变量的区别，但是所有的变量都必须是定距以上层次变量。
 
 
 
3 因子分析
 
因子分析与主成分分析常被研究者混用，因为二者的功能都是通过对变量间的相关分析，以达到简化数据功能。但不同的是，主成分分析是在找出变量间最佳线性组合（linear combination）的主成分，以说明变量间最多的变异量；至于因子分析，则在于找出变量间共同的潜在结构（latent structure）或因子，以估计每一个变量在各因子上的负荷量（loading）。进行因子分析时，并无自变量和因变量的区分，但是所有变量都必须是定距以上层次变量。
 
 
 
4 典型相关
 
典型相关可视为积差相关或多元回归分析的扩展，主要功能在分析两个变量间的相关。进行多元回归分析的目的，是在分析一个或多个自变量与一个因变量间的关系，而典型相关中因变量也可以是多个；也就是说，典型相关的目的在于通过计算得到两个变量线性组合的加权系数。以使（maximum）两个变量间的相关达到最大化。进行典型相关时，并无自变量和因变量的区分，但是所有变量都必须是定距以上层次变量。
 
 
 
5 聚类分析
 
聚类分析的主要功能在进行分类（classification），当研究者有观测值时，常会根据观测值的相似性或差异性进行分类，以形成几个性质不同的类别，简化解释的工作。也就是说，聚类分析根据对变量进行测量的观察值进行分类，以达到组内同质、组间异值的目的。其次，聚类分析完成后，通常可以进行判别分析，以识别分类的效度。当然，在某些时候也可以对变量进行分类（此功能类似因子分析，因此多采用因子分析解决问题）。进行聚类分析时，并无自变量和因变量的区分，但是所有变量都必须是定距以上层次变量。
 
 
 
6 判别分析
 
判别分析是多变量分析中应用相当广泛的统计方法，它可以用来对样本进行分类的工作；也可以用来了解不同类别样本在某些变量上的差异情形；同时也可以根据不同类别的样本在某些变量的实际表现，用来预测新的样本属于某一类别的概率。因此，在行为科学中，常见的研究者单独使用判别分析，建立判别函数（discriminant function），以对新样本进行预测；或是多变量方差分析的检验值达到显著性水平后，比较不同组别样本在因变量平均数的差异情形；或是聚类分析后，检验聚类分析的正确性。进行判别分析时，自变量是定距以上层次变量，至于因变量通常是离散变量。
 
 
 
7 多维量表分析
 
多维量表分析基本上也是一种分类的统计方法，他在市场上普遍被应用。当研究者想要解释一群受试者（例如消费者）对一组客体（例如商品）在某些变量上相似性的测量中所包含的信息，此时多维量表分析就是一个相当适用的方法。研究者只要将这一组客体在变量上的测量值转化成多维度的几何表征，就能够将这些客体有效地显示在这个几何空间中，达到分类的目的，同时也可以进一步解释这些几何表征所代表的潜在结构或意义。进行多维量表分析时，并无自变量和因变量的区分，同时变量可以是等距以上变量，也可以是定类或定序变量。
 
 
 
8 线性结构方程
 
线性结构方程是一个相当具有变通与弹性的统计方法，随着研究者对变量间关系界定的差异，LISREL的常见名称包括协方差结构分析，潜变量分析、线性结构模型或验证性因子分析。LISREL可视为多元回归分析与因子分析两个方法论的整合模型，让研究者可以探讨变量间的线性关系（回归分析），并对可测量显变量与不可测量的潜变量见（因子分析）的因果模型作假设检验。
 
 
 
9 逻辑斯蒂回归分析
 
逻辑斯蒂回归可视为传统多元回归分析的一个特列。它和多元回归分析一样，都具有解释自变量与因变量之间的关系，并可进行预测。所不同的是在进行多元回归分析时，包括自变量与因变量都必须是定距以上层次变量；但在进行逻辑斯蒂回归分析时，自变量仍是定距以上层次变量，因变量则是二分的定类变量或多分定类变量或定序变量。
 
 
 
10 对数线性方程
 
在基本统计学中，当研究者面对探讨两个定类或定序变量间关系的研究问题时，都是以卡方检验来进行假设检验。当问题的性质是探讨两个定类变量间是否独立或是关联强度时，是以卡方独立性检验来进行假设检验。进行卡方独立性检验时，研究者必须将样本在两个定类变量上的反应，建立二维列联表（contingency table），以进一步根据列联表中各单元格（cell）的次数反应，进行显著性检验。但当研究者面对三个或三个以上的定类变量时，所建立的多元列联表间变量关联的分析，卡方独立性检验将无法解决这样的问题，此时适合的方法就是对数线性模型。利用对数线性模型来解决多元列联表的问题的目的，主要就在于探讨构成列联表的多个定类变量间的关系，进而在精简原则下构建拟合的解释模型，并根据所建立的模型估计单元格参数值，以了解各变量效果对单元格次数的影响。
 
 
 
11 Logit对数线性模型
 
在对数线性模型中，多个定类变量间是互为因果的关系（即相关关系），并无自变量与因变量的区分，研究目的在于探讨变量间的关联强度和性质。但有时研究者会面临变量间有自变量和因变量的区分的情境。在基本统计学中，当研究者面对的问题性质是两个定类变量间有自变量和因变量的区别，目的在于探讨两个变量间的因果关系时，多是以卡方齐性检验来进行假设检验。但自变量个数在两个以上时，卡方齐性检验就不再适用，而必须改用logit对数线性模型方法来对数据进行分析。Logit对数线性模型的功能与多元回归分析相当类似，都可以用来探讨与解释因变量与自变量间的关系，但不同的是，多元回归分析的变量都是定距以上层次变量，通常以最小二乘法进行模型估计与检验；logit对数线性模型的变量都是定类变量，通常以最大似然估计法进行模型估计与检验。
 
 
 
·END·
 

转载自：http://www.360doc.com/content/18/0309/10/52857582_735589820.shtml

在实际应用中，一个实验的指标往往受到多个因素的影响。
 
例如饮料的销量有可能受到销售地区或者饮料颜色的影响。在方差分析中，若把饮料的颜色看做影响销量的因素A，把销售地区看做影响因素B。同时对因素A和因素B进行分析，就称为双因素方差分析。
 
a	b	c
a1	b1	20
a1	b2	22
a1	b3	24
a1	b4	16
a1	b5	26
a2	b1	12
a2	b2	10
a2	b3	14
a2	b4	4
a2	b5	22
a3	b1	20
a3	b2	20
a3	b3	18
a3	b4	8
a3	b5	16
a4	b1	10
a4	b2	12
a4	b3	18
a4	b4	6
a4	b5	20
a5	b1	14
a5	b2	6
a5	b3	10
a5	b4	18
a5	b5	10

 
from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm

formula = 'c~ a + b '
anova_results = anova_lm(ols(formula,df).fit())
print(anova_results)

 
            df  sum_sq  mean_sq         F    PR(>F)
a          4.0  335.36    83.84  3.874307  0.021886
b          4.0  199.36    49.84  2.303142  0.103195
Residual  16.0  346.24    21.64       NaN       NaN

 
检验的结论：
 
因素A的p值0.021886<0.05,拒绝原假设，说明饮料颜色对销量有显著影响；而因素B的p值0.103195>0.05,不能拒绝原假设，因此没有充分的理由说明销售地区对销量有显著影响。
 
然而，我们知道了颜色对销量有显著影响，那么是哪种颜色呢？
 使用tukey方法对颜色进行多重比较
 
from statsmodels.stats.multicomp import pairwise_tukeyhsd
print(pairwise_tukeyhsd(df['c'], df['a']))
 
Multiple Comparison of Means - Tukey HSD,FWER=0.05
==============================================
group1 group2 meandiff  lower    upper  reject
----------------------------------------------
  1      2      -9.2   -19.0855  0.6855 False 
  1      3      -5.2   -15.0855  4.6855 False 
  1      4      -8.4   -18.2855  1.4855 False 
  1      5     -10.0   -19.8855 -0.1145  True 
  2      3      4.0    -5.8855  13.8855 False 
  2      4      0.8    -9.0855  10.6855 False 
  2      5      -0.8   -10.6855  9.0855 False 
  3      4      -3.2   -13.0855  6.6855 False 
  3      5      -4.8   -14.6855  5.0855 False 
  4      5      -1.6   -11.4855  8.2855 False 
----------------------------------------------
 
结果说明：1和5的reject=True，说明这两种颜色有显著性差异