基于Rstudio的时间序列预测和计算（第⼀节）
引⾔:
时间序列，⼜称动态数列。顾名思义，是指将同⼀统计指标的数值按其发⽣的时间先后顺序排列⽽成的数列。涉及应⽤的领域⾮常⼴泛，譬如机器学习，统计学，数学建模等等。这篇⽂章将会以R中⾃带的康涅狄格州纽⿊⽂地区从1912年到1971年每年的平均⽓温数据集(nhtemp)为例，进⾏时间序列的预测和计算。
1.数据的转化
在我们进⾏对时间序列的预测和计算之前，我们要清楚：我们⼿⾥的数据是来⾃数据集的，或者是通过不同的统计平台获取的，他们还不是对象，需要转化，R⾥的与时间序列有关的很多库函数，如HoltWinters()函数,ets()函数等对参数的要求都是有严格的标准的，必须是时序对象才能使这些函数正确运⾏，那么转化的代码如下:
nhtemp1<-ts(nhtemp,frequency=12,start=c(1912,2),end=c(1912,12)) %对nhtemp数据集进转化
上⾯这⾏代码就把nhtemp转化成了⼀个时序对象，并把其储存在nhtemp1中，frequency，start，end都是参数，根据⾃⼰的需要进⾏增删:
frequency：顾名思义是频率，周期的意思。直观地来说就是你使⽤时间序列时你⾃⼰定义的后⾯要进⾏差分，统计等等动作的周期。start，开始，我们需要使⽤定义向量的函数c()来定义序列起始位置的⽇期，上⾯这个意思就是⽤建⽴向量的⽅式定义⽬标数据的开头的时间是1912年2⽉。
end，结束，⽤法和start⼀样。值得注意的是，start和end可以⽤于对数据集取⼦集，⼗分好⽤。
2.基于Rstudio的各种操作
作为⼀门⾼级可视化分析的语⾔，⼯具。很明显在我们读取完数据后，很常见的操作就是绘图，这也是R语⾔相较于其他数学建模⼯具的⼀个很鲜明的特点：⾼级可视化统计分析(⼩编的体会)。
2.1对时间序列对象的绘图
我们可以使⽤下⾯这些命令，对前⾯的数据集进⾏绘图(两者在效果上是等效的)：
plot(nhtemp1) %对已经转化的数据，就已经是时间序列对象直接执⾏画图
plot.ts(nhtemp)_ %使⽤plot.ts()对未经转化的数据集进⾏⼀个制图
在R语⾔中，有着⾮常多的plot.___()之类的函数，它们旨在对各种各类的对象进⾏⼀个制图，不可否认，在程序开发的过程中⼗分简便.。效果如下：
2.2对时序对象进⾏简单平滑处理
在该⼩节开始前，我们要清楚：为什么要对时序对象进⾏平滑处理呢？
因为数据集他本⾝是⼀个真实或者仿真的数据集合，⼀般不存在规律性。⽽且数据⼀般是离散型，因此，在时序数据集会有很显著的随机或误差部分，这对我们后续的统计和数值分析的合理性和准确性会产⽣很多不必要的影响，⽽平滑处理就可以解决、避开这些不合理的“数据波动”。
简单移动平均是画出平滑曲线最简易的⽅法，其中⽐较简单易操作的⼀种移动⽅法就是：每个数据点都可使⽤这⼀点和前后两个点的平均值来表⽰，其名⽈：居中移动平均。他的数学表达如下:
其中St就是时间点的平滑值,k=2q+1就是每次⽤来平均的观测值个数(⼤家可以结合数轴来想象)。使⽤该⽅法的代价就是我们将会失去时序集⾥⾯最后的(k-1)/2个观测值。
R中有⼏个函数都能完成这种简单平均移动，包括TTR包中的SMA()函数，zoo包中的rollmean()函数，forecast包中的ma()函数。我们这⾥的代码将使⽤ma()函数对数据集进⾏平滑处理。代码如下:
library(forecast) %载⼊forecast包，⼀般是使⽤library()函数载⼊，若不载⼊，可能⽆法使⽤包⾥专有的函数
opar<-adonly=TRUE) %该语句是⽤来更改当前变量环境， par(opar)⽤来还原默认变量环境。
par(mfrow=c(2,2)) %par()函数就是打开⼀个绘图设备，mfrow=c(2,2)意思就是产⽣4张图，两⾏两列
ylim<-c(min(nhtemp),max(nhtemp)) %以数据集nhtemp的最⼤值和最⼩值为边界绘制后⾯的的四张图的阈值
plot(nhtemp,main="中北⼤学数学建模") %画出第⼀张图
plot(ma(nhtemp,3),main="中北⼤学数学建模，光滑⽔平为3"，ylim=ylim) %画出第⼆张图
plot(ma(nhtemp,5),main="中北⼤学数学建模，光滑⽔平为5",ylim=ylim) %画出第三张图
plot(ma(nhtemp,7),main="中北⼤学数学建模，光滑⽔平为7",ylim=ylim) %画出第四张图
在程序段⾥⾯ma()函数，第⼀个参数是要平滑处理的数据集，第⼆个参数是指定图像的光滑⽔平。从上⾯四张图我们可以看出：随着光滑⽔平K越来越⼤，图像会变的越来越来光滑，因此我们以后在处理平滑处理数据集的问题上找到最能画出数据中规律的K，避免过于平滑或⽋平滑。由于相关科学理论还没有解决怎么公式⾼效化地选取K的值，因此我们只能多多尝试，看出他们的区别。
2.3季节性分解
对于时间间隔⼤于1的时序对象，我们需了解的就不仅仅是总体趋势了，此时，我们需要通过季节性分解帮助我们探究季节性波动以及总体趋势。
⾸先我们要知道的是：存在季节性因素的时间序列数据(如⽉度，季度数据等等)可以被分解为趋势因⼦、季节性因⼦和随机因⼦。趋势因⼦能捕捉到长期变化；季节性因⼦能捕捉到⼀年内的周期性变化;⽽随机误差因⼦则能捕捉到那些不能被趋势或季节效应解释的变化。
接着,⼩编要先向⼤家普及⼀下加法模型和乘法模型：就是将形成时间序列变动的四类构成因素，按照它们的影响⽅式不同，设定的两种组合模型:
乘法模型：Y = T·S·C·I
加法模型：Y = T+S+C+I
式中：Y：时间序列的指标数值 T：长期趋势成分 S：季节变动成分 C：循环变动成分 I：不规则变动成分
季节因素的表述的不同在于：
乘法模型是假定四个因素对现象的发展的影响是相互作⽤的，以长期趋势成分的绝对量为基础，其余量均以⽐率表⽰，加法模型是假定四个因素的影响是相互独⽴的，每个成分均以绝对量表⽰。
那么，我们接下来分解时序数据就要⽤到上⾯提到的两种模型。
对于加法模型，各种因⼦之和应等于对应的时序值，即:
其中时刻t的观测值即这⼀时刻的趋势值、季节效应以及随机影响之和。
⽽乘法模型则将时间序列表⽰为:
即趋势项、季节项和随机影响相乘。
2.3.1加法模型和乘法模型之间的转化(基于R语⾔)
在进⾏这⼀步时我们必须知道，加法模型和乘法模型有什么不同，当然不考虑形式上(幼⼉园同学都知道公式形式不⼀样)。根据⼩编的导师所⾔，两者最⼤的区别可概括为：加法模型与变量是⽆关的，乘法模型与变量是有关的。我在这⾥可以举个例⼦：假设我们有⼀个时序，记录了10年摩托车的⽉销量。在加法模型中，11⽉和12⽉(圣诞节)的销量⼀般会增加500，⽽1⽉(⼀般是销售淡季)的销量则会减少200。此时季节性波动和加法模型是⽆关的！在乘法模型中，11⽉和12⽉的销售量则会增加20%，1⽉的销量减少10%。即季节性的波动量和当时的销售是成⽐例的。这也使得在很多时候，乘法模型⽐加法模型更贴近实际。
好的，在我们知道了两者⽐较明显的区别之后，我们就要知道怎么两者之间怎么通过R语⾔达到转化的⽬的:
从数学思想来说，要进⾏乘法模型和加法模型的转化，我们很容易就想到进⾏对数转化的⽅法(对数函数有相乘到相加的转化性质)，式⼦如下:
通过对数函数的基本性质，我们就可以简单地对两种模型进⾏⼀种转化。
⽽在实际的算法编写中，两者的转化⼯作是这样实现的:
plot(nhtemp) %画出原始的时间序列
Lnhtemp<-log(nhtemp) %做对数化，数据⼤⼩会发⽣变化，但⽤于统计的⽐例关系不变
plot(Lnhtemp,ylab="log(nhtemp)") %画出对数化后的时间序列图像
fit<-stl(Lnhtemp,s.windows="period") %使⽤stl()函数对时间序列的三个分量进⾏分解
plot(fit)
fit$time.ries
exp(fit$time.ries) %对fit的time.ries分量进⾏指数化，就能回到对数前的数据⼤⼩
这段程序主要是想告诉⼤家如何对乘法模型和加法模型进⾏转化，其实适⽤季节性模型的时间序列使⽤加法模型和乘法模型都是可以的。使⽤stl()函数对时间序列对象进⾏分解主要是因为该函数只能处理加法模型，但这种条件现在看来其实不是很⼤的⼀种限制，因为乘法模型总可以通过对数变换转换成加法模型。
2.4对序列是否平稳进⾏判断
在时间序列这个科学概念中，我们要进⾏统计模型的拟合，简单地差分，增强时间序列对象在统计分析上的可⽤性，就必须要知道我们得到的时间序列是否平稳，若序列不是平稳的，我们要进⾏差分，差分的次数也需要确定，若序列是平稳的，那么，该时间序列我们就可以直接使⽤来进⾏后⾯统计的处理，很庆幸的是，基于R语⾔我们可以仅仅通过ndiffs()函数来判断我们的⽬标序列是否是平稳的，情况如下:
上⾯这段代码是什么意思呢？在使⽤ndiffs()函数的时候要记得导⼊forecast包进⾏备⽤，这是当你使⽤⼀些函数时必须完成的动作，如果不载⼊，系统就会说没有这个函数，谨记！
接着我们使⽤ndiffs()函数，函数的对象是时间序列对象，我们使⽤的是R库⾃备的nhtemp数据集，不需要担⼼对象的类型问题。我们可以看到该函数返回的值为1，就说明该时间序列不是平稳的，返回的值就是该序列需要差分的次数。
那么，到这⾥⼤家可能已经明⽩了这个函数的⽤法:
1.若返回值为0，序列是平稳的，根本不需要差分！
2.若返回值为其他值，序列是不平稳的，返回的值就是需要差分的次数。
