python - 在 python/scikit/numpy 中替代 r 的指数平滑状态空间模型

Question

在 R 中，我们有一个很好的预测模型，例如:

ets(y, model="ZZZ", damped=NULL, alpha=NULL, beta=NULL, gamma=NULL, 

phi=NULL, additive.only=FALSE, lambda=NULL, 

lower=c(rep(0.0001,3), 0.8), upper=c(rep(0.9999,3),0.98), 

opt.crit=c("lik","amse","mse","sigma","mae"), nmse=3, 

bounds=c("both","usual","admissible"), ic=c("aicc","aic","bic"),

restrict=TRUE, allow.multiplicative.trend=FALSE, use.initial.values=FALSE, ...)

在此方法中，如果我们分配任何变量，它会自动获取季节类型、趋势和错误类型，如 model="ZZZ"/"AMA"/"MMZ"并且某些因素会自动调整以获得准确的结果。

• 在 python 中我们是否有任何类似于 ets 的东西 Pandas /numpy/scipy/scikit?

  根据我的研究:
  Ewma pandas 中类似，但我们需要将所有参数硬编码为固定参数。
  在 Holtwinter 中，我们需要为所有趋势和季节类型编写详细的方法。

• 所以我们有任何现成的函数来代替它 数据框作为输入并提供预测值，无需编写 我们自己有任何参数的内部函数吗？

• 任何微调回归模型 scikit/statsmodels？

Answer 1

经过一番搜索，我没有找到任何看起来真正有希望作为 python 的 ets 替代品的东西。不过也有一些尝试:StatsModels和 pycast's Forecasting methods ，您可以检查它们是否适合您的需求。

您可以用来解决缺少的实现的一个选项是使用 subprocess 从 python 运行 R 脚本。模块。关于如何做到这一点有一篇非常好的文章here .

为了以后做:

1. 您需要创建一个 R 脚本(例如 my_forecast.R)，它将 计算(使用 ets)并打印文件或 stdout 上的预测(使用 cat() 命令)，以便在脚本之后使用它们 运行。

2. 您可以从 python 脚本运行 R 脚本，如下所示:

   import subprocess
   
   # You need to define the command that will run the Rscript from the subprocess
   command = 'Rscript'
   path2script = 'path/to/my_forecast.R'
   cmd = [command, path2script]
   
   # Option 1: If your script prints to a file
   subprocess.run(cmd)
   f = open('path/to/created/file', 'r')
   (...Do stuff from here...)
   
   # Option 2: If your script prints to stdout
   forecasts = subprocess.check_output(cmd, universal_newlines=True)
   (...Do stuff from here...)
   

   您还可以在 cmd 中添加参数，Rscript 将使用这些参数作为命令行参数，如下所示:

   args = [arg0, arg1, ...]
   
   cmd = [command, path2script] + args 
   Then pass cmd to the subprocess
   

---

编辑:

我找到了有关 Holt-Winters 预测的一系列示例文章:part1 , part2和 part3 .除了这些文章中易于理解的分析之外，Gregory Trubetskoy(作者)还提供了他开发的代码:

初步趋势:

def initial_trend(series, slen):
    sum = 0.0
    for i in range(slen):
        sum += float(series[i+slen] - series[i]) / slen
    return sum / slen

# >>> initial_trend(series, 12)
# -0.7847222222222222

初始季节性组件:

def initial_seasonal_components(series, slen):
    seasonals = {}
    season_averages = []
    n_seasons = int(len(series)/slen)
    # compute season averages
    for j in range(n_seasons):
        season_averages.append(sum(series[slen*j:slen*j+slen])/float(slen))
    # compute initial values
    for i in range(slen):
        sum_of_vals_over_avg = 0.0
        for j in range(n_seasons):
            sum_of_vals_over_avg += series[slen*j+i]-season_averages[j]
        seasonals[i] = sum_of_vals_over_avg/n_seasons
    return seasonals

# >>> initial_seasonal_components(series, 12)
# {0: -7.4305555555555545, 1: -15.097222222222221, 2: -7.263888888888888,
#  3: -5.097222222222222,  4: 3.402777777777778,   5: 8.069444444444445,  
#  6: 16.569444444444446,  7: 9.736111111111112,   8: -0.7638888888888887,
#  9: 1.902777777777778,  10: -3.263888888888889, 11: -0.7638888888888887}

最后是算法:

def triple_exponential_smoothing(series, slen, alpha, beta, gamma, n_preds):
    result = []
    seasonals = initial_seasonal_components(series, slen)
    for i in range(len(series)+n_preds):
        if i == 0: # initial values
            smooth = series[0]
            trend = initial_trend(series, slen)
            result.append(series[0])
            continue
        if i >= len(series): # we are forecasting
            m = i - len(series) + 1
            result.append((smooth + m*trend) + seasonals[i%slen])
        else:
            val = series[i]
            last_smooth, smooth = smooth, alpha*(val-seasonals[i%slen]) + (1-alpha)*(smooth+trend)
            trend = beta * (smooth-last_smooth) + (1-beta)*trend
            seasonals[i%slen] = gamma*(val-smooth) + (1-gamma)*seasonals[i%slen]
            result.append(smooth+trend+seasonals[i%slen])
    return result

# # forecast 24 points (i.e. two seasons)
# >>> triple_exponential_smoothing(series, 12, 0.716, 0.029, 0.993, 24)
# [30, 20.34449316666667, 28.410051892109554, 30.438122252647577, 39.466817731253066, ...

您可以将它们放在一个文件中，例如:holtwinters.py 在具有以下结构的文件夹中:

forecast_folder
|
└── __init__.py
|
└── holtwinters.py

从这里开始，这是一个 python 模块，您可以将其放置在您想要的每个项目结构中，并在该项目中的任何位置使用它，只需导入它。