如何模拟 ARIMA 的分布值以绘制成扇形图
How to simulate distribution values from an ARIMA in order to plot as a fanchart
我有以下时间序列的 CPI 数据，我正在寻找创建扇形图（类似于 https://journal.r-project.org/archive/2015-1/abel.pdf 或 ggplot2 中的英格兰银行示例，如果可能的话） .
到目前为止，我已经根据我的时间序列创建了一个 ARIMA 模型。我正在寻找有关如何从我的模型中模拟随机变量分布并将其绘制为扇形图的解决方案。我希望为分发提前模拟 10 个周期。
这是我的数据集的可复制品 cpi
structure(list(Date = structure(c(1356998400, 1359676800, 1362096000,
1364774400, 1367366400, 1370044800, 1372636800, 1375315200, 1377993600,
1380585600, 1383264000, 1385856000, 1388534400, 1391212800, 1393632000,
1396310400, 1398902400, 1401580800, 1404172800, 1406851200, 1409529600,
1412121600, 1414800000, 1417392000, 1420070400, 1422748800, 1425168000,
1427846400, 1430438400, 1433116800, 1435708800, 1438387200, 1441065600,
1443657600, 1446336000, 1448928000, 1451606400, 1454284800, 1456790400,
1459468800, 1462060800, 1464739200, 1467331200, 1470009600, 1472688000,
1475280000, 1477958400, 1480550400, 1483228800, 1485907200, 1488326400,
1491004800, 1493596800, 1496275200, 1498867200, 1501545600, 1504224000,
1506816000, 1509494400, 1512086400, 1514764800, 1517443200, 1519862400,
1522540800, 1525132800, 1527811200, 1530403200, 1533081600, 1535760000,
1538352000, 1541030400, 1543622400, 1546300800, 1548979200, 1551398400,
1554076800, 1556668800, 1559347200, 1561939200, 1564617600, 1567296000,
1569888000, 1572566400, 1575158400, 1577836800, 1580515200, 1583020800,
1585699200, 1588291200, 1590969600, 1593561600), class = c("POSIXct",
"POSIXt"), tzone = "UTC"), CPI = c(100.943613610327, 101.355726290109,
101.920519704091, 102.251765014058, 102.399483334481, 102.654230611209,
103.366370423635, 103.771996583604, 104.069828647932, 104.475897454947,
104.745585890252, 104.9, 105.877675706645, 106.600613244374,
107.25658797107, 108.285287342243, 108.607710827378, 108.935592526775,
109.11670321665, 109.390661099815, 109.563232156331, 109.694215435852,
109.939646273932, 109.754097918499, 110.601049654351, 110.415206179718,
110.905507883552, 111.45837834832, 111.873469766967, 112.253828314821,
112.699336213665, 113.056054221625, 113.204653466884, 113.387164759728,
113.581282843726, 113.810860009533, 116.506784014018, 117.199721025597,
118.107968739773, 118.823678758349, 119.420709143437, 119.808600479962,
120.575551335206, 120.774779709305, 121.014544917053, 121.61732414169,
121.917354377998, 122.116542025261, 126.058371342546, 126.285551233707,
126.43426615261, 126.763103151148, 126.92061331762, 127.095652703716,
127.146439944094, 127.257270861715, 127.754395868046, 127.897364611267,
128.227889139291, 128.426778898969, 130.540032633942, 130.730222134177,
130.87769195147, 131.302356289165, 131.797387843531, 132.126557217198,
132.823218725753, 132.868685232286, 133.870800057958, 134.439906096246,
135.351580975176, 135.040382301698, 136.620612224767, 136.503608878263,
136.763944144826, 137.24925661824, 137.169191683167, 137.331600194512,
137.656945057261, 137.792027588476, 137.792027588476, 138.493686354623,
138.681976535356, 138.535078801086, 139.421769773802, 139.848223614133,
139.983926150073, 139.504431667605, 139.994961370897, 140.280481556844,
140.529583177439)), row.names = c(NA, -91L), class = c("tbl_df",
"tbl", "data.frame"))
到目前为止，这是我的模型的代码
# Load Packages
library(pacman)
pacman::p_load(tseries, tidyverse, urca, forecast, tbl2xts)
# Create a log transformation for CPI and convert from tibble to time series format
cpi.ts <- cpi %>%
mutate(CPI = log(CPI)) %>%
tbl_xts()
# Test for a unit root using an ADF test
adf.cpi.ts <- ur.df(cpi.ts, type = "none", selectlags = "AIC")
summary(adf.cpi.ts)
# Create an ARIMA Model using cpi.ts
arima <- auto.arima(cpi.ts)
这是 arima
的结果
ARIMA(0,1,0) with drift
Coefficients:
drift
0.0037
s.e. 0.0005
sigma^2 estimated as 2.255e-05: log likelihood=354.77
AIC=-705.54 AICc=-705.4 BIC=-700.54
我可以使用 arima.sim 函数来做这件事吗（如果可以，我该怎么做呢？）。理想情况下，我正在寻找类似于下图的最终解决方案（不过，如果我能找到 ggplot2 解决方案就更好了。
TIA
这里有两个问题——如何从模型中模拟未来值，以及如何将预测（或模拟）绘制成扇形图。两者都可以使用寓言包来完成。
library(tidyverse)
library(tsibble)
library(fable)
# Create tsibble object
cpi <- cpi %>%
mutate(Date = yearmonth(Date)) %>%
as_tsibble(index=Date)
# Fit ARIMA model to log data
fit <- cpi %>%
model(arima = ARIMA(log(CPI)))
# Simulated future sample paths
fit %>%
generate(times=20, h="1 year") %>%
autoplot(.sim) + autolayer(cpi, CPI) +
ylab("CPI") +
theme(legend.position="none")
# Fan plot
fit %>%
forecast(h="1 year") %>%
autoplot(cpi, level=seq(10,90,by=10)) +
theme(legend.position="none")
由 reprex package (v0.3.0)
于 2020-08-19 创建
我有以下时间序列的 CPI 数据，我正在寻找创建扇形图（类似于 https://journal.r-project.org/archive/2015-1/abel.pdf 或 ggplot2 中的英格兰银行示例，如果可能的话） .
到目前为止，我已经根据我的时间序列创建了一个 ARIMA 模型。我正在寻找有关如何从我的模型中模拟随机变量分布并将其绘制为扇形图的解决方案。我希望为分发提前模拟 10 个周期。
这是我的数据集的可复制品 cpi
structure(list(Date = structure(c(1356998400, 1359676800, 1362096000,
1364774400, 1367366400, 1370044800, 1372636800, 1375315200, 1377993600,
1380585600, 1383264000, 1385856000, 1388534400, 1391212800, 1393632000,
1396310400, 1398902400, 1401580800, 1404172800, 1406851200, 1409529600,
1412121600, 1414800000, 1417392000, 1420070400, 1422748800, 1425168000,
1427846400, 1430438400, 1433116800, 1435708800, 1438387200, 1441065600,
1443657600, 1446336000, 1448928000, 1451606400, 1454284800, 1456790400,
1459468800, 1462060800, 1464739200, 1467331200, 1470009600, 1472688000,
1475280000, 1477958400, 1480550400, 1483228800, 1485907200, 1488326400,
1491004800, 1493596800, 1496275200, 1498867200, 1501545600, 1504224000,
1506816000, 1509494400, 1512086400, 1514764800, 1517443200, 1519862400,
1522540800, 1525132800, 1527811200, 1530403200, 1533081600, 1535760000,
1538352000, 1541030400, 1543622400, 1546300800, 1548979200, 1551398400,
1554076800, 1556668800, 1559347200, 1561939200, 1564617600, 1567296000,
1569888000, 1572566400, 1575158400, 1577836800, 1580515200, 1583020800,
1585699200, 1588291200, 1590969600, 1593561600), class = c("POSIXct",
"POSIXt"), tzone = "UTC"), CPI = c(100.943613610327, 101.355726290109,
101.920519704091, 102.251765014058, 102.399483334481, 102.654230611209,
103.366370423635, 103.771996583604, 104.069828647932, 104.475897454947,
104.745585890252, 104.9, 105.877675706645, 106.600613244374,
107.25658797107, 108.285287342243, 108.607710827378, 108.935592526775,
109.11670321665, 109.390661099815, 109.563232156331, 109.694215435852,
109.939646273932, 109.754097918499, 110.601049654351, 110.415206179718,
110.905507883552, 111.45837834832, 111.873469766967, 112.253828314821,
112.699336213665, 113.056054221625, 113.204653466884, 113.387164759728,
113.581282843726, 113.810860009533, 116.506784014018, 117.199721025597,
118.107968739773, 118.823678758349, 119.420709143437, 119.808600479962,
120.575551335206, 120.774779709305, 121.014544917053, 121.61732414169,
121.917354377998, 122.116542025261, 126.058371342546, 126.285551233707,
126.43426615261, 126.763103151148, 126.92061331762, 127.095652703716,
127.146439944094, 127.257270861715, 127.754395868046, 127.897364611267,
128.227889139291, 128.426778898969, 130.540032633942, 130.730222134177,
130.87769195147, 131.302356289165, 131.797387843531, 132.126557217198,
132.823218725753, 132.868685232286, 133.870800057958, 134.439906096246,
135.351580975176, 135.040382301698, 136.620612224767, 136.503608878263,
136.763944144826, 137.24925661824, 137.169191683167, 137.331600194512,
137.656945057261, 137.792027588476, 137.792027588476, 138.493686354623,
138.681976535356, 138.535078801086, 139.421769773802, 139.848223614133,
139.983926150073, 139.504431667605, 139.994961370897, 140.280481556844,
140.529583177439)), row.names = c(NA, -91L), class = c("tbl_df",
"tbl", "data.frame"))
到目前为止，这是我的模型的代码
# Load Packages
library(pacman)
pacman::p_load(tseries, tidyverse, urca, forecast, tbl2xts)
# Create a log transformation for CPI and convert from tibble to time series format
cpi.ts <- cpi %>%
mutate(CPI = log(CPI)) %>%
tbl_xts()
# Test for a unit root using an ADF test
adf.cpi.ts <- ur.df(cpi.ts, type = "none", selectlags = "AIC")
summary(adf.cpi.ts)
# Create an ARIMA Model using cpi.ts
arima <- auto.arima(cpi.ts)
这是 arima
ARIMA(0,1,0) with drift
Coefficients:
drift
0.0037
s.e. 0.0005
sigma^2 estimated as 2.255e-05: log likelihood=354.77
AIC=-705.54 AICc=-705.4 BIC=-700.54
我可以使用 arima.sim 函数来做这件事吗（如果可以，我该怎么做呢？）。理想情况下，我正在寻找类似于下图的最终解决方案（不过，如果我能找到 ggplot2 解决方案就更好了。
TIA
这里有两个问题——如何从模型中模拟未来值，以及如何将预测（或模拟）绘制成扇形图。两者都可以使用寓言包来完成。
library(tidyverse)
library(tsibble)
library(fable)
# Create tsibble object
cpi <- cpi %>%
mutate(Date = yearmonth(Date)) %>%
as_tsibble(index=Date)
# Fit ARIMA model to log data
fit <- cpi %>%
model(arima = ARIMA(log(CPI)))
# Simulated future sample paths
fit %>%
generate(times=20, h="1 year") %>%
autoplot(.sim) + autolayer(cpi, CPI) +
ylab("CPI") +
theme(legend.position="none")
# Fan plot
fit %>%
forecast(h="1 year") %>%
autoplot(cpi, level=seq(10,90,by=10)) +
theme(legend.position="none")
由 reprex package (v0.3.0)
于 2020-08-19 创建