Granger因果与误差修正模型:Python实现协整分析与VECM全流程
发布时间:2026/9/18 19:00:10 锦皓数字建站

简介这份doc格式的专题资料完整梳理了Granger因果关系与误差修正模型ECM的分析方法面向经济学、金融学及能源经济领域的研究生和科研人员也适合需要处理时间序列数据的实证分析者。文档以中国电力与经济增长关系为案例系统讲解平稳性检验ADF与PP单位根检验、结构断点分析、Granger因果检验、Johansen协整检验以及ECM模型的建立步骤配有模型公式、检验统计量和判定逻辑有助于读者快速掌握从单位根检验到协整建模的完整分析流程。压缩包共1个文件为doc类型大小784KB内容结构紧凑便于在Word中阅读和标注。目前已有156人学习下载适合作为课程论文、课题研究或论文复现的参考资料。1. Granger 和 ECM 不是两套独立方法而是一条因果分析链路两个带趋势的时间序列直接做 OLS 回归 R² 能到 0.95t 值全部显著换个样本区间结论就翻脸——Granger 因果和 ECM误差修正模型这套方法组合就是为这类伪回归问题设计的。Granger 检验回答x 的过去值能否提升对 y 的预测能力它不涉及经济学意义上的因果ECM 则把长期均衡关系和短期偏离修正放进同一个方程两者单独用都会出问题序列不平稳时 Granger 检验的 F 统计量没有标准分布序列不存在协整时 ECM 写出来还是伪回归。所以 2021-2022 年前后主流的计量分析流程是一条固定链路ADF 单位根检验 → 协整检验 → Granger 因果与 ECM/VECM 估计。下面按这条链路推进用 Python statsmodels 给出可复现代码把每一步的参数设置和判读标准讲清楚适合处理金融、宏观或任何带趋势数据的数据分析者。2. 建模前的数据体检平稳性检验与滞后阶选择2.1 ADF 单位根检验先分清 I(0) 和 I(1)Granger 因果检验的核心前提是变量平稳。对带单位根的序列直接做水平 VAR统计量的渐近分布与标准情形不同p 值不可信ECM 的长期均衡项也只有在单整阶数明确时才有意义。所以第一步永远是单位根检验最常用 ADF辅以 KPSS 做交叉验证。import pandas as pd from statsmodels.tsa.stattools import adfuller # 演示数据两列月度序列 y 和 x2015-2022 年约 96 个观测 df pd.read_csv(macro_series.csv, parse_dates[date], index_coldate) y df[y] # 默认带常数项(regressionc)autolag 按 AIC 自动选滞后 res adfuller(y, autolagAIC, maxlag12) print(fADF 统计量: {res[0]:.4f}) print(fp 值: {res[1]:.4f}) print(f使用滞后阶: {res[2]}) print(f临界值: { {k: round(v, 4) for k, v in res[4].items()} })adfuller 返回顺序是统计量、p 值、实际采用的滞后阶、样本数、临界值字典和最优信息准则。判断先看 p 值p 0.05 拒绝存在单位根的原假设序列平稳否则不能拒绝序列非平稳。maxlag 按数据频率给月度数据取 12季度数据取 4 到 8太小会漏掉自相关太大损失自由度样本只有 96 个观测时尤其明显。regression 参数有三个选择c 带常数项默认、ct 带常数和趋势项、n 都不带。序列有明显趋势时用 ct 更稳妥不带趋势项时 ADF 检验功效会明显下降甚至把漂移项误判成随机游走。提示如果 ADF 在 c 和 ct 两种设定下结论相反通常以 ct 为准。趋势项没控制住时单位根检验不可信。只用 ADF 不够它对某些结构性断点序列功效偏低。补一个 KPSS 检验原假设与 ADF 相反是序列平稳from statsmodels.tsa.stattools import kpss kpss_res kpss(y, regressionc, nlagsauto) print(fKPSS 统计量: {kpss_res[0]:.4f}, p 值: {kpss_res[1]:.4f})ADF 不拒绝单位根、KPSS 拒绝平稳I(1) 的证据就一致了两个检验都拒绝说明序列是 I(0)两个都不拒绝多半是样本量太小或序列里有断点需要先处理异常段再重新检验。对 x 重复同样流程后把结论记下来y 和 x 是否都是 I(1)决定了后面能不能走协整路线。2.2 滞后阶AIC/BIC 与残差自相关一起看滞后阶数同时影响 ADF、Granger 检验和 VECM 三处设定值得单独定一次。先在水平 VAR 上跑信息准则from statsmodels.tsa.api import VAR for lag in range(1, 7): var_model VAR(df[[y, x]]).fit(lagslag) ic var_model.info_criteria print(f滞后 {lag}: AIC{ic[aic]:.2f}, BIC{ic[bic]:.2f}, fHQIC{ic[hqic]:.2f}) # 对备选阶数做残差白噪声检验p 值要大于 0.05 var_final VAR(df[[y, x]]).fit(lags3) print(var_final.test_whiteness(6).summary())下表是一组典型输出示意值实际以你自己数据为准滞后阶AICBICHQIC残差 Ljung-Box p 值1-820.4-791.2-808.70.0032-835.1-798.6-820.50.0713-842.8-800.2-826.30.1864-841.5-805.3-828.90.212信息准则给了一个候选区间AIC 和 HQIC 都指向 3 阶BIC 更节俭指向 4 阶。此时不要只按 AIC 最小选——还要看残差2 阶残差的 Ljung-Box p 值只有 0.071说明仍接近有自相关到 3 阶上升到 0.186残差白噪声成立。常见做法是选 3 阶然后额外把 4 阶的估计结果作为稳健性对照。这个水平 VAR 的最优滞后阶记为 p后面 VECM 的差分滞后项要取 p−1statsmodels 的 k_ar_diff 参数就是这个含义别在这里省事。2.3 单整阶数不一致时ECM 直接不成立协整的定义要求所有变量单整阶数一致最常见的组合是两个 I(1) 序列。如果检验结果是 y 是 I(1)、x 是 I(0)两者之间不存在可估计的长期均衡关系常规做法是把 x 以差分形式放进回归不要硬套 ECM。如果出现 I(2)先做二阶差分再回到常规流程。这里有一个高频误解很多人先对两个序列各差分一次再做 OLS以为这样平稳了但差分回归丢失了水平信息——如果 y 和 x 真的协整正确建模用的是水平残差不是差分残差。正确顺序永远是先确认 I(1)再做协整检验最后依据协整结论决定模型形式。做完这一章的数据体检得到两个明确结论y 和 x 都是 I(1)VAR 最优滞后 p 3。接下来才轮到协整检验。3. 协整检验先证明长期均衡存在再做 ECM3.1 Engle-Granger 两步法及其临界值陷阱Engle-Granger 两步法思路直白第一步y 对 x 做水平回归得到残差第二步检验残差是否平稳。残差是 I(0) 说明 y 和 x 的线性组合消除了单位根存在长期均衡。这里最关键的一点是残差是估计出来的OLS 有让残差看起来更平稳的倾向直接用 ADF 的高斯临界值会高估协整关系。import statsmodels.api as sm from statsmodels.tsa.stattools import adfuller # 第一步带常数项的水平回归 X_level sm.add_constant(x) long_run sm.OLS(y, X_level).fit() resid long_run.resid # 第二步直接对残差跑 ADF——注意这个 p 值只能作为参考 adf_resid adfuller(resid, autolagAIC, maxlag12) print(f残差 ADF 统计量: {adf_resid[0]:.4f}, p 值: {adf_resid[1]:.4f})这段代码里 adfuller 的 p 值不能当最终结论。正确做法是用 statsmodels 的 coint()它内置 Engle-Granger 专用临界值MacKinnon 响应面而不是标准 ADF 的临界值表from statsmodels.tsa.stattools import coint t_stat, p_value, crit coint(y, x, autolagAIC, maxlag12) print(fEG 协整检验 t 统计量: {t_stat:.4f}, p 值: {p_value:.4f}) print(f临界值 1% / 5% / 10%: {[round(c, 4) for c in crit]})coint() 的趋势参数 trend 默认 c对应协整关系含常数项。p 值小于 0.05 拒绝无协整原假设说明 y 和 x 之间存在长期均衡。EG 两步法适合两个变量、最多一个协整向量变量多于两个或可能存在多条协整关系时要用 Johansen 检验。下面这张表总结了两个方法的边界选型时直接对照维度Engle-Granger 两步法Johansen 检验适用变量数两个两个及以上协整向量数量最多一个无法检验多条能确定协整秩 r长期参数估计第二步残差 ADF 需要专用临界值特征向量直接给出主要局限第一步估计误差会污染第二步标准误小样本下倾向高估协整秩3.2 Johansen 检验迹统计量与最大特征值统计量Johansen 在 VAR 框架下同时估计所有可能的协整向量输出两类统计量迹统计量trace检验协整秩 ≤ r对r最大特征值统计量检验秩 r对秩 r1。statsmodels 的实现from statsmodels.tsa.vector_ar.vecm import coint_johansen # det_order0 表示协整关系内含常数项k_ar_diff 取 VAR 最优阶减 1 joh_data df[[y, x]].dropna() joh coint_johansen(joh_data, det_order0, k_ar_diff2) print(迹统计量:, joh.trace_stat.round(4)) print(迹统计量临界值(90%,95%,99%):\n, joh.trace_stat_crit_vals.round(4)) print(最大特征值统计量:, joh.max_eig_stat.round(4)) print(最大特征值统计量临界值:\n, joh.max_eig_stat_crit_vals.round(4))判读从第一行开始trace_stat[0] 大于 95% 临界值则拒绝不存在协整再看 trace_stat[1] 是否拒绝至多一条协整关系。两个变量的系统r 1 是最常见结论也就是存在一条长期均衡关系。示意判读如下具体临界值以 statsmodels 输出为准原假设迹统计量5% 临界值结论r 0无协整32.7415.41拒绝存在协整r ≤ 1至多一条3.023.84不能拒绝止于 r 1两个统计量冲突时一般以迹统计量为准但小样本下 Johansen 倾向高估协整秩样本量小于 80 时要多依赖最大特征值统计量的保守结论。k_ar_diff 的取值必须来自 2.2 节定出的水平 VAR 阶数p 3 时这里填 2。如果明明定了 p 却在这里填了 3等于多估了一阶检验结果会系统性偏移。3.3 协整向量归一化与长期参数Johansen 输出的特征向量evec不唯一需要归一化。对 r 1、两个变量的情况把 y 的系数归一为 1另一个系数就是长期关系 y θ·x c 里的 θbeta joh.evec[:, 0] # 第一个协整向量对应最大特征值 beta_norm beta / beta[0] # 归一化令 y 的系数为 1 theta -beta_norm[1] # 长期系数 print(fJohansen 长期系数 θ {theta:.4f}) # 与 EG 两步法的斜率对比两者应该接近 print(fEG 长期系数 θ {long_run.params[x]:.4f})两个来源的 θ 应该接近。如果差异大先检查 k_ar_diff 和 det_order 是否一致再看数据里有没有结构性断点比如 2020-2021 年那一段异常波动会使两步法的长期参数都受影响。θ 的经济含义是长期均衡比例x 变动 1 单位y 在长期中向 θ 单位调整。注意 θ 本身不构成因果证据它只描述长期共变关系——因果方向留给下一步 ECM 的调整系数去识别。4. Granger 因果与 ECM 估计短期预测力与长期调整分开读4.1 用 grangercausalitytests 跑短期 Granger 因果先看基础工具。grangercausalitytests 一次输出多个滞后阶的结果默认同时算 F 检验、卡方检验和似然比检验报告时通常只写 F 检验那一行。from statsmodels.tsa.stattools import grangercausalitytests # 重要原假设是第二列不 Granger 引起第一列 # 想检验 x → y就把 y 放第一列、x 放第二列 gc_data df[[y, x]].dropna() results grangercausalitytests(gc_data, maxlag3, verboseFalse) for lag, res in results.items(): f_stat, f_pval, df_num, df_den res[0][ssr_ftest] print(f滞后 {lag}: F{f_stat:.4f}, p{f_pval:.4f})参数顺序是最容易错的地方把 x、y 放反检验的就是y 是否 Granger 引起 x结论完全反过来。maxlag 取 2.2 节定出的 VAR 最优水平滞后阶 p不要只跑一阶因果可能滞后多期才显现。这里有个比参数顺序更深的坑grangercausalitytests 对 I(1) 序列直接用水平值跑F 统计量不再服从标准分布p 值偏小。数据确认协整后正确的短期 Granger 检验应该放到 VECM 里去见 4.3 节。4.2 单方程 ECM误差修正项的构造与解读协整关系存在时误差修正模型把长期均衡和短期动态放在一个方程里Δy_t c β·Δx_t γ·ECM_{t−1} ε_t其中 ECM_{t−1} y_{t−1} − α − θ·x_{t−1} 是上一期对长期均衡的偏离α 是长期回归的常数项。Python 里按 Engle-Granger 两步法估计# 第一步的长期回归已经得到 long_run 和 theta ecm long_run.resid # y_t - (α θ·x_t) ecm_lag ecm.shift(1) # 取上一期的均衡偏离 # 构建差分样本并去掉缺失值 dy y.diff() dx x.diff() ecm_df pd.concat([dy, dx, ecm_lag], axis1, keys[dy, dx, ecm_lag]).dropna() X_ecm sm.add_constant(ecm_df[[dx, ecm_lag]]) ecm_model sm.OLS(ecm_df[dy], X_ecm).fit() print(ecm_model.summary())三个系数的读法如下表系数对应变量含义判读标准c常数项差分回归漂移项一般不解释βdx短期乘子x 当期变动 1 单位对 y 当期变动的影响γecm_lag调整速度必须显著为负绝对值越大回拉越快γ −0.15 表示每期修复上一期偏离的 15%半衰期 ln(2)/0.15 ≈ 4.6 期γ 为正说明系统发散长期关系存疑γ 不显著说明 y 这一侧不承担向均衡调整的压力调整可能全部发生在 x 那边。遇到 γ 不显著把被解释变量和解释变量对调用同样的长期关系重新估计一次往往另一侧的调整系数就显著了——这种不对称性本身就是 Granger 因果方向的证据承担调整压力的是被动方。4.3 VECM双向调整与短期 Granger 因果的正式检验两个变量的调整系数都显著或者变量多于两个时单方程 ECM 有效率损失用 VECM 同时估计所有方程from statsmodels.tsa.vector_ar.vecm import VECM # k_ar_diff2 对应水平 VAR 3 阶deterministicci 常数在协整关系内 vecm_model VECM(df[[y, x]], k_ar_diff2, coint_rank1, deterministicci) vecm_res vecm_model.fit() print(协整向量 β:, vecm_res.beta) print(调整系数 α:, vecm_res.alpha) print(α 的 p 值:, vecm_res.pvalues_alpha)vecm_res.beta 是 (2, 1) 的协整向量vecm_res.alpha 是 (2, 1) 的调整系数矩阵两行分别对应 y 方程和 x 方程。x 方程的 α 不显著时称 x 弱外生说明 x 不响应均衡偏离此时单方程 ECM 是有效率的直接报 4.2 节的结果就够了两个 α 都显著说明调整是双向的必须报完整 VECM。coint_rank 来自 Johansen 检验结论不要凭空设置k_ar_diff 仍是 p−1。如果你用 EViews 复现对应菜单是 Estimate VAR → Vector Error Correction统计量与此处一致。VECM 下的 Granger 因果分两段短期因果看差分滞后项的联合显著性长期因果看误差修正项对应的调整系数显著性。# 短期x 的差分滞后项对 y 方程是否联合显著 gc_res vecm_res.test_granger_causality(caused0, causing1, kindf) print(f短期 Granger 因果: F{gc_res.stat:.4f}, p{gc_res.pvalue:.4f}) # 长期y 方程调整系数 α_y 的显著性 alpha_y vecm_res.alpha[0, 0] p_alpha_y vecm_res.pvalues_alpha[0, 0] print(f长期调整: α_y{alpha_y:.4f}, p{p_alpha_y:.4f})caused0 指 y 方程causing1 指 x 方程原假设是 x 不 Granger 引起 y。短期 p 值显著说明 x 的短期变动对 y 有预测力α_y 显著说明 x 通过长期均衡渠道对 y 发挥作用。实际数据里短期显著、长期不显著和长期显著、短期不显著都常见分别对应短期冲击传导和长期均衡纠偏两种机制报告时务必分开表述不要合成一句存在 Granger 因果。5. 报告落地与稳健性验证半衰期、残差诊断和三个翻车点5.1 调整速度、半衰期与长期关系一起解读拿到 ECM 估计结果后报给读者的不只是显著与否还要让调整速度可感知。半衰期 h ln(2)/|γ| 是最直观的指标月度数据 h 是月数季度数据是季度数不要混着说。一组典型解读估计结果数值解读γ −0.15p 0.008每期修复 15%半衰期约 4.6 期回拉较快两年内基本收敛γ −0.05p 0.21方向正确但太小调整弱且不显著y 不承担纠偏γ 0.08p 0.04正且显著系统发散检查设定或数据断点注意两步法里 γ 的标准误没有计入第一步估计长期参数的误差通常偏小正式报告要注明这一点需要用 Johansen 支撑长期参数、用 VECM 的 α 和标准误做最终结论时直接引用 vecm_res.stderr_alpha。5.2 稳健性检查换滞后、换确定性设定、残差诊断ECM 的稳健性检查有固定套路滞后阶在 p−1 和 p1 之间切换Johansen 的 det_order 在 0 和 1 之间切换删掉样本两端异常段比如 2020 年前后的断点重跑对所有模型做残差诊断。诊断代码可以直接固化到脚本里from statsmodels.stats.diagnostic import acorr_ljungbox, het_arch from scipy import stats resid_ecm ecm_model.resid print(残差自相关 Ljung-Box:, acorr_ljungbox( resid_ecm, lags[6, 12], return_dfTrue)) print(残差正态性 JB: p , stats.jarque_bera(resid_ecm)[1]) print(残差 ARCH: p , het_arch(resid_ecm, nlags5)[1])残差存在自相关说明滞后阶不足存在 ARCH 效应说明条件方差有聚集JB 正态性在金融数据里经常拒绝只要前两项通过正态性稍差可以接受用稳健标准误兜底即可。5.3 三个高频翻车点用标准 ADF 临界值判协整残差得到伪协整结论。判断协整只能用 EG 专用临界值coint() 内部处理或 Johansen 分布表。确认协整后还在差分序列上单独跑 Granger 检验把长期因果丢掉。长期信息在误差修正项里不在差分项里。把 EG 第一步回归的 R² 当卖点。协整回归的 R² 再高也不说明任何问题能说明问题的只有残差平稳性和 γ或 α的符号与显著性。报告正文的固定格式先给 ADF 表说明 I(1)再给 EG/Johansen 检验说明协整秩为 1最后给 ECM 表或 VECM 表系数表里必须同时给出长期参数 θ 和调整系数 γ或 α缺一个读者都没法判断长期关系和纠偏机制是否成立。按这个顺序写完审稿人或技术评审能直接复现每一步。本文还有配套的精品资源点击获取
锦
锦皓数字建站
深耕本土企业品牌数字化升级,专注原创端正雅致商务官网,从视觉设计到稳定运维全程保驾护航。