Python中含多组间因素的重复测量ANOVA实现方法问询

针对你需要的1个组内因素（time-point）+2个组间因素（Group、Sex）的混合ANOVA，以下是两种可行的Python实现方案：

方案一：用statsmodels的混合效应模型+嵌套模型比较 #

statsmodels的LinearMixedModel可以拟合包含多因素交互的混合效应模型，通过对比嵌套模型的拟合差异，就能得到对应ANOVA效应的检验结果（F值、p值）。

代码示例 #

import pandas as pd
import statsmodels.api as sm
from statsmodels.formula.api import mixedlm
from statsmodels.stats.anova import anova_lm

# 假设数据结构：包含被试ID(subject)、组间因素(Group/Sex)、组内因素(time_point)、因变量(dv)
data = pd.read_csv("your_data.csv")

# 拟合全模型：包含所有主效应和交互效应
full_model = mixedlm(
    formula="dv ~ Group * Sex * time_point",
    data=data,
    groups=data["subject"],
    re_formula="~time_point"  # 考虑被试内time_point的随机斜率和截距
)
full_result = full_model.fit()

# 检验三向交互效应：对比全模型和去掉三向交互的模型
model_no_3way = mixedlm(
    formula="dv ~ Group * Sex + Group * time_point + Sex * time_point",
    data=data,
    groups=data["subject"],
    re_formula="~time_point"
)
result_no_3way = model_no_3way.fit()
anova_3way = anova_lm(full_result, result_no_3way)
print("三向交互效应检验结果：")
print(anova_3way)

# 检验Group×Sex组间交互效应
model_no_group_sex = mixedlm(
    formula="dv ~ Group + Sex + Group * time_point + Sex * time_point",
    data=data,
    groups=data["subject"],
    re_formula="~time_point"
)
result_no_group_sex = model_no_group_sex.fit()
anova_group_sex = anova_lm(full_result, result_no_group_sex)
print("\nGroup×Sex交互效应检验结果：")
print(anova_group_sex)

# 同理可检验主效应（如Group、Sex、time_point）和其他交互效应

原理说明 #

混合ANOVA的各效应本质对应混合效应模型中的固定效应项。通过比较包含目标效应和不包含目标效应的两个嵌套模型的拟合差异（似然比检验转化为F检验），就能得到该效应的统计显著性。

方案二：用pymer4库直接拟合并输出ANOVA结果 #

pymer4是基于R语言lme4库的Python封装，专门简化混合效应模型的操作，支持直接输出多因素混合ANOVA的结果，无需手动对比嵌套模型。

代码示例 #

import pandas as pd
from pymer4.models import Lmer

# 加载数据
data = pd.read_csv("your_data.csv")

# 拟合混合模型：指定随机效应为被试的截距和time_point斜率
model = Lmer("dv ~ Group * Sex * time_point + (1 + time_point|subject)", data=data)
model.fit()

# 直接输出ANOVA风格的效应检验结果
print(model.anova())

优势 #

pymer4的anova()方法会自动计算所有固定效应（主效应、交互效应）的F值、p值，输出格式和传统混合ANOVA结果一致，操作更简洁。

内容的提问来源于stack exchange，提问作者haidan