0 library
library(ggplot2)
library(gridExtra)
library(stats)
1 MA(1)模型的ACVF,ACF,PACF
1.1生成高斯白噪声
set.seed("1234")
n=250
a=rnorm(n)
1.2生成含250个值的序列,分布如图
u=0
Theta=0.5
z=array(NA,n)
z[1]=u+a[1]
for(t in 2:n)
z[t]=u+a[t]-Theta*a[t-1]
plot(z,type = "l")
1.3手工计算ACVF,ACF估计值
Ez=mean(z)
Gamma=array(NA,n-1)
Rou=array(NA,n-1)
sum=0
for(t in 1:n)
sum=sum+(z[t]-Ez)*(z[t]-Ez)
Gamma0=1/n*sum
for(k in 1:(n-1)){
sum=0
for(t in 1:(n-k))
sum=sum+(z[t]-Ez)*(z[t+k]-Ez)
Gamma[k]=1/n*sum
Rou[k]=Gamma[k]/Gamma0
}
No=0:(n-1)
data=data.frame(no = No, gamma=c(Gamma0,Gamma),rou = c(1,Rou))
p <-ggplot(data)
p1=p+geom_bar(width=0.1,mapping=aes(x=no, y=rou),stat = "identity")+
xlab("NO")+ ylab("ACF")+ xlim(-1,24)+geom_hline(yintercept = 0)
p2=p+geom_bar(width=0.1,mapping = aes(x=no, y=gamma),stat ="identity")+
xlab("NO")+ ylab("ACVF")+ xlim(-1,24)+geom_hline(yintercept = 0)
grid.arrange(p1, p2, nrow=1, ncol=2)
1.4 acf函数生成
acf(z,type = "correlation")
acf(z,type = "covariance")
acf(z,type = "partial")
2.1 AR(1)模型的标准谱密度函数
f <- function(w,b) (1/(2*pi))*((1-b^2)/(1+b^2-2*b*cos(w)))
# 创建数据点
w <- seq(0, pi, by=0.01)
y1 <- f(w,0.7)
y2 <- f(w,0.5)
y3 <- f(w,0.3)
df <- data.frame(w, y1, y2, y3)
# 用ggplot2来画图
g <- ggplot(df)
g1 <- g + geom_line(col='red',mapping=aes(w,y1),linetype = 1)+
geom_line(col='blue',mapping=aes(w,y2),linetype = 2)+
geom_line(col='black',mapping=aes(w,y3),linetype = 3)+
geom_hline(yintercept = 0)+geom_vline(xintercept = 0)+
ylab("f(w)")+
ggtitle("f(w)")
print(g1)
2.2谱密度估计
三种常用方法:
1.(最小二乘)平滑周期图
对有限样本,异方差性影响平滑的有效性
2.(最小二乘)平滑对数周期图
3.局部似然估计(推荐这个,因为前两个不是有效的估计)
## Examples from Venables & Ripley
## spec.pgram
##Allowed methods are "pgram" (the default) and "ar".
##这是啥方法
par(mfrow = c(2,2))
spectrum(lh)
spectrum(lh, spans = 3)
spectrum(lh, spans = c(3,3))
spectrum(lh, spans = c(3,5))
spectrum(ldeaths)
spectrum(ldeaths, spans = c(3,3))
spectrum(ldeaths, spans = c(3,5))
spectrum(ldeaths, spans = c(5,7), log = "dB", ci = 0.8)
# for multivariate examples see the help for spec.pgram
## spec.ar
spectrum(lh, method = "ar")
spectrum(ldeaths, method = "ar")
关于R新增的了解
1.```{r,warning=FALSE,dev=‘svg’} 输出矢量图
2.ggplot2 的使用
3.使用ggplot2时,常用的par不能用了,可以使用gridExtra包里的grid.arrage分页
4.acf,spectrum函数的使用
5.ggplot2的坐标轴画法
6.set.seed(“1234”) 设置伪随机数起点,使每次模拟结果相同
另外:结构断点与谱分析http://www.xzbu.com/3/view-4489847.htm