时间序列分析详解
时间序列是按时间顺序排列的、随时间变化且相互关联的数据序列。
分析时间序 列的方法构成数据分析的一个重要领域,即时间序列分析。 时间序列根据所研究的依据不同,可有不同的分类。
- 按所研究的对象的多少分,有一元时间序列和多元时间序列。
- 按时间的连续性可将时间序列分为离散时间序列和连续时间序列两种。
- 按序列的统计特性分,有平稳时间序列和非平稳时间序列。如果一个时间序列 的概率分布与时间t无关,则称该序列为严格的(狭义的)平稳时间序列。
如果序列的 一、二阶矩存在,而且对任意时刻t满足:
(1)均值为常数
(2)协方差为时间间隔τ的函数。
则称该序列为宽平稳时间序列,也叫广义平稳时间序列。我们以后所研究的时间序列主 要是宽平稳时间序列。
4. 按时间序列的分布规律来分,有高斯型时间序列和非高斯型时间序列。
时间序列方法概述
使用时间序列进行预测
简单平均法
简单平均法就是预测的值为之前过去所有值的平均.当然这不会很准确,但这种预测方法在某些情况下效果是最好的。
移动平均法
我们经常会遇到这种数据集,比如价格或销售额某段时间大幅上升或下降。如果我们这时用之前的简单平均法,就得使用所有先前数据的平均值,但在这里使用之前的所有数据是说不通的,因为用开始阶段的价格值会大幅影响接下来日期的预测值。因此,我们只取最近几个时期的价格平均值。很明显这里的逻辑是只有最近的值最要紧。这种用某些窗口期计算平均值的预测方法就叫移动平均法。
指数平滑法
在做时序预测时,一个显然的思路是:认为离着预测点越近的点,作用越大。比如我这个月体重100斤,去年某个月120斤,显然对于预测下个月体重而言,这个月的数据影响力更大些。假设随着时间变化权重以指数方式下降——最近为0.8,然后0.82,0.83…,最终年代久远的数据权重将接近于0。将权重按照指数级进行衰减,这就是指数平滑法的基本思想。
指数平滑法有几种不同形式:一次指数平滑法针对没有趋势和季节性的序列,二次指数平滑法针对有趋势但没有季节性的序列,三次指数平滑法针对有趋势也有季节性的序列。“
所有的指数平滑法都要更新上一时间步长的计算结果,并使用当前时间步长的数据中包含的新信息。它们通过”混合“新信息和旧信息来实现,而相关的新旧信息的权重由一个可调整的参数来控制。
ARMA模型
AR(Auto Regressive Model)自回归模型是线性时间序列分析模型中最简单的模型。通过自身前面部分的数据与后面部分的数据之间的相关关系(自相关)来建立回归方程,从而可以进行预测或者分析
MA(Moving Average Model)移动平均模型通过将一段时间序列中白噪声(误差)进行加权和,可以得到移动平均方程。如下模型为q阶移动平均过程,表示为MA(q)。
ARMA(Auto Regressive and Moving Average Model)自回归移动平均模型是与自回归和移动平均模型两部分组成。所以可以表示为ARMA(p, q)。p是自回归阶数,q是移动平均阶数。
ARIMA模型
由于序列可能并不平稳但可能在差分后平稳,所以先差分再ARMA就变成了ARIMA
ARIMA(Auto Regressive Integrate Moving Average Model)差分自回归移动平均模型是在ARMA模型的基础上进行改造的,ARMA模型是针对t期值进行建模的,而ARIMA是针对t期与t-d期之间差值进行建模,我们把这种不同期之间做差称为差分,这里的d是几就是几阶差分。ARIMA模型也是基于平稳的时间序列的或者差分化后是稳定的,另外前面的几种模型都可以看作ARIMA的某种特殊形式。表示为ARIMA(p, d, q)。p为自回归阶数,q为移动平均阶数,d为时间成为平稳时所做的差分次数,也就是Integrate单词的在这里的意思。
代码:
clc,clear;
close all;t = 1:240;
y=[2006 1138 993 1027 1155 1175 1299 1786 1781 1570 1387 1635 2795 1679 1329 1585 1406 1260 1814 1556 2143 1963 1555 1217 1738 1261 1141 1099 1274 1325 1197 1596 1923 1789 1606 2174 3089 1746 1571 1545 1593 1480 1925 1290 2637 1907 1922 1660 2355 1104 1093 1085 1125 1257 1055 1547 1927 1881 1490 2156 3147 1706 1599 1544 1738 1747 2380 2049 2713 2543 2139 2105 2378 1461 1295 1313 1360 1582 1425 1656 1822 1602 1404 1837 2589 1052 956 895 972 1124 1620 1280 1492 1587 1330 1219 2064 1347 1325 1273 1461 1481 1511 2000 2137 1643 1464 1706 2236 1065 994 907 948 997 1788 1394 1787 1435 1225 1221 1900 1259 1357 509 557 1375 1548 1802 1995 1720 1468 1625 2515 1371 1439 1642 1344 1247 1737 1434 1782 1595 1348 1083 1864 1367 1299 1301 1308 1382 1401 1837 2045 1790 1409 1179 1852 796 780 602 584 558 831 639 943 960 859 775 976 581 563 547 635 664 684 828 981 889 1303 1490 2246 1335 1270 1281 1210 1078 1664 1337 1632 1454 1333 1146 1721 1059 963 1000 1107 1132 1086 1399 1575 1507 1469 1650 2305 1328 1323 1314 1317 1443 1963 2201 2577 1900 1606 1390 1859 1009 1001 940 981 1239 1023 1284 1510 1497 1448 1677 2188 1499 1313 1351 1237 1454 2031 1615 2230 1899 1574 1666];
figure
plot( t, y )adftest(y)
%% ACF和PACF
figure
subplot(211),autocorr( y );
subplot(212),parcorr( y );
figure
dy = diff( y );
subplot(211),autocorr( dy );
subplot(212),parcorr( dy );
Y=[2006 1138 993 1027 1155 1175 1299 1786 1781 1570 1387 1635 2795 1679 1329 1585 1406 1260 1814 1556 2143 1963 1555 1217 1738 1261 1141 1099 1274 1325 1197 1596 1923 1789 1606 2174 3089 1746 1571 1545 1593 1480 1925 1290 2637 1907 1922 1660 2355 1104 1093 1085 1125 1257 1055 1547 1927 1881 1490 2156 3147 1706 1599 1544 1738 1747 2380 2049 2713 2543 2139 2105 2378 1461 1295 1313 1360 1582 1425 1656 1822 1602 1404 1837 2589 1052 956 895 972 1124 1620 1280 1492 1587 1330 1219 2064 1347 1325 1273 1461 1481 1511 2000 2137 1643 1464 1706 2236 1065 994 907 948 997 1788 1394 1787 1435 1225 1221 1900 1259 1357 509 557 1375 1548 1802 1995 1720 1468 1625 2515 1371 1439 1642 1344 1247 1737 1434 1782 1595 1348 1083 1864 1367 1299 1301 1308 1382 1401 1837 2045 1790 1409 1179 1852 796 780 602 584 558 831 639 943 960 859 775 976 581 563 547 635 664 684 828 981 889 1303 1490 2246 1335 1270 1281 1210 1078 1664 1337 1632 1454 1333 1146 1721 1059 963 1000 1107 1132 1086 1399 1575 1507 1469 1650 2305 1328 1323 1314 1317 1443 1963 2201 2577 1900 1606 1390 1859 1009 1001 940 981 1239 1023 1284 1510 1497 1448 1677 2188 1499 1313 1351 1237 1454 2031 1615 2230 1899 1574 1666];
Y=Y(1:168)';
x=1:168;
% 通过AIC,BIC等准则暴力选定阶数
%[AR_Order,MA_Order] = ARMA_Order_Select(Y,max_ar,max_ma,0);
AR_Order=6;
MA_Order=3;
Mdl = arima(6,0,3);
EstMdl = estimate(Mdl,Y);
[res,~,logL] = infer(EstMdl,Y); %res即残差stdr = res/sqrt(EstMdl.Variance);
figure('Name','残差检验')
subplot(2,3,1)
plot(stdr)
title('Standardized Residuals')
subplot(2,3,2)
histogram(stdr,10)
title('Standardized Residuals')
subplot(2,3,3)
autocorr(stdr)
subplot(2,3,4)
parcorr(stdr)
subplot(2,3,5)
qqplot(stdr)
% 5.预测
step = 24; %预测步数为24
[forData,YMSE] = forecast(EstMdl,step,'Y',Y); %matlab2018及以下版本写为Predict_Y(i+1) = forecast(EstMdl,1,'Y0',Y(1:i));
lower = forData - 1.96*sqrt(YMSE); %95置信区间下限
upper = forData + 1.96*sqrt(YMSE); %95置信区间上限
subplot(2,3,6);
plot(Y,'Color',[.7,.7,.7]);hold on;
h1 = plot(length(Y):length(Y)+step,[Y(end);lower],'r:','LineWidth',2);subplot(2,3,6);
plot(length(Y):length(Y)+step,[Y(end);upper],'r:','LineWidth',2),hold on;
h2 = plot(length(Y):length(Y)+step,[Y(end);forData],'k','LineWidth',2);subplot(2,3,6);
legend([h1 h2],'95% 置信区间','预测值','Location','NorthWest')
title('Forecast')
拓展:
SARIMA模型
SARIMA(Seasonal Autoregressive Integrated Moving Average)模型在非稳 的数据上多了一些对有周期性特征数据的处理。其主要原理是来源于外界对 数据周期性的观察,得到季节的长度(s)、季节自回归的阶数(P,取值由 PACF判断)、季节移动平均的阶数(Q,取值由ACF判断)季节差分的阶数(D, 一般就是0或1),再建立关于季节性的模型,将该模型和ARIMA模型结合,就 是SARIMA模型。
SARIMA季节性自回归移动平均模型模型在ARIMA模型的基础上添加了季节性的影响,结构参数有七个:SARIMA(p,d,q)(P,D,Q,s)
其中p,d,q分别为之前ARIMA模型中我们所说的p:趋势的自回归阶数。d:趋势差分阶数。q:趋势的移动平均阶数。
P:季节性自回归阶数。
D:季节性差分阶数。
Q:季节性移动平均阶数。
s:单个季节性周期的时间步长数。
AIC-BIC准则
实际上我们仍然会得到针对训练数据的多个较为合理的模型,比如混合高斯模型中不同的协方差矩阵和不同 的分量数会产生不同的结果,如何在这些模型中选择最优的模型呢? 模型选择问题就是需要在模型复杂性和模型对数据集描述能力(即似然函数)之间寻找一个平衡,而常用的 模型选择方法就是BIC和AIC,模型选择中AIC和BIC值均是越小越好。在BIC和AIC描述过程中,当训练数据 足够多时,可以不断提高模型精度,即是似然函数会越大,模型对数据集的描述能力越强,而复杂性方面往 往体现在参数数量,即参数越少越好,通过下式就可以理解AIC和BIC公式的原理,实际运用中BIC的运用较 多。