时间序列作业

 

第二次上机实验：模拟MA模型和ARMA模型

 作业内容请访问：http://www.doc88.com/p-209946619868.html

参考解答：

#2--5上机结果n<-100;     x1<-0;a2<-0;a3<-0;x<-0;      for(i in -50:n)          {a1<-rnorm(1);          x1<-a1+0.65*a2+0.24*a3;a3<-a2;a2<-a1;                    if(i>0)x[i]<-x1}        xmean(x)#-0.007158113var(x)#1.532393fivenum(x)#-2.33380162 -0.96968524  0.01153527  0.78820663  2.37949770plot(x,type="b",xlab="i",ylab=expression(X[i]),col=2)#see picture1par(mfrow=c(1,2))acf(x)pacf(x)#see picture 2#6set.seed(100)n<-100;     x1<-0;a2<-0;a3<-0;x<-0;      for(i in -50:n)          {a1<-rnorm(1);          x1<-a1-0.65*a2-0.24*a3;a3<-a2;a2<-a1;                    if(i>0)x[i]<-x1}        xmean(x)#0.009997361var(x)# 1.745586fivenum(x)#-3.54820575 -0.91206753 -0.02810196  0.85317906  3.06819255par(mfrow=c(1,1))plot(x,type="b",xlab="i",ylab=expression(X[i]),col=2)#see picture3par(mfrow=c(2,1))acf(x)pacf(x)#see picture 4#7set.seed(100)n<-100;     x1<-0;a2<-0;a3<-0;x<-0;      for(i in -50:n)          {a1<-rnorm(1);          x1<-a1-0.65*a2+0.24*a3;a3<-a2;a2<-a1;                    if(i>0)x[i]<-x1}        x;mean(x)#-0.04177703var(x)# 1.864563fivenum(x)#-3.46235063 -1.03775369 -0.06152364  0.93842732  2.96942047par(mfrow=c(1,1))plot(x,type="b",xlab="i",ylab=expression(X[i]),col=2)#see picture 5par(mfrow=c(2,1))acf(x)pacf(x)#see picture 6#8set.seed(100)x<-arima.sim(n = 100, list( ma = c(0.65, -0.24)))mean(x)#-0.03727247var(x)# 1.401805fivenum(x)#-3.54417536 -0.69506541 -0.04359503  0.66448787  3.33484756par(mfrow=c(1,1))plot(x,type="b",xlab="i",ylab=expression(X[i]),col=2)#see picture 7par(mfrow=c(2,1))acf(x)pacf(x)#see picture 8#8 use function in R directlyset.seed(100)x<-arima.sim(n = 100, list(ar=c( ma = c(0.65, -0.24)))mean(x)#-0.03727247var(x)# 1.401805fivenum(x)#-3.54417536 -0.69506541 -0.04359503  0.66448787  3.33484756par(mfrow=c(1,1))plot(x,type="b",xlab="i",ylab=expression(X[i]),col=2)#see picture 7par(mfrow=c(2,1))acf(x)pacf(x)#see picture 8#9set.seed(100)x<-arima.sim(n = 100, list(ar=c(-0.75,-0.5), ma = c(0.3, 0.4)))+5mean(x)#4.98var(x)# 1.692439fivenum(x)# 1.658756 4.142209 4.996332 5.941509 8.064161par(mfrow=c(1,1))plot(x,type="b",xlab="i",ylab=expression(X[i]),col=2)#see picture 9par(mfrow=c(2,1))acf(x)pacf(x)#see picture 10#10回到graph窗口观察各种序列图形的异同。#acf均是2步截尾的（ARMA不是）pacf均是拖尾的，arma的acf也是拖尾的#acf衰减情况不同，有指数衰减的，也有震荡衰减的#11 模拟其他随机序列，体会不同类型随机序列的样本自相关函数和样本偏自相关函数的差异，尤其是拖尾性与截尾性的区分。require(graphics)x1<-arima.sim(n = 63, list(ar = c(0.8897, -0.4858), ma = c(-0.2279, 0.2488)),          sd = sqrt(0.1796))# mildly long-tailedx2<-arima.sim(n = 63, list(ar=c(0.8897, -0.4858), ma=c(-0.2279, 0.2488)),          rand.gen = function(n, ...) sqrt(0.1796) * rt(n, df = 5))par(mfrow=c(2,1))acf(x1)pacf(x1)par(mfrow=c(2,1))acf(x2)pacf(x2)#截尾性与拖尾性其实也没有那么好判断，这与acf的衰减幅度有关。如果衰减幅度特别大的时候，用ARMA拟合效果也不错#思考：如何利用样本自相关函数和样本偏自相关函数识别序列类型？#看截尾性与拖尾性，也可以使用广义自相关函数eacf去判断