[时间序列分析][2]--趋势和(季节)因子
一句话概括今天要讲的内容:
对于有趋势的时间序列的数据,我们可以先把趋势给拟合出来,然后在对残差进行分析,看残差是否是稳定的,是否是白噪声序列。
我们对于趋势和因子举两个例子,分别来看一下:
一.趋势
下面这组数据是 “上海证券交易所1991年1月 - 2001 年10月每月末上证指数序列”
{130.44, 133.47, 120.19, 113.94, 114.83, 137.56, 143.8, 178.43, \180.92, 218.6, 259.6, 292.75, 313.24, 364.66, 381.24, 445.38, \1234.71, 1191.19, 1052.07, 823.27, 702.32, 507.25, 724.6, 780.39, \1198.48, 1339.88, 925.91, 1358.78, 935.48, 1007.05, 881.07, 895.68, \890.27, 814.82, 984.93, 833.8, 770.25, 770.98, 704.46, 592.56, \556.26, 469.29, 333.92, 785.33, 791.15, 654.98, 683.59, 647.87, \562.59, 549.26, 646.92, 579.93, 700.51, 630.58, 695.55, 723.87, \722.43, 717.32, 641.13, 555.29, 537.34, 552.93, 556.39, 681.16, \643.65, 804.25, 822.48, 809.94, 875.52, 976.71, 1032.95, 917.02, \964.74, 1040.27, 1234.62, 1393.75, 1285.18, 1250.27, 1189.76, \1221.06, 1097.38, 1180.39, 1139.63, 1194.1, 1222.91, 1206.53, \1243.01, 1343.44, 1411.2, 1339.2, 1316.91, 1150.22, 1242.09, 1217.31, \1247.42, 1146.7, 1134.67, 1090.09, 1158.05, 1120.92, 1279.32, \1689.42, 1601.45, 1627.12, 1570.7, 1504.56, 1434.97, 1366.58, 1535., \1714.58, 1800.22, 1836.32, 1894.55, 1928.11, 2023.54, 2021.2, \1910.16, 1961.29, 2070.61, 2073.48, 2065.61, 1959.18, 2112.78, \2119.18, 2213.18, 2218.03, 1961.29, 2070.61, 2073.48, 2065.61}看一下散点图
从散点图上我们可以看出这个形状类似于 a^x^2 + b*x + c 这个的形式,我们下面采用非线性拟合,来看一下
sol = FindFit[data, a*x^2 + b*x + c, {a, b, c}, x]p1 = Plot[a^x^2 + b*x + c /. sol, {x, 1, 130}, PlotStyle -> {Red}];把大概的趋势拟合出来了。
我们来做一下平稳性的检验,分别对原始数据和拟合后数据的残差(之后我们就是对于残差做分析了)
原始数据的平稳性
drawdata = data;ListPlot[CorrelationFunction[drawdata, {50}], Filling -> Axis, FillingStyle -> Directive[Thickness[.005], Green, Dashed], PlotRange -> All]可见不是平稳的,其实看原始数据的散点图我们就能知道他是不平稳的了。
拟合后的残差
先画出残差的图
rad = data - (a*#^2 + b*# + c /. sol & /@ Range[Length[data]]);ListPlot[rad]
其实可以看到残差还是具有一定的规律的,我们来检验一下平稳性
ListPlot[CorrelationFunction[rad, {50}], Filling -> Axis, FillingStyle -> Directive[Thickness[.005], Green, Dashed], PlotRange -> All]
(这个图好像还是有点问题的)
二.季节因子
接下来我们看一个有季节因子的例子
下面数据来源:1993-2000年中国社会消费品零售总额序列
我们看一下数据的样子
最上面一行是年份,左边一列是月标
我直接把处理好的数据贴在下面,方便直接使用
{ {977.5, 892.5, 942.3, 941.3, 962.2, 1005.7, 963.8, 959.8, 1023.3, 1051.1, 1102., 1415.5}, {1192.2, 1162.7, 1167.5, 1170.4, 1213.7, 1281.1, 1251.5, 1286., 1396.2, 1444.1, 1553.8, 1932.2}, {1602.2, 1491.5, 1533.3, 1548.7, 1585.4, 1639.7, 1623.6, 1637.1, 1756., 1818., 1935.2, 2389.5}, {1909.1, 1911.2, 1860.1, 1854.8, 1898.3, 1966., 1888.7, 1916.4, 2083.5, 2148.3, 2290.1, 2848.6}, {2288.5, 2213.5, 2130.9, 2100.5, 2108.2, 2164.7, 2102.5, 2104.4, 2239.6, 2348., 2454.9, 2881.7}, {2549.5, 2306.4, 2279.7, 2252.7, 2265.2, 2326., 2286.1, 2314.6, 2443.1, 2536., 2652.2, 3131.4}, {2662.1, 2538.4, 2403.1, 2356.8, 2364., 2428.8, 2380.3, 2410.9, 2604.3, 2743.9, 2781.5, 3405.7}, {2774.7, 2805., 2627., 2572., 2637., 2645., 2597., 2636., 2854., 3029., 3108., 3680.}}画出散点图
可以看到很明显的趋势,和周期性质。
1.首先计算出月因子-- 拿当月的平均/总的平均
md = Mean[Flatten@data];mouthdata = Table[Mean[data[[All, i]]]/md, {i, 1, 12}];Grid[Table[{ToString[i] <> "月", mouthdata[[i]]}, {i, 1, 12}], Frame -> All]
画出月因子指数图
从图中可以看到12月份比其他月份都要高,和一些其他的什么什么之类的,大家自己看看吧,反正意思大家都懂的。
2.
消除月因子的办法就是拿每个月的数据去除他对应的月因子,看下面代码p1 = ListPlot[Flatten[data]/Flatten[Table[mouthdata, {i, 1, 8}]]]画出来的图如下所示
3.拟合消除月因子的直线
把月因子消除掉之后,我们就要来拟合总体的趋势了
lm = LinearModelFit[ Flatten[data]/Flatten[Table[mouthdata, {i, 1, 8}]], x, x];p2 = Plot[Normal[lm], {x, 0, 100}, PlotStyle -> {Red, Dashed}];Show[p1, p2]
然后求一下残差看一下(其实好像没什么用,即时残差不好我也不想改模型了....当然,如果想要认真做的话就要分析一下了)
4.
最后一步了.....为什么感觉写了好久....
p1 = ListLinePlot[ Table[lm[i]*mouthdata[[If[Mod[i, 12] == 0, 12, Mod[i, 12]]]], {i, 1, 96}], PlotStyle -> {Red}, PlotRange -> All];Show[p, p1, PlotRange -> All]
到这里就把预测都做好了,红色的为预测的线,可以看到和蓝色的点重合还是可以的,我们来看一下残差...我都不想看了...但是我们还是要看一下
ListPlot[Flatten[data] - Table[lm[i]*mouthdata[[If[Mod[i, 12] == 0, 12, Mod[i, 12]]]], {i, 1, 96}]]
恩,形状是和上面的是一样的。
那么,到这里这一节就结束了。今天就主要讲了对于有趋势的数据应该如何来处理
希望大家多多留言,给我提一下建议,谢谢!
2017/3/26