R语言常用上机命令分功能整理 - 时间序列分析为主

2020-04-14 01:52

第一讲 应用实例

? R的基本界面是一个交互式命令窗口,命令提示符是一个大于号,命令的结果马上

显示在命令下面。

? S命令主要有两种形式:表达式或赋值运算(用’<-’或者’=’表示)。在命令提示符后

键入一个表达式表示计算此表达式并显示结果。赋值运算把赋值号右边的值计算出来赋给左边的变量。

? 可以用向上光标键来找回以前运行的命令再次运行或修改后再运行。 ? S是区分大小写的,所以x和X是不同的名字。

我们用一些例子来看R软件的特点。假设我们已经进入了R的交互式窗口。如果没有打开的图形窗口,在R中,用:> x11() 可以打开一个作图窗口。然后,输入以下语句: x1 = 0:100

x2 = x1*2*pi/100 y = sin(x2)

plot(x2,y,type=\

这些语句可以绘制正弦曲线图。其中,“=”是赋值运算符。0:100表示一个从0到100 的等差数列向量。第二个语句可以看出,我们可以对向量直接进行四则运算,计算得到的x2 是向量x1的所有元素乘以常数2*pi/100的结果。从第三个语句可看到函数可以以向量为输入,并可以输出一个向量,结果向量y的每一个分量是自变量x2的每一个分量的正弦函数值。

plot(x2,y, type=\画图练习\好好练\轴\轴')

有关作图命令plot的详细介绍可以在R中输入help(plot)

数学函数

abs,sqrt:绝对值,平方根 log, log10, log2 , exp:对数与指数函数 sin,cos,tan,asin,acos,atan,atan2:三角函数 sinh,cosh,tanh,asinh,acosh,atanh:双曲函数 简单统计量

sum, mean, var, sd, min, max, range, median, IQR(四分位间距)等为统计量,sort,order,rank与排序有关,其它还有ave,fivenum,mad,quantile,stem等。

下面我们看一看S的统计功能: > marks <- c(10, 6, 4, 7, 8) > mean(marks) > sd(marks) > min(marks) > max(marks)

第一个语句输入若干数据到一个向量,函c()用来把数据组合为一个向量。后面用了几个函数来计算数据的均值、标准差、最小值、最大值。

可以把若干行命令保存在一个文本文件中,然后用source函数来运行整个文件: > source(\ 注意字符串中的反斜杠。

例:计算6, 4, 7, 8,10的均值和标准差,把若干行命令保存在一个文本文件(比如C:\\1.R)中,然后用source函数来运行整个文件。 a<- c(10, 6, 4, 7, 8) b<-mean(a) c<-sd(a)

source(\

时间序列数据的输入 使用函数ts

ts(1:10, frequency = 4, start = c(1959, 2))

print( ts(1:10, frequency = 7, start = c(12, 2)), calendar = TRUE) a<-ts(1:10, frequency = 4, start = c(1959, 2)) plot(a)

将外部数据读入R read.csv

默认header = TRUE,也就是第一行是标签,不是数据。 read.table

默认header = FALSE

将R中的数据输出 write

write.table write.csv

第二讲

1. 绘制时序图、自相关图 例题2.1

d=scan(\

sha=ts(d,start=1964,freq=1)

plot.ts(sha) #绘制时序图

acf(sha,22) #绘制自相关图,滞后期数22 pacf(sha,22) #绘制偏自相关图,滞后期数22 corr=acf(sha,22) #保存相关系数

cov=acf(sha,22,type = \ #保存协方差 图的保存,单击选中图,在菜单栏选中“文件”,再选“另存为”。 同时显示多个图:用x11()命令生成一个空白图,再输入作图命令。

2. 同时绘制两组数据的时序图 d=read.csv(\double=ts(d,start=1964,freq=1)

plot(double, plot.type = \ #两组数据两个图 plot(double, plot.type = \ #两组数据一个图

plot(double, plot.type = \设置每组数据图的颜色、曲线类型)

3.产生服从正态分布的随机观察值

例题2.4 随机产生1000白噪声序列观察值 d=rnorm(1000,0,1) #个数1000 均值0 方差1 plot.ts(d)

4.纯随机性检验 例题2.3续

d=scan(\

temp=ts(d,freq=1,start=c(1949))

Box.test(temp, type=\

5.差分计算 x=1:10 y=diff(x) k步差分

?kxt?xt?xt?k 加入参数 lag=k

如计算x的3步差分为 y=diff(x, lag = 3)

pp?1p?1?x??x??xt?1加入参数differences = p ttp阶差分

2?xt??xt??xt?1

如2阶差分

y=diff(x,differences = 2)

第三讲 例题3.1

plot.ts(arima.sim(n = 100, list(ar = 0.8))) #模拟AR(1)模型,并作时序图。

plot.ts(arima.sim(n = 100, list(ar = -1.1))) #非平稳,无法得到时序图。

plot.ts(arima.sim(n = 100, list(ar = c(1,-0.5)))) plot.ts(arima.sim(n = 100, list(ar = c(1,0.5))))

例题3.5

acf(arima.sim(n = 100, list(ar = 0.8))) acf (arima.sim(n = 100, list(ar = -1.1))) acf (arima.sim(n = 100, list(ar = c(1,-0.5)))) acf (arima.sim(n = 100, list(ar = c(1,0.5))))

例题3.7

arima.sim(n = 1000, list(ar = 0.5, ma = -0.8))

acf(arima.sim(n = 1000, list(ar = 0.5, ma = -0.8)),20) pacf(arima.sim(n = 1000, list(ar = 0.5, ma = -0.8)),20)

例题2.5

d=scan(\ #导入数据

prop=ts(d,start=1950,freq=1) #转化为时间序列数据 plot(prop) #作时序图

acf(prop,12) #作自相关图,拖尾

pacf(prop,12) #作偏自相关图,1阶截尾 Box.test(prop, type=\ #纯随机性检验,p值小于5%,序列为非白噪声 Box.test(prop, type=\arima(prop, order = c(1,0,0),method=\

#用AR(1)模型拟合,如参数method=\,估计方法为条件最小二乘法,用条件最小二乘法时,不显示AIC。

arima(prop, order = c(1,0,0),method=\用AR(1)模型拟合,不含截距项。 tsdiag(arima(prop, order = c(1,0,0),method=\ #对估计进行诊断,判断残差是否为白噪声

summary(arima(prop, order = c(1,0,0),method=\a=arima(prop, order = c(1,0,0),method=\r=a$residuals#用r来保存残差

Box.test(r,type=\对残差进行纯随机性检验

predict(arima(prop, order = c(1,0,0)), n.ahead =5) #预测未来5期 prop.fore = predict(arima(prop, order = c(1,0,0)), n.ahead =5) #将未来5期预测值保存在prop.fore变量中 U = prop.fore$pred + 1.96* prop.fore$se

L = prop.fore$pred – 1.96* prop.fore$se#算出95%置信区间

ts.plot(prop, prop.fore$pred,col=1:2)#作时序图,含预测。 lines(U, col=\

lines(L, col=\在时序图中作出95%置信区间

例题3.9

d=scan(\x=diff(d)

arima(x, order = c(1,0,1),method=\

tsdiag(arima(x, order = c(1,0,1),method=\

第一点:

对于第三讲中的例2.5,运行命令arima(prop, order = c(1,0,0),method=\之后,显示: Call:

arima(x = prop, order = c(1, 0, 0), method = \Coefficients:

ar1 intercept 0.6914 81.5509 s.e. 0.0989 1.7453

sigma^2 estimated as 15.51: log likelihood = -137.02, aic = 280.05

注意:intercept下面的81.5509是均值,而不是截距!虽然intercept是截距的意思,这里如果用mean会更好。(the mean and the intercept are the same only when there is no AR term,均值和截距是相同的,只有在没有AR项的时候)

如果想得到截距,利用公式计算。int=(1-0.6914)*81.5509= 25.16661。课本P81的例2.5续中的截距25.17是正确的。

第二点:

如需计算参数的t统计量值和p值,利用下面的公式。 ar的t统计量值=0.6914/0.0989= 6.9909

(注:数值与课本略有不同,因为课本用sas算的se= 0.1029,R计算的se=0.0989) p值=pt(6.9909,df=48,lower.tail = F)*2

pt()为求t分布求p值的函数,6.99为t统计量的绝对值,df为自由度=数据个数-参数个数,lower.tail = F表示所求p值为P[T > t],如不加入这个参数表示所求p值为P[T <=t]。 乘2表示p值是双侧的(课本上的p值由sas算出,是双侧的)

均值的t统计量值和p值同理。

在时间序列中对参数显著性的要求与回归模型不同,我们更多的是考察模型整体的好坏,而不是参数。所以,R中的arima拟合结果中没有给出参数的t统计量值和p值,如果题目没有特别要求,一般不需要手动计算。

第三点:

修正第三讲中的错误:


R语言常用上机命令分功能整理 - 时间序列分析为主.doc 将本文的Word文档下载到电脑 下载失败或者文档不完整,请联系客服人员解决!

下一篇:草坪、时令草花及水生植物

相关阅读
本类排行
× 注册会员免费下载(下载后可以自由复制和排版)

马上注册会员

注:下载文档有可能“只有目录或者内容不全”等情况,请下载之前注意辨别,如果您已付费且无法下载或内容有问题,请联系我们协助你处理。
微信: QQ: