R语言 时间序列有关各种函数总结

一、基础篇

1、安装程序包

es('程序包名字')

例：es('TSA') //安装TSA程序包

注：直接输入es() 可弹出程序包窗口，若未选择镜像，会先弹出镜像窗口，选择镜像后再弹出程序包窗口

2、加载程序包

library(程序包名字)

例：library(TSA) 加载TSA程序包

注：程序包名字前后不需要打引号

3、加载程序包内数据

data(数据名)

例：data(larain) //加载larain数据

注：加载之后才可以使用此数据，否则会提示数据不存在

4、构建时间序列

时间序列名字=ts(数据名,freq=单位时间内的数据个数,start=第一个数据对应的时间)

例：a=ts(larain,freq=12,start=2000) //将larain转化成从2000年1月开始的月度数据

注：freq=4表示季度数据，12表示月度数据；start可以等于一个数，表示年份，即数据是从该年第一月开始的，也可以等于c(年份,月份)，表示从某年某月开始

5、从TXT文档读取单变量数据

数据名=scan('文件名字.txt')

例：a=scan('') //读取file1中的变量，并命名为a

注：读取的文档需放到工作目录之下

6、从TXT文档读取多变量数据

数据名=('文件名字.txt',header=T/F)

attach(数据名字)

例：a= ('', header=F) //将file2中数据命名为a，系统自带列标V1,V2……

attach(a) //打开a，之后可以直接调用a中各变量

V1 //输出a中第一列数据

注：读取的文档需放到工作目录之下，header=T表示文档中第一行是各列数据的名字，第二行开始是数据，header=F表示系统会自动给各列数据命名，调用数据之前必须先打开它

7、从CSV格式的表格文件中读取数据

数据名=('文件名字.csv',header=T/F)

attach(数据名)

例：a=('', header=T) //将file3中数据命名为a，第一行为各列数据的名字

attach(a) //打开a，之后可以直接调用a中各变量

X //输出a中第一列数据（假如表中第一列第一行为X）

注：读取的文件需放到工作目录之下，header参数含义同上，调用数据之前也必须先打开

8、导出数据

(需要导出的数据名,file ="生成的文件名.txt",sep="数据分隔符号",

=T/F, =T/F, quote =T/F)

例： (larain, file ="", =F, =F, quote =F)

//将larain导出到file4文档中，没有分隔符，行标不导出，列标也不导出，数据不加引号

注：生成的文件在工作目录下，sep参数在多列数据导出时，可以用到来间隔各列数据，如

果是单列数据，可以省略。表示行标（数据的序号）是否导出，表示列标（变量名）是否导出。当数据不是数值型时，有时可能需要加引号，此时可令quote =T。

二、确定趋势篇

1、画出数据的时间序列图

plot(数据名,type='数据显示方式代码')

例：plot(hours,type='o') //绘制hours的时间序列图，显示方式为点线相连

注：type常见的有o为点线相连、p为只有点、l为只有线，h为垂直线表示每个数据

2、拟合线性趋势模型

模型名字=lm(数据名~time(数据名))

3、拟合二次时间趋势模型

模型名字=lm(数据名~time(数据名)+I(time(数据名)^2))

4、拟合季节均值模型

模型名字=lm(数据名~season(数据名)-1)

5、拟合余弦趋势模型

模型名字=lm(数据名~harmonic(数据名,1))

6、拟合季节均值+二次时间趋势模型

模型名字=lm(数据名~ season(数据名)+time(数据名)+ I(time(数据名)^2))

注：输出结果中第一行为1月份系数，2~12月的系数均需要加上1月份系数才是各月实际系数，但是时间t是二次时间t^2的系数不需要加。

7、输出拟合结果

summary(模型名字)

8、保存模型的残差序列

残差序列名=rstudent(模型名字) //数据类型为时间序列

注：这里保存的是学生化残差，R中还有其它残差，residual()是原始残差（也可以用resid()），rstandard()是标准残差，但并不是任何时候都可以调用成功。

9、对数据进行游程检验（检验独立性）

runs(数据名)

注：游程检验的原假设为检验对象是独立的，因此一般来说，当p值小于0.05时，拒绝原假设，认为检验对象不是独立的，当p值大于0.05时，不能拒接原假设，认为检验对象是独立的。

9、对数据进行SW正态性检验

(数据名)

注：SW正态性检验的原假设为检验对象具有正态性，因此一般来说，当p值小于0.05时，拒绝原假设，认为检验对象不具有正态性，当p值大于0.05时，不能拒接原假设，认为检验对象具有正态性。

10、绘制样本自相关图判断独立性

acf(数据名)

注：一般来说，当各滞后上的样本自相关系数大致在两条虚线范围之内时，认为数据的自相关性很弱，即可认为数据是独立的，反之不独立

11、绘制直方图判断正态性

hist(数据名)

注：直方图呈现中间大两端小的情况时，可认为数据可能服从正态分布，但还需进一步验证

12、绘制正态图（qq图）判断正态性

qqnorm (数据名)

qqline(数据名) //在正态图上加一条穿过第一分位数和第三分位数的直线帮助判断

注：当正态图近似呈一条直线时，可认为数据可能服从正态分布；当图形有明显偏离时，认为数据可能不服从正态分布，但还需进一步验证

13、绘制模型残差和拟合值的对比图

拟合序列名=fitted(模型名字) //从模型中读取拟合值，数据类型为时间序列

残差序列名=rstudent(模型名字)

plot(x=(拟合序列名),y=残差序列名) //将拟合值转化成向量数据

三、ARIMA模型篇

1、模拟一个ARIMA序列

数据名=(model=list(ma=-c(1,,q ),ar=c(1,,p),n=100)+均值

2、对时间序列差分

diff(数据名,lag=滞后,diff=差分次数)

例：diff(data1,lag=2,diff=1) //表示滞后为2的一次差分

3、求幂变换的最优值

(数据名, lambda = seq(最小值,最大值,间隔))

注：lambda参数可以用来放大95%置信区间，可以更加清楚地确定最优值

4、求ARMA模型的理论自相关系数

ARMAacf(ma=-c(1,,q),ar=c (1,,p),=最大滞后)

5、绘制样本自相关（ACF）图

acf(数据名,plot=T/F,xaxp=c(最小值,最大值,刻度个数))

注：plot=F会显示每个滞后下的样本自相关系数，但不显示图形；xaxp是横坐标刻度参数；若滞后q阶后呈现截尾现象，则合适的模型可能为MA(q)

6、绘制样本偏自相关（PACF）图

pacf(数据名)

注：若滞后p阶后呈现截尾现象，则合适的模型可能为AR(p)

7、计算样本扩展自相关函数（EACF）

eacf(数据名)

注：找到一个角∠，其内部尽可能地都是O，可由该角的顶点所处位置找到合适的模型，其横坐标为MA阶数，纵坐标为AR阶数

8、ADF单位根检验判断平稳性

library(tseries) //需要安装并加载tseries程序包

(数据名)

注：ADF单位根检验的原假设为检验对象为非平稳的，因此一般来说，当p值小于0.05时，拒绝原假设，认为检验对象平稳，当p值大于0.05时，不能拒接原假设，认为检验对象非平稳。

9、最优子集ARMA法

结果=armasubsets(y=数据名,nar=AR最大阶,nma=MA最大阶)

plot(结果)

10、对适合AR模型的时间序列进行参数估计

ar(数据名,=最大AR阶数,AIC=F/T,method='估计方法')

注：AIC=F，则模型的阶数为最大AR阶数，AIC=T，则在0-最大AR阶数之间，选择最小AIC来确定模型的阶数；估计方法主要有三种，yw表示Yule-Walker方程，也就是矩估计法，ols表示最小二乘法，mle表示极大似然估计法

11、对适合ARIMA模型的时间序列进行参数估计

模型名=airma(数据名,order=c(p,d,q),method='估计方法',fixed=c(NA/0))

模型名 //显示结果

例：model1=arima(a,order=c(3,0,0),fixed=c(NA,0,NA,NA)) //设定AR(2)的系数为0

注：估计方法主要有两种：CSS表示条件平方和法，即最小二乘法，ML表示极大似然估计，默认是ML；fixed参数可以在一些系数不显著时，对它们设定为0再进行参数估计，顺序是AR、MA、截距项；截距项的估计值为均值估计值，不是0

12、获得模型的残差值

residuals(模型名) 或者resid(模型名) //得到原始残差

rstandard(模型名) //得到标准残差

注：学生化残差在这里不适用了

13、LB检验判断独立性

(模型名，lag=最大滞后)

注：LB检验的原假设为检验对象为不相关的，因此一般来说，当p值小于0.05时，拒绝原假设，认为检验对象是相关的，即不独立，当p值大于0.05时，不能拒接原假设，认为检验对象是不相关的，可看作是独立的

14、综合诊断（标准残差图、样本自相关、LB检验）

tsdiag(模型名,gof=最大滞后,l=F/T)

注：l=T，会忽略数据最初的几个滞后，适用于有周期性或季节性的数据

15、得到时间序列的子序列

子序列名=window(原序列名,start=起点,end=终点)

注：起点默认为1，重点默认为原序列长度

16、利用估计模型进行预测

模型名=arima(……) //先进行参数估计,内置参数见11

预测结果名=plot(模型名,=预测个数) //预测并画图

预测序列=预测结果名$pred //将预测值保存在预测序列中

预测序列 //显示预测值

cbind(真实序列,预测序列) //可对比真实值与预测值