R中Rob J Hyndman的Arima与Python中SARIMAX结果复现问题

我尝试用Python的SARIMAX函数复现R语言中Rob J Hyndman的Arima函数计算结果，但同一数据集下两者系数无法匹配。以下是原代码对比及修复方案：

原代码对比 #

R代码 #

md1 = Arima(y, xreg = as.matrix(df[,c('x1', 'x2'), drop = FALSE]), order = c(1,0,1), include.mean = TRUE)

原Python代码（存在语法/参数遗漏问题） #

md2 = SARIMAX(y, exog = df[np.array(['x1', 'x2'])], trend='c', mle_regression=True)

核心问题与修复方案 #

• 遗漏ARMA阶数参数：R中明确指定了order=c(1,0,1)（AR(1)、差分0、MA(1)），但原Python代码未设置order=(1,0,1)，这是结构不匹配的核心原因。
• 外生变量格式不一致：R将xreg转为矩阵，Python中需确保exog是二维数组，避免索引错误，可直接用df[['x1','x2']]或转为numpy矩阵。
• 优化器差异：R的Arima默认用L-BFGS优化器，而Python SARIMAX默认用BFGS，指定method='lbfgs'可对齐优化逻辑，减少收敛结果差异。

修正后的Python代码 #

import numpy as np
from statsmodels.tsa.statespace.sarimax import SARIMAX

# 处理外生变量，确保格式正确
exog_data = df[['x1', 'x2']].values
# 对齐R的模型参数设置
md2 = SARIMAX(
    y,
    exog=exog_data,
    trend='c',
    order=(1, 0, 1),
    mle_regression=True,
    method='lbfgs'
)
# 拟合模型（关闭迭代信息输出）
results = md2.fit(disp=False)
# 查看结果
print(results.summary())

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