matlab估计arma残差,pythonARIMA时间序列
1 时间序列与时间序列分析
在⽣产和科学研究中，对某⼀个或者⼀组变量 x(t) 进⾏观察测量，将在⼀系列时刻 t1,t2,…,tn 所得到的离散数字组成的序列集合，称之为时间序列。
时间序列分析是根据系统观察得到的时间序列数据，通过曲线拟合和参数估计来建⽴数学模型的理论和⽅法。时间序列分析常⽤于国民宏观经济控制、市场潜⼒预测、⽓象预测、农作物害⾍灾害预报等各个⽅⾯。
2 时间序列建模基本步骤
获取被观测系统时间序列数据；
对数据绘图，观测是否为平稳时间序列；对于⾮平稳时间序列要先进⾏d阶差分运算，化为平稳时间序列；
经过第⼆步处理，已经得到平稳时间序列。要对平稳时间序列分别求得其⾃相关系数ACF 和偏⾃相关系数PACF，通过对⾃相关图和偏⾃相关图的分析，得到最佳的阶层 p和阶数 q
由以上得到的d、q、p ，得到ARIMA模型。然后开始对得到的模型进⾏模型检验。
3 ARIMA实战解剖
原理⼤概清楚，实践却还是会有诸多问题。相⽐较R语⾔，Python在做时间序列分析的资料相对少很多。下⾯就通过Python语⾔详细解析后三个步骤的实现过程。
⽂中使⽤到这些基础库： pandas,numpy,scipy,matplotlib,statsmodels。 对其调⽤如下
from __future__ import print_function
import pandas as pd
import numpy as np竹笋做法
from scipy import stats
import matplotlib.pyplot as plt
import statsmodels.api as sm
aphics.api import qqplot
3.1 获取数据
这⾥我们使⽤⼀个具有周期性的测试数据，进⾏分析。
数据如下：软装设计师
dta=[10930,10318,10595,10972,7706,6756,9092,10551,9722,10913,11151,8186,6422,
6337,11649,11652,10310,12043,7937,6476,9662,9570,9981,9331,9449,6773,6304,9355,
10477,10148,10395,11261,8713,7299,10424,10795,11069,11602,11427,9095,7707,10767,
12136,12812,12006,12528,10329,7818,11719,11683,12603,11495,13670,11337,10232,
13261,13230,15535,16837,19598,14823,11622,19391,18177,19994,14723,15694,13248,
9543,12872,13101,15053,12619,13749,10228,9725,14729,12518,14564,15085,14722,
会议记录模板范文
11999,9390,13481,14795,15845,15271,14686,11054,10395]
dta=np.array(dta,dtype=np.float) //这⾥要转下数据类型，不然运⾏会报错
dta=pd.Series(dta)
dta.index = pd.Index(sm.tsa.datetools.dates_from_range('2001','2090')) //应该是2090，不是2100
dta.plot(figsize=(12,8))
plt.show()
3.2 时间序列的差分d
ARIMA 模型对时间序列的要求是平稳型。因此，当你得到⼀个⾮平稳的时间序列时，⾸先要做的即是做时间序列的差分，直到得到⼀个平稳时间序列。如果你对时间序列做d次差分才能得到⼀个平稳序列，那么可以使⽤ARIMA(p,d,q)模型，其中d是差分次数。
fig = plt.figure(figsize=(12,8))
ax1= fig.add_subplot(111)
diff1 = dta.diff(1)
diff1.plot(ax=ax1)
plt.show()
⼀阶差分的时间序列的均值和⽅差已经基本平稳，不过我们还是可以⽐较⼀下⼆阶差分的效果
fig = plt.figure(figsize=(12,8))
ax2= fig.add_subplot(111)
diff2 = dta.diff(2)
有请diff2.plot(ax=ax2)
plt.show()
可以看出⼆阶差分后的时间序列与⼀阶差分相差不⼤，并且⼆者随着时间推移，时间序列的均值和⽅差保持不变。因此可以将差分次数d设置为1。
其实还有针对平稳的检验，叫“ADF单位根平稳型检验”，以后再更。
3.3 合适的p,q
现在我们已经得到⼀个平稳的时间序列，接来下就是选择合适的ARIMA模型，即ARIMA模型中合适的p,q。
第⼀步我们要先检查平稳时间序列的⾃相关图和偏⾃相关图。
dta= dta.diff(1)#我们已经知道要使⽤⼀阶差分的时间序列，之前判断差分的程序可以注释掉 //原⽂有错误应该是diff1= dta.diff(1)，⽽⾮dta= dta.diff(1)
fig = plt.figure(figsize=(12,8))
ax1=fig.add_subplot(211)
fig = sm.graphics.tsa.plot_acf(dta,lags=40,ax=ax1)
ax2 = fig.add_subplot(212)
fig = sm.graphics.tsa.plot_pacf(dta,lags=40,ax=ax2)
其中lags 表⽰滞后的阶数，以上分别得到acf 图和pacf 图
通过两图观察得到：
音乐枕* ⾃相关图显⽰滞后有三个阶超出了置信边界；
* 偏相关图显⽰在滞后1⾄7阶(lags 1,2,…，7)时的偏⾃相关系数超出了置信边界，从lag 7之后偏⾃相关系数值缩⼩⾄0
则有以下模型可以供选择：
1. ARMA(0,1)模型：即⾃相关图在滞后1阶之后缩⼩为0，且偏⾃相关缩⼩⾄0，则是⼀个阶数q=1的移动平均模型；
2. ARMA(7,0)模型：即偏⾃相关图在滞后7阶之后缩⼩为0，且⾃相关缩⼩⾄0，则是⼀个阶层p=7的⾃回归模型； //原⽂错写为3
3. ARMA(7,1)模型：即使得⾃相关和偏⾃相关都缩⼩⾄零。则是⼀个混合模型。
4. …还可以有其他供选择的模型
现在有以上这么多可供选择的模型，我们通常采⽤ARMA模型的AIC法则。我们知道：增加⾃由参数的数⽬提⾼了拟合的优良性，AIC⿎励数据拟合的优良性但是尽量避免出现过度拟合(Overfitting)的情况。所以优先考虑的模型应是AIC值最⼩的那⼀个。⾚池信息准则的⽅法是寻找可以最好地解释数据但包含最少⾃由参数的模型。不仅仅包括AIC准则，⽬前选择模型常⽤如下准则：
* AIC=-2 ln(L) + 2 k 中⽂名字：⾚池信息量 akaike information criterion
* BIC=-2 ln(L) + ln(n)*k 中⽂名字：贝叶斯信息量 bayesian information criterion放松大脑
* HQ=-2 ln(L) + ln(ln(n))*k hannan-quinn criterion
构造这些统计量所遵循的统计思想是⼀致的，就是在考虑拟合残差的同时，依⾃变量个数施加“惩罚”。但要注意的是，这些准则不能说明某⼀个模型的精确度，也即是说，对于三个模型Ａ，Ｂ，Ｃ，我们能够判断出Ｃ模型是最好的，但不能保证Ｃ模型能够很好地刻画数据，因为有可能三个模型都是糟糕的。
arma_mod20 = sm.tsa.ARMA(dta,(7,0)).fit()
print(arma_mod20.aic,arma_mod20.bic,arma_mod20.hqic)
arma_mod30 = sm.tsa.ARMA(dta,(0,1)).fit()
print(arma_mod30.aic,arma_mod30.bic,arma_mod30.hqic)
arma_mod40 = sm.tsa.ARMA(dta,(7,1)).fit()
print(arma_mod40.aic,arma_mod40.bic,arma_mod40.hqic)
arma_mod50 = sm.tsa.ARMA(dta,(8,0)).fit()
print(arma_mod50.aic,arma_mod50.bic,arma_mod50.hqic)
arma_mod70 = sm.tsa.ARMA(dta,(7,0)).fit()
print(arma_mod70.aic,arma_mod70.bic,arma_mod70.hqic)
arma_mod30 = sm.tsa.ARMA(dta,(0,1)).fit()
print(arma_mod30.aic,arma_mod30.bic,arma_mod30.hqic)
arma_mod71 = sm.tsa.ARMA(dta,(7,1)).fit()
print(arma_mod71.aic,arma_mod71.bic,arma_mod71.hqic)
arma_mod80 = sm.tsa.ARMA(dta,(8,0)).fit()
print(arma_mod80.aic,arma_mod80.bic,arma_mod80.hqic)
把原⽂的mod名字改了，以qp阶数命名更直观，不过我跑出来的结果和原⽂不太⼀样
1619.19185018 1641.69013721 1628.26448199
皂角米的功效
1657.21729729 1664.71672631 1660.2415079
1605.68656094 1630.68465765 1615.76726295
1597.93598102 1622.93407772 1608.01668303
这样的话应该是ARMA(8，0)模型拟合效果最好。
可以看到ARMA(7,0)的aic，bic，hqic均最⼩，因此是最佳模型。
3.4 模型检验
在指数平滑模型下，观察ARIMA模型的残差是否是平均值为0且⽅差为常数的正态分布(服从零均值、⽅差不变的正态分布)，同时也要观察连续残差是否(⾃)相关。
3.4.1 我们对ARMA(7,0)模型所产⽣的残差做⾃相关图
使⽤ARMA(8,0)
接下来检验下残差序列：
resid = sid //原⽂把这个变量赋值语句漏了，所以⽼是出错
fig = plt.figure(figsize=(12,8))
ax1 = fig.add_subplot(211)
fig = sm.graphics.tsa.plot_acf(resid.values.squeeze(), lags=40, ax=ax1)
ax2 = fig.add_subplot(212)一曲红绡不知数
fig = sm.graphics.tsa.plot_pacf(resid, lags=40, ax=ax2)
残差的ACF和PACF图，可以看到序列残差基本为⽩噪声
