数学建模论文-基于SAS分析的山西省饮料销售量季节模型时序预测研究.doc

基于SAS分析的山西省饮料销售量季节模型时序预测研究 一、 前言 本文以山西省的软饮料销售量为研究对象，着重应用 SAS软件，运用时间序列分析方法，通过研究各季度饮料销售量，预测未来市场走向，并以此作为事前信息，及时调整企业生产安排，为企业抢占市场份额，获取高额利润提供及时快速信息。在国家食品及零售业“十二五”规划支持和居民收入水平不断提高的影响下，我国饮料行业产销量呈现明显增长态势。根据ChinaVenture投中集团旗下金融数据产品CVSource 统计显示，2006至2011 年的6 年间，我国饮料市场产量从2006 年4220万吨增至2011年11762万吨，年复合增长率为22.75%； 饮料销量从2006年4196万吨增至2011年11622万吨，年复合增长率为22.59%。历年数据中软饮料产销量基本持平，行业并未出现库存积压现象，虽在2010年后产销量增速有所放缓，但行业整体规模增长态势明显，并将在未来几年持续该趋势（见图1）。图1 2006年至2011年中国饮料产销量规模趋势图饮料（soft drink) 是酒精含量低于0.5%( 质量比)的天然的或人工配制的饮料。又称清凉饮料、无醇饮料。所含酒精限指溶解香精、香料、色素等用的乙醇溶剂或乳酸饮料生产过程的副产物。饮料市场可以分为饮料、酒精饮料、固体饮料三类，本文中所定义的软饮料主要分为果蔬汁、瓶装水、碳酸饮料、功能饮料、即饮茶饮料、即饮咖啡饮料六大类（见图2）。图2 软饮料分类图【关键词】SAS软件 时间序列分析 预测 软饮料 二、数据来源 CNKI中国统计年鉴数据库：/kns55/Dig/dig.aspx#数据跨越2004年到最近的2011年第四季度（共8年32个数值）山西省饮料销售量季节数据（如下表）： 表1 数据表季度年份第1季度第2季度第3季度第4季度200422.626.334.943.8200528.138.438.584.6200661.667.973.649.7200754.586.467.156.6200872.388.587.191.9200987.6105.8121.6106.42010103.9132.3151.9144.62011119.5138.1153.1181.9单位：万吨 数据来源：CNKI中国统计年鉴数据库 三、数据分析绘制时序图，在SAS软件输入代码，data a; /*建立一个名字为a的临时数据集 */ input beverage; /* 变量取名为“beverage”将会以行的方式被读取，系统自动将它视为数值型变量，自动读取任意长度的变量值*/ t=intnx(quarter,1jan2004d,_n_-1); /*等时间间隔是quarter（季度）；依照时间是1jan2004d；依次自动生成时间数据 */ format t yyq4.; /*“t”变量的输出格式是字符长度为4的年度月度数据 */ cards; /* 开始输入数据*/ 22.626.334.943.828.138.438.584.661.667.973.649.754.586.467.156.672.388.587.191.987.6105.8121.6106.4103.9132.3151.9144.6119.5138.1153.1181.9 ; proc gplot data=a; /*对临时数据集a中的数据绘图*/plot beverage*t; /*输出以t为横坐标，beverage为纵坐标的曲线*/symbol c=black i=join v=star; /*曲线的颜色是黑色，观察值之间线性连接，各点的形状为星号*/ run; /*可以运行 */输出结果如图3：图3 原始序列时序图由图3可以看出该序列不是平稳序列，故对其一阶差分，在SAS软件输入代码，data b; input beverage; t=intnx(quarter,1jan2004d,_n_-1); format t yyq4.; dif=dif(beverage); /*对beverage进行1阶差分，差分后的序列值赋值给dif*/cards; 22.626.334.943.828.138.438.584.661.667.973.649.754.586.467.156.672.388.587.191.987.6105.8121.6106.4103.9132.3151.9144.6119.5138.1153.1181.9 ; proc gplot data=b; plot dif*t; /*输出以t为横坐标，dif为纵坐标的曲线*/symbol v=circle i=join c=black; /*曲线各点的形状为圆圈，观察值之间线性连接，曲线的颜色是黑色*/ run;输出结果如图4：图4 一阶差分后时序图由图4 可以看出该序列平稳性不显著并且具有很强的季节效应，其周期为4，因此再对其进行4步差分，在SAS软件输入代码，data c; input beverage; t=intnx(quarter,1jan2004d,_n_-1); format t yyq4.; dif1=dif4(dif(beverage); /*对beverage进行1阶差分后，再进行4步差分，并将差分后的序列值赋值给变量dif1*/ cards; 22.626.334.943.828.138.438.584.661.667.973.649.754.586.467.156.672.388.587.191.987.6105.8121.6106.4103.9132.3151.9144.6119.5138.1153.1181.9 ; proc gplot data=c; plot dif1*t; /*输出以t为横坐标，dif1为纵坐标的曲线*/symbol v=circle i=join c=black; /*各点的形状为圆圈，观察值之间线性连接，曲线的颜色是黑色*/ run;输出结果如图5：图5 一阶4步差分后时序图由图3可以看出，该序列没有显著的不平稳性，认为一阶4歩差分后的序列平稳。做平稳性检验，在SAS软件输入代码，proc arima data=c; /*对临时数据集c中的数据进行ARIMA分析*/identify var=dif1 nlag=9; /*对进行了1阶4步差分后的dif1进行识别,延迟期数到第9阶*/run;生成自相关图（图6）、偏自相关图（图7）以及纯随机检验结果（图8）：图6 自相关系数图图7 偏自相关系数图图8 纯随机检验由图8可以看出延迟6阶的检验P值为0.0152小于给定的显著性水平，因此，拒绝原假设，认为该序列为非白噪声序列。所以，对该序列建模是有意义的。四、模型的建立及模型的检验（一）、模型的识别优化观察图7样本偏自相关系数图可知滞后2阶和4阶偏自相关系数都落在两倍标准差外。因此，尝试对一阶4步差分后的序列拟合AR（2 4）模型。在SAS软件输入代码，对模型进行估计。具体代码如下：proc arima data=c; estimate p=(2 4); /*用最小二乘估计方法（系统自动默认）对AR（2 4）模型进行参数估计*/forecast lead=5 id=time out=out; /*对序列进行5期预测，指定变量的序列观测值标识时间的id变量，并把数据集储存到out中*/run：（二）、模型参数的估计由以上代码，输出如下的模型的估计值：图9由图9可以看到参数的估计方法为条件最小二乘估计法；由图可以得出MU、AR1,1和AR1,2的估计值分别为0.87795、-0.43585和-0.54998；图10 拟合优度统计量表通过拟合优度统计量表可以看出相关统计量，这些统计量可以帮助比较该模型和其他模型的优劣。其中“Constant Estimate”表示的为均值项MU和自回归参数的函数；“Variance Estimate”表示残差序列的方差；“Std Error Estimate”代表方差估计值的平方根；AIC和SBC函数值的大小分别为238.932和242.8195；“Numbers of Residuals”表示的是残差个数，本例残差个数为27个。图11 参数估计值的相关系数表通过参数估计值的相关系数表可以帮助我们了解参数相关性可能影响结果的程度。从该表可以发现，任何两参数估计值的相关性都不高。图12 残差序列检验值表通过残差序列检验值表来检验残差序列是否为白噪声序列，从而检验模型的显著性。由表可以看出延迟6、12和18期的P值都明显大于0.05，认为残差序列为白噪声序列，虽然延迟24期的P值都明显小于0.05，但我们仍可以认为模型拟合较为良好。 图13 拟合模型参数值图13输出的是拟合模型的具体形式。其中，均值的估计值为0.877949；在本图下一部分显示的是自相关因子。得到的模型表达式如下：五、模型的预测图14 序列预测值 由图14可以看到由代码：“forecast lead=5 id=time out=out;”语句输出5期的预测值，其未来5期预测值分别为：172.9778，182.9771，195.1994，209.6372和194.7717。proc gplot data=out; /*对临时数据集out中的数据绘图*/plot x*time=1 forecast*time=2 /overlay; /*第一条以time为横坐标，x为纵坐标；第二条以time为横坐标，forecast为纵坐标；overlay 选项指令系统将这两条时序线绘制在同一张图中，同时显示*/symbol1 c=black i=none v=star; /*曲线1的颜色是黑色，观察值之间不做人任何连接，各点的形状为星号 */symbol2 c=red i=join v=none; /*曲线2的颜色是红色，观察值之间线性连接，各点不使用特别图形标注*/run;为了观察模型的拟合效果，本文将原序列图、序列拟合图、预测值95%置信下限和上限图画在一起，如图15：图15 序列拟合效果图在图15中，星号表示原序列值，曲线