python 实现可重复双因素方差分析（3）

案例：

一家超市连锁店进行一项研究，以确定超市所在的位置和竞争者的数量对其销售额是否有显著影响。

取显著性水平， 检验：

（1）竞争者的数量对销售额是否有显著影响。

（2）超市的位置对销售额是否显著影响。

（3）竞争者的数量和超市的位置对销售额是否有交互影响。

# 导入相关包
import pandas as pd
import numpy as np
import math
import scipy
from scipy import stats 
import matplotlib.pyplot as plt
from statsmodels.formula.api import ols
from statsmodels.stats.anova import anova_lm

# 自定义函数
def level_n(data, row_name, col_name):
    r_level = list(data[row_name].unique())
    c_level = list(data[col_name].unique())
    tuple_ = [(i,j) for i in r_level for j in c_level]
    row_n = []
    for i, j in tuple_:
        p = data.query("{}=='{}'".format(row_name, i)).query("{}=='{}'".format(col_name, j)).shape[0]
        row_n.append(p)
    m = max(row_n) 
    return m

def SST(Y):
    sst = sum(np.power(Y - np.mean(Y), 2))
    return sst

def SSA(data, row_name, col_name, y_name):
    n = len(data[col_name].unique())
    m = level_n(data, row_name, col_name)
    
    total_avg = np.mean(data[y_name])
    df = data.groupby([row_name]).agg(['mean', 'count'])
    df = df[y_name]
    
    ssa = n * m * sum((np.power(df["mean"] - total_avg, 2)))
    return ssa

def SSRC(data, row_name, col_name, y_name):
    
    m = level_n(data, row_name, col_name)
    
    df_row = data.groupby([row_name]).agg(['mean'])
    row_dict = dict(df_row[y_name]["mean"])
    
    df_col = data.groupby([col_name]).agg(['mean'])
    col_dict = dict(df_col[y_name]["mean"])
    
    df_rc  = data.groupby([row_name, col_name]).agg(['mean'])
    rc_dict = dict(df_rc[y_name]["mean"])
    
    total_avg = np.mean(data[y_name])
    
    sum = 0
    for i in rc_dict:
        rc_avg= rc_dict[i]
        r_avg = row_dict[i[0]]
        c_avg = col_dict[i[1]]
        sum += np.power(rc_avg - r_avg - c_avg + total_avg, 2)
        
    ssrc = m * sum   
        
    return ssrc

def two_way_anova(data, row_name, col_name, y_name, alpha=0.05):
    """可重复双因素方差分析"""
    
    n = len(data)                          # 总观测值数
    k = len(data[row_name].unique())       # 行变量的水平个数
    r = len(data[col_name].unique())       # 列变量的水平个数
    m = level_n(data, row_name, col_name)  # 行变量中每个水平的行数
    
    sst = SST(data[y_name])                         # 总平方和
    ssr = SSA(data, row_name, col_name, y_name)     # 行变量平方和
    ssc = SSA(data, col_name, row_name, y_name)     # 列变量平方和
    ssrc = SSRC(data, row_name, col_name, y_name)   # 交互作用平方和
    sse = sst - ssr - ssc - ssrc                    # 误差平方和
    
    msr = ssr / (k-1)
    msc = ssc / (r-1)
    msrc = ssrc / ((k-1)*(r-1))
    mse = sse / (k*r*(m-1))
    
    Fr = msr / mse
    Fc = msc / mse
    Frc = msrc / mse
    
    pfr = scipy.stats.f.sf(Fr, k-1, k*r*(m-1))
    pfc = scipy.stats.f.sf(Fc, r-1, k*r*(m-1))
    pfrc = scipy.stats.f.sf(Frc, (k-1)*(r-1), k*r*(m-1))
    
    Far = scipy.stats.f.isf(alpha, dfn=k-1, dfd=k*r*(m-1))   # 行 F临界值
    Fac = scipy.stats.f.isf(alpha, dfn=r-1, dfd=k*r*(m-1))   # 列 F临界值
    Farc = scipy.stats.f.isf(alpha, dfn=(k-1)*(r-1), dfd=k*r*(m-1))   # 交互 F临界值
    
    r_square = (ssr+ssc+ssrc) / sst
    
    table = pd.DataFrame({'差异源':[row_name, col_name, '交互作用', '内部(误差)', '总计'],
                          '平方和SS':[ssr, ssc, ssrc, sse, sst],
                          '自由度df':[k-1, r-1, (k-1)*(r-1), k*r*(m-1), n-1],
                          '均方MS':[msr, msc, msrc, mse, '_'],
                          'F值':[Fr, Fc, Frc, '_', '_'],
                          'P值':[pfr, pfc, pfrc, '_', '_'],
                          'F临界值':[Far, Fac, Farc, '_', '_'],
                          'R^2':[r_square, '_', '_', '_', '_']})
    
    return table

# 导入数据
df = pd.read_excel("E:xx业务数据.xlsx", sheet_name='source_05')

# 输出方差分析结果
two_way_anova(df, '超市位置', '竞争者数量', '销售额', alpha=0.01)

根据以上方差分析结果解释如下：

（1）竞争者数量：P-value=1.57485e-05 <  (或 F值>F临界值=4.71805），拒绝原假设。表明竞争者的数量对销售额有显著影响。

（2）超市位置：P-value=9.17534e-08    <  (或 F值>F临界值=5.61359），拒绝原假设。表明超市的位置对销售额有显著影响。

（3）交互作用：P-value=0.0160501       >  (或 F值>F临界值=3.66672），不拒绝原假设。没有证据表明竞争者的数量和超市的位置对销售额有交互影响。