Python含固定效应虚拟变量的回归性能远逊于Stata的技术求助

阿华AIGC实验室

2026-5-21

解决Python中带组固定效应的OLS回归性能瓶颈问题

我之前也踩过这个坑——当给回归模型加入组固定效应后，Python的statsmodels OLS直接慢到跑不动，而Stata用areg或者reghdfe几秒就出结果。这根本不是Python不行，而是两者处理固定效应的底层逻辑差了十万八千里，咱们来一步步解决这个问题：

为什么差距这么大？

Stata的areg（或者更常用的reghdfe）用的是组内去均值法：先把每个组里的所有变量（因变量+自变量）都减去组内均值，相当于把固定效应的影响直接“吸收”掉，之后只需要对去均值后的变量做普通回归。这种方法的计算复杂度极低，尤其是当组数多但每组观测数不多时，效率碾压直接生成虚拟变量。
而如果你在Python里直接用smf.ols('income ~ height + Number_children + C(group_id)', data=humans)，本质是把每个组都转换成一个虚拟变量加入模型。假设你有1万个组，那设计矩阵会凭空多出1万个列，矩阵求逆的复杂度直接飙升到O((K+G)^3)（K是自变量数，G是组数），20万条数据加上这么多变量，不卡才怪。

具体解决方案

方案1：用linearmodels的PanelOLS（最推荐）

linearmodels是专门为面板数据设计的库，处理固定效应的逻辑和Stata完全一致，性能拉满。先安装库：

pip install linearmodels

然后写代码：

import linearmodels as plm

# 先给数据设置面板索引：组ID + 观测标识（如果是时间序列面板就用时间）
humans_panel = humans.set_index(['group_id', 'obs_id'])

# 用EntityEffects指定组固定效应，还可以加TimeEffects处理时间固定效应
model = plm.PanelOLS.from_formula(
    'income ~ height + Number_children + EntityEffects',
    data=humans_panel
)
# 还可以指定聚类标准误，和Stata的cluster()参数对应
results = model.fit(cov_type='clustered', cluster_entity=True)
print(results.summary())

方案2：手动实现组内去均值（无需额外库）

如果不想装新库，自己写去均值逻辑也能解决问题：

import statsmodels.api as sm
import pandas as pd

# 对所有变量按组去均值
demeaned_data = humans.groupby('group_id').transform(lambda x: x - x.mean())

# 注意去均值后不需要截距项，因为固定效应已经被吸收了
X = sm.add_constant(demeaned_data[['height', 'Number_children']], has_constant='add')
y = demeaned_data['income']

model = sm.OLS(y, X).fit()
print(model.summary())

⚠️ 注意：如果有组只有1条观测，去均值后变量会变成0，建议提前过滤掉这类组，避免影响结果。

方案3：用pyreghdfe复刻Stata的reghdfe

如果你习惯Stata的reghdfe语法，可以用Python版的pyreghdfe，用法几乎一模一样：

pip install pyreghdfe

代码示例：

from pyreghdfe import regress

results = regress(
    data=humans,
    formula='income ~ height + Number_children + absorb(group_id)',
    cluster=['group_id']  # 支持多维度聚类
)
print(results.summary())