数学建模B作业：非参数统计、灰色系统、时间序列分析010

1、2014年数学建模B作业：非参、灰色、时间序列分析非参数统计-1 某制造商想要比较两种不同的生产方法所花费的生产时间是否有差异。随机地选取了11个工人，每一个工人都分别使用两种不同的生产方法来完成一项相同的任务，在样本中的每一个工人都做了观察。数据见表，试用Wilcoxon秩和检验这两种方法有无差异？工人编号1234567891011方法110.29.69.210.69.910.210.610.011.210.710.6方法29.59.88.810.110.39.310.510.010.610.29.8解：提出原假设，这两组方法没有显著性差异，用配对实验的符号检验法，相应代码如下：data e

2、x; input x1 x2; y=x1-x2; cards; 10.29.59.69.89.28.810.610.19.910.310.29.310.610.5101011.210.610.710.210.69.8;proc univariate; var y; run; 运行结果如下： 从结果中可以看出，sign统计量为3，其显著性为0.1094，大于0.05，故接受原假设，认为这两组方法没有显著性差异。-2 为培训大学生志愿者为社区服务，设计了4种培训方案，记作为A,B,C,D.将报名的30名大学生随机地分为4组，分别接受不同培训。训练一周后，按规定的要求考试，评定的成绩如下，试用非参数

3、检验方法检验这四种培训方案的有效性是否存在显著差异？培训方案A60,75,62,76,73,98,86培训方案B72,52,68,82,74,64,87培训方案C61,85,78,66,70,59,69,79培训方案D63,58,65,71,84,77,80,89解：提出原假设，这四种培训方案方法没有显著性差异，相应代码如下：data ex;do a=1 to 4;input n;do i=1 to n;input x;output;end;end; cards;7 60 75 62 76 73 98 867 72 52 68 82 74 64 878 61 85 78 66 70 59 69

4、 798 63 58 65 71 84 77 80 89;proc npar1way wilcoxon;class a;var x;run;运行结果如下： 从结果中可以看出，Chi-Square统计量为0.5537，其显著性为0.9069，大于0.05，故接受原假设，认为四种培训方案方法没有显著性差异。-3 双胞胎智力的相关分析某研究所对10对双胞胎儿童的智力进行调查，试计算其Pearson、Spearman和Kendall相关系数并对其进行相关性检验。双胞胎编号先出生儿童X后出生儿童Y19.07.8216.619.3316.220.1411.37.1516.213.067.14.877.88

5、.984.07.4911.210.0101.31.5解：求其Pearson，Spearman和Kendall相关系数，代码如下：DATA new;INPUT x y;CARDS;9.07.816.619.316.220.111.37.116.213.07.14.87.88.94.07.411.210.01.31.5;PROC CORR pearson spearman kendall;VAR x y;RUN;结果如下： Pearson Correlation Coefficients, N = 10 Prob > |r| under H0: Rho=0 x y x 1.00000 0.8

6、8081 0.0008 y 0.88081 1.00000 0.0008 Spearman Correlation Coefficients, N = 10 Prob > |r| under H0: Rho=0 x y x 1.00000 0.82067 0.0036 y 0.82067 1.00000 0.0036 Kendall Tau b Correlation Coefficients, N = 10 Prob > |r| under H0: Rho=0 x y x 1.00000 0.67420 0.0071 y 0.67420 1.00000 0.0071可见，x与y的

7、Pearson相关系数为0.88081，概率为0.0008，达到极显著水平；Spearman相关系数为0.82067，概率为0.0036，达到极显著水平；Kendall相关系数0.67420，概率为0.0071达到极显著水平；故，x与y显著相关。灰色系统作业：-4 陕西省农业总产值数据如下：已知下表数据年份1985198619871888198919901991199219931994总产值62.958.861.487.2104.9124.8110.7129.0155.3219.03请建立灰色系统GM（1，1）模型，并预测1995-1997三年的农业总产值。解：有原始时间1985-1994序列

8、，对生成1-AGO序列另外可得Yn见表：、1-AGO序列、Ynk1234567891062.958.861.487.2104.9124.8110.7129155.3219.0362.9121.7183.1270.3375.2500610.7739.78951114.03Yn58.861.487.2104.9124.8110.7129155.3219.03利用MATLAB编程得：function X,c,error1,error2=example9_11()%利用MATLAB编程预测2003年中国蔬菜产量，%并对预测结果做残差检验和后验差检验，程序如下：X0=62.958.861.487.210

9、4.9124.8110.7129.0155.3219.03;k=3;X,c,error1,error2=GM11(X0,k)plot(1985:1994,X0,'g*-')hold on plot(1985:1997,X)%function X,c,error1,error2=GM11(X0,k)% 建立函数X,c,error1,error2=example9_3_2_3(X0,k)% 其中X0为输入序列，k为预测长度，% X为预测输出序列，c为后验差检验数，error1为残差，error2为相对误差format long;n=length(X0);X1=;X1(1)=X0(1

10、);for i=2:n X1(i)=X1(i-1)+X0(i); %计算累加生成序列endfor i=1:n-1B(i,1)=-0.5*(X1(i)+X1(i+1); %计算B，Yn B(i,2)=1; Y(i)=X0(i+1);endalpha=(B'*B)(-1)*B'*Y' %做最小二乘估计a=alpha(1,1);b=alpha(2,1);d=b/a; %计算时间响应函数参数c=X1(1)-d;X2(1)=X0(1);X(1)=X0(1);for i=1:n-1 X2(i+1)=c*exp(-a*i)+d; X(i+1)=X2(i+1)-X2(i); %计算预测

11、序列endfor i=(n+1):(n+k) X2(i)=c*exp(-a*(i-1)+d; %计算预测序列 X(i)=X2(i)-X2(i-1);endfor i=1:n error(i)=X(i)-X0(i); error1(i)=abs(error(i); %计算残差 error2(i)=error1(i)/X0(i); %计算相对误差endc=std(error1)/std(X0); %计算后验差检验数运行结果见表格：年份1985198619871888198919901991199219931994199519961997总产值62.958.861.487.2104.9124.811

12、0.7129155.3219.03预测值62.958.5832668.1568779.2949992.25329107.3292124.8688145.2748169.0154196.6357228.7697266.155309.6498残差00.2166.7567.90512.6417.4714.1616.2713.7122.39相对误差00.0030.1100.0900.1200.1390.1270.1260.0880.102画出预测与实际值变化曲线，如图所示：预测与实际值变化曲线实验模型以及结果检验：由表与图的结果可见，预测值与实际值偏离不大，其后验残差检验数C=0.1475小于0.35

13、，所以模型精度为优。时间序列分析作业-5 某车站1993-1997年各月的列车运行数量数据如下表，试用时间序列建立合适的模型。并预测1998年1月的数值1196.81181.31222.61229.31221.5 1148.41250.21174.41234.51209.71206.51204.01234.11146.01304.9 1221.91244.11194.41281.51277.31238.91267.51200.91245.51249.9 1220.11267.41182.31221.71178.11261.61274.51196.41222.61174.7 1212.61215

14、.01191.01179.01224.01183.01288.01274.01218.01263.0 1205.01210.01243.01266.01200.01306.01209.01248.01208.01231.0 1244.01296.01221.01287.01191.0 解：（1） 首先进行平稳性检验：data a;/*a为数据名*/input lieche;/*lieche为变量名*/month=intnx('month','1jan1993'd,_n_-1);/*intnx间隔取时间变量*/format month date.;/*月按？*/c

15、ards;1196.81181.31222.61229.31221.5 1148.41250.21174.41234.51209.71206.51204.01234.11146.01304.9 1221.91244.11194.41281.51277.31238.91267.51200.91245.51249.9 1220.11267.41182.31221.71178.11261.61274.51196.41222.61174.7 1212.61215.01191.01179.01224.01183.01288.01274.01218.01263.0 1205.01210.01243.012

16、66.01200.01306.01209.01248.01208.01231.0 1244.01296.01221.01287.01191.0 ;run ;proc gplot;/*画图*/plot lieche*month;/*纵轴为lieche,横轴为mouth*/symbol v=square i=join c = red;/*图形特征，v表示点的形状，i表示图形连线的情况，c代表颜色*/proc arima data = a;/*调用arima模块*/identify var=lieche nlag = 22;/*延迟阶数为22阶*/run;运行得自相关图：由此自相关图可看出，自相关系

17、数很快的衰减向0，且始终控制在2倍范围内，可以认为该序列为平稳序列。时序图：由图可知，此车站列车运行数量数据在一个常数值附近随机波动，而且波动范围有界，无明显趋势及周期特征，基本可以视序列为平稳序列。（2）进行随机性检验：选取结果中The ARIMA Procedure部分：由于统计量P值均大于0.05，则认为在0.05的显著水平下，无法拒绝原假设，即不能显著拒绝序列为纯随机序列的假定，因而认为此车站列车运行数量为纯随机波动序列，各序列之间没有任何行相关关系，即为无记忆序列，也就是说，该车站列车运行数量前后两年并无大的联系，也就是实说，我们很难根据历史信息预测未来年份此车站列车运行数量，故，该

18、平稳序列不值得继续分析下去，对该序列分析到此结束。-6 对我国1952-1994年的社会消费品零售总额数据建立合适的时间序列模型，并预测1995-1997年的数据。社会消费品零售总额1952262.7328.8356.11955364.0424.0441.6481.2556.51960595.4537.7543.7544.8572.71965590.1632.8679.1649.2698.21970728.8776.9853.5917.7967.419751046.41099.01174.31264.91476.019801794.02002.52181.52426.12899.2198538

19、01.44374.05115.06534.67074.219907250.38245.79704.812462.116264.7解：（1）首先进行平稳性检验：data a;/*a为数据名*/input xf;/*xf为变量名*/year=intnx('year','1jan1952'd,_n_-1);/*intnx间隔取时间变量*/format year year4.; /*年按四位数显示*/cards;262.7328.8356.1364.0424.0441.6481.2556.5595.4537.7543.7544.8572.7590.1632.8679.1

20、649.2698.2728.8776.9853.5917.7967.41046.41099.01174.31264.91476.01794.02002.52181.52426.12899.23801.44374.05115.06534.67074.27250.38245.79704.812462.116264.7;run ;proc gplot;/*画图*/plot xf*year;symbol v=square i=join c = red;/*图形特征，v表示点的形状，i表示图形连线的情况，c代表颜色*/proc arima data = a;/*调用arima模块*/identify v

21、ar=xf nlag = 22;/*延迟阶数为22阶*/run;首先分析时序图：由时序图可得，该时间序列显著递增，初步判断此序列不平稳。再分析自相关图：由自相关图中，自相关系数从正数缓慢递减为到零后，又不断在负值范围内增大，该序列自相关系数并未较快的衰减为零，因此该序列并非为平稳时间序列。（2） 随机性检验：选取结果中The ARIMA Procedure部分：从运行结果得出，次统计量的P值均小于0.0001，则认为在0.05的显著水平下拒绝原假设，可以认为此序列为非随机序列。这说明我们可以根据历史是信息预测未来年份我国的社会消费品零售总额。（3） 模型选取原序列自相关系数拖尾，偏自相关系数一

22、阶截尾，根据ARMA模型相关性特征表，应该选取AR(1)模型。 首先对其进行一阶差分:data a;/*a为数据名*/input xf;/*lieche为变量名*/year=intnx('year','1jan1952'd,_n_-1);/*intnx间隔取时间变量*/format year year4.; /*年按四位数显示*/dif1=dif(xf);cards;262.7328.8356.1364.0424.0441.6481.2556.5595.4537.7543.7544.8572.7590.1632.8679.1649.2698.2728.8776.

23、9853.5917.7967.41046.41099.01174.31264.91476.01794.02002.52181.52426.12899.23801.44374.05115.06534.67074.27250.38245.79704.812462.116264.7;run ;proc gplot;/*画图*/plot dif1*year ;symbol v=square i=join c = red;/*图形特征，v表示点的形状，i表示图形连线的情况，c代表颜色*/proc arima data = a;/*调用arima模块*/identify var=xf ;proc arima data = a;/*调用arima模块*/identify var=dif1 ;run;得到时序图如下：又该图可以简单看出差分后，数据在某个数据间波动，范围有界，无明显趋势及周期性特征，初步判断一阶差分后序列平稳。（4） 模型建立：选取AR(1)摸型data a;/*aΪÊý¾ÝÃû*/input xf;year=intnx('year',&#