R 机器学习预测时间序列模型
机器学习在时间序列数据上应用
随着疫情的变化,急性传染病数据经常会随时间变化,我们通过对每天传染病的记录,就形成了时间序列数据,周期可以是天,周,月,年。目前我们经常会用到ARIMA来预测疾病在未来的变化趋势。
image.png
但是随着机器学习的广泛应用,在时间序列上,也可以采用机器学习发方法去预测,结果比传统的ARIMA
EST
更加快速,简洁,准确。
这次将要介绍关于的时间序列预测的Modeltime
包,旨在加快模型评估,选择和预测的速度。modeltime
通过将tidymodels
机器学习软件包生态系统集成到简化的工作流中以进行tidyverse
预测来实现此目的。modeltime
结合了机器学习模型,经典模型和自动化模型等。
主要优点:
简化数据建模预测流程,包括数据建模,评估,预测及输出
- 预测的系统工作流程。
modeltime_table(
),modeltime_calibrate()
和modeltime_refit()
- 结合Tidymodels以期加入机器学习算法。如
XGBoost
,GLMnet
,Stan
,Random Forest
等 - 改进传统时间序列模型。如
arima_boost()
,prophet_boost()
1.数据
我们选取bike_sharing_daily时间序列数据集,其中包括自行车每日的使用数据。
这里只需要日期与当日的使用量“date”
and “value”
。然后可以简单绘制一下。
注意这里的时间序列是tibble格式。
install.packages("modeltime")
devtools::install_github("business-science/timetk")
library(tidyverse)
library(tidymodels)
library(modeltime)
library(timetk)
library(lubridate)
# data
bike_transactions_tbl <- bike_sharing_daily %>%
select(dteday, cnt) %>%
set_names(c("date", "value"))
bike_transactions_tbl
# plot
bike_transactions_tbl %>%
plot_time_series(date, value, .interactive = FALSE)
2.模型构建
下面继续对数据进行处理,在建模之前,我们需要将数据分成train
与test
。
使用time_series_split()
来分割我们的数据,assess = "3 months"
来确定后三个月为test数据集,cumulative = TRUE
指定前面部分为train。
# Split your time series into training and testing sets
splits <- bike_transactions_tbl %>%
time_series_split(assess = "3 months", cumulative = TRUE)
splits
<Analysis/Assess/Total>
<641/90/731>
# plot
splits %>%
tk_time_series_cv_plan() %>%
plot_time_series_cv_plan(date, value, .interactive = FALSE)
image.png
image.png
2.1经典ARIMA与Prophet
因为model time可以简化数据建模流程,包括简化参数的设置,自动模型通常是已经包含了自动化的建模方法。包括“自动ARIMA”和“自动ETS”功能以及“ Prophet”算法。主要包含三个参数设置:
- Model Spec: 指定预测模型种类(e.g. arima_reg(), prophet_reg())
- Engine: 指定模型
- Fit Model: 加载trian数据
结下来我们尝试建立ARIMA模型与Prophet模型
# ARIMA
model_fit_arima <- arima_reg() %>%
set_engine("auto_arima") %>%
fit(value ~ date, training(splits))
model_fit_arima
parsnip model object
Fit time: 426ms
Series: outcome
ARIMA(0,1,3) with drift
Coefficients:
ma1 ma2 ma3 drift
-0.6106 -0.1868 -0.0673 9.3169
s.e. 0.0396 0.0466 0.0398 4.6225
sigma^2 estimated as 730568: log likelihood=-5227.22
AIC=10464.44 AICc=10464.53 BIC=10486.74
# Prophet
model_fit_prophet <- prophet_reg() %>%
set_engine("prophet", yearly.seasonality = TRUE) %>%
fit(value ~ date, training(splits))
model_fit_prophet
parsnip model object
Fit time: 455ms
PROPHET Model
- growth: 'linear'
- n.changepoints: 25
- seasonality.mode: 'additive'
- extra_regressors: 0
2.2机器学习算法
机器学习模型前面设置比自动化经典模型更为复杂。通常在进行机器学习建模之前,对数据进行预处理,称之为workflow
一般过程如下:
- 创建预处理配方 Preprocessing Recipe
- 创建模型规格 Model Specifications
- 使用工作流将模型规格和预处理相结合,并拟合模型 Fit Model
首先,我将使用配recipe()
创建预处理数据的先前步骤。该过程使用“日期”列创建了我要建模的45个新的列。这些列包含了时间序列的详细信息及傅立叶变化的数据。
recipe_spec <- recipe(value ~ date, training(splits)) %>%
step_timeseries_signature(date) %>%
step_rm(contains("am.pm"), contains("hour"), contains("minute"),
contains("second"), contains("xts")) %>%
step_fourier(date, period = 365, K = 5) %>%
step_dummy(all_nominal())
recipe_spec %>% prep() %>% juice()
有了recipe,我们可以建立机器学习模型了。
为什么需要recipe是因为在tidymodel里面,设置了建立机器学习模型的一套准则,感兴趣可以去:
机器学习模型
这里我们新建了glmnet与RF模型。
主要参数包括:
- Start with a workflow()
- Add a Model Spec: add_model(model_spec_glmnet)
- Add Preprocessing: add_recipe(recipe_spec %>% step_rm(date))
- Fit the Workflow: fit(training(splits))
# Elastic NET model
model_spec_glmnet <- linear_reg(penalty = 0.01, mixture = 0.5) %>%
set_engine("glmnet")
workflow_fit_glmnet <- workflow() %>%
add_model(model_spec_glmnet) %>%
add_recipe(recipe_spec %>% step_rm(date)) %>%
fit(training(splits))
# Random Forest
model_spec_rf <- rand_forest(trees = 500, min_n = 50) %>%
set_engine("randomForest")
workflow_fit_rf <- workflow() %>%
add_model(model_spec_rf) %>%
add_recipe(recipe_spec %>% step_rm(date)) %>%
fit(training(splits))
2.3机器学习结合传统算法
Prophet Boost算法将Prophet与XGBoost结合使用,从而获得了两全其美的效果(即Prophet Automation + Machine Learning)
- 首先使用prophet对单个时间序列建模
- 使用通过预处理配方提供的回归数据(生成的45个新列),并使用XGBoost模型对prophet残差进行回归
model_spec_prophet_boost <- prophet_boost() %>%
set_engine("prophet_xgboost", yearly.seasonality = TRUE)
workflow_fit_prophet_boost <- workflow() %>%
add_model(model_spec_prophet_boost) %>%
add_recipe(recipe_spec) %>%
fit(training(splits))
## [05:22:38] WARNING: amalgamation/../src/learner.cc:480:
## Parameters: { validation } might not be used.
##
## This may not be accurate due to some parameters are only used in language bindings but
## passed down to XGBoost core. Or some parameters are not used but slip through this
## verification. Please open an issue if you find above cases.
workflow_fit_prophet_boost
3.模型调整
image.pngModeltime
工作流程旨在加速模型评估和选择。
现在我们有了几个时间序列模型,让我们对其进行分析,并通过模型时间工作流程预测未来变化趋势。
Modeltime
使用ID来定位我们之前建立的模型,以帮助我们识别模型。
让我们将模型添加到modeltime_table()
中。
[图片上传中...(image.png-ccae46-1593703307979-0)]
model_table <- modeltime_table(
model_fit_arima,
model_fit_prophet,
workflow_fit_glmnet,
workflow_fit_rf,
workflow_fit_prophet_boost
)
model_table
## # Modeltime Table
## # A tibble: 5 x 3
## .model_id .model .model_desc
## <int> <list> <chr>
## 1 1 <fit[+]> ARIMA(0,1,3) WITH DRIFT
## 2 2 <fit[+]> PROPHET
## 3 3 <workflow> GLMNET
## 4 4 <workflow> RANDOMFOREST
## 5 5 <workflow> PROPHET W/ XGBOOST ERRORS
模型校准用于量化误差并估计置信区间。
我们将使用modeltime_calibrate()
函数对test数据进行模型校准。
将生成两个新列(“ .type”和“ .calibration_data”),其中最重要的是“ .calibration_data”。
包括l测试集的实际值,拟合值和残差。
calibration_table <- model_table %>%
modeltime_calibrate(testing(splits))
calibration_table
## # Modeltime Table
## # A tibble: 5 x 5
## .model_id .model .model_desc .type .calibration_data
## <int> <list> <chr> <chr> <list>
## 1 1 <fit[+]> ARIMA(0,1,3) WITH DRIFT Test <tibble [90 × 4]>
## 2 2 <fit[+]> PROPHET Test <tibble [90 × 4]>
## 3 3 <workflow> GLMNET Test <tibble [90 × 4]>
## 4 4 <workflow> RANDOMFOREST Test <tibble [90 × 4]>
## 5 5 <workflow> PROPHET W/ XGBOOST ERRORS Test <tibble [90 × 4]>
4.模型预测
利用校准后的数据,我们可以对test数据进行预测。
使用modeltime_forecast()
生成测试集的预测数据。
使用plot_modeltime_forecast()
绘制预测结果。
calibration_table %>%
modeltime_forecast(actual_data = bike_transactions_tbl) %>%
plot_modeltime_forecast(.interactive = FALSE)
image.png
5.模型评估
接下来,计算预测值与实际值的差异来评估模型好坏。
使用modeltime_accuracy()
生成样本外准确性指标数据。
使用table_modeltime_accuracy()
生成预测结果表格
calibration_table %>%
modeltime_accuracy() %>%
table_modeltime_accuracy(.interactive = FALSE)
image.png
后续还有对最优模型进行预测,使用全部数据作为training,预测未来1年的变化趋势。
参考
Modeltime GitHub Page - Give it a Star if you like it!
Timetk Documentation - Data wrangling, visualization, and preprocessing for time series.
Tidymodels.org - The tidymodels framework is a collection of packages for modeling and machine learning using tidyverse principles.
Introducing Modeltime