#线性时间序列模型
#苹果公司 2007 年到 2017 年股票日收盘价 library(quantmod) AAPL <- getSymbols("AAPL", src="yahoo", auto.assign=FALSE) chartSeries(AAPL, type="line",major.ticks="years", minor.ticks=FALSE, TA=NULL,subset="2003/2017",theme="white", name="Apple") #股价序列呈现缓慢的、非单调的上升趋势,局部又有短暂的波动。
#可口可乐公司每季度发布的每股盈利数据 library(readr);library(xts);library(quantmod);library(lubridate) da <- read_table("q-ko-earns8309.txt",col_types=cols( .default = col_double(), pends = col_date("%Y%m%d"),anntime = col_date("%Y%m%d"))) ko.Rqtr <- xts(da[["value"]], da[["pends"]])#xts前变量后日期 chartSeries(ko.Rqtr, type="line", TA=NULL, major.ticks="years", minor.ticks=FALSE, theme="white", name="Coca Kola Quarterly Return") #序列仍体现出缓慢的、非单调的上升趋势, #又有明显的每年的周期变化(称为季节性),还有短期的波动。
#用不同颜色标出不同季节,plot() 作图,x是时间(连续的数值型序列) tmp.x <- year(index(ko.Rqtr)) + (quarter(index(ko.Rqtr))-1)/4 tmp.y <- c(coredata(ko.Rqtr)) plot(tmp.x, tmp.y, type="l", col="gray") cpal <- c("green", "red", "yellow", "black") points(tmp.x, tmp.y, pch=16, col=cpal[quarter(index(ko.Rqtr))]) legend(locator(2),cex=0.7,xjust=0,yjust=1,bty="n", text.width=40,inset=5, pch=16, col=cpal, legend=c("Spring", "Summer", "Autumn", "Winter")) #每年一般冬季和春季最低,夏季最高,秋季介于夏季和冬季之间。
#标普 500 指数月对数收益率 d <- read_table("m-ibmsp-2611.txt",col_types=cols( .default=col_double(),date=col_date(format="%Y%m%d"))) sp.rmon <- xts(log(1 + d[["sp"]]), d[["date"]]) chartSeries(sp.rmon, type="line", TA=NULL,theme="white", major.ticks="auto", minor.ticks=FALSE, name="S&P 500 Monthly Returns") #收益率在 0 上下波动,除了个别时候基本在某个波动范围之内。 #由弱平稳性,可以对未来的标普 500 收益率预测如下: #均值在 0 左右,上下幅度在 ±0.2 之间。
#美国国债 3 月期和 6 月期周利率 d <- read_table2("w-tb3ms.txt",col_types=cols(.default=col_double())) x1 <- xts(d[["rate"]], make_date(d[["year"]], d[["mon"]], d[["day"]])) d <- read_table2("w-tb6ms.txt",col_types=cols(.default=col_double())) x2 <- xts(d[["rate"]], make_date(d[["year"]], d[["mon"]], d[["day"]])) tb36ms <- merge(x1, x2)#或用cbind,效果一致 names(tb36ms) <- c("tb3ms", "tb6ms")#给列命名 rm(d, x1, x2)#只要tb36ms
plot(tb36ms, type="l", grid.ticks.on="years") plot(tb36ms, type="l", subset="2004") plot(tb36ms, type="l", subset="1980")
#平稳性
#相关系数 d<-read.table('m-ibmsp6709.txt',header=T) ibm=d$ibm;sp5=d$sp; cor(sp5,ibm)#等同于 cor(d[,"sp"], d[,"ibm"]) cor(sp5,ibm,method='spearman') cor(sp5,ibm,method='kendall')
# acf(x) 可以估计时间序列 x 的自相关函数并对其前面若干项画图
d <- read_table2("m-dec12910.txt",col_types=cols( .default=col_double(),date=col_date(format="%Y%m%d"))) #19670131这样的被date=col_date(format="%Y%m%d")改成1967-01-31 dec <- xts(as.matrix(d[,-1]), d$date) #变量去掉日期列,日期列放在xts函数的第二项 tclass(dec) <- "yearmon"#修改行名(日期)格式 #'indexClass<-' is deprecated. Use 'tclass<-' instead.
d10 <- ts( coredata(dec)[,"dec10"], start=c(1967,1), frequency=12) plot(d10, main="CRSP Lower 10% Mothly Returns")
#用 acf() 作时间序列的自相关函数图: acf(d10) acf(d10,plot=FALSE) plot(acf(d10,plot=FALSE))#效果与acf一致
#ACF 图中横轴上下两条水平下 #是在独立同分布白噪声假设下的加减两倍标准差 #如果独立同分布白噪声假设成立, #每个 ˆρk 有 95% 以上的概率落入这两条线之间。 # ACF图 k = 0 处总对应 ˆρ0 = 1 。#因为k=0时ρk=ρ0=1
#上图的 ˆ ρ 1 和ˆ ρ 12 都超出了界限 #(因为是月度数目,横轴的单位是 1/12 为一个时间点)。 #从此图可以认为此投资组合的收益率不是白噪声。
win.graph(4,4,4) par(mfcol=c(2,1)) plot(d10) title(main='simple return') acf(d10,lag=24,main=' Autocorrelation function')