时间序列模型参数估计
1 理论基础
1.1 矩估计
1.1.1 AR 模型
矩估计法参数估计的思路：
即从样本中依次求中r k 然后求其对应的参数Φk 值
方差：
1.1.2 MA 模型
对于MA 模型采用矩估计是比较不精确的，所以这里不予讨论
1.1.3 ARMA （1，1）
矩估计法参数估计的思路：
方差：
1.2 最小二乘估计
1.2.1 AR 模型
最小二乘参数估计的思路：
对于AR （P ）而言也可以得到类似矩估计得到的方程，即最小二乘与矩估计得到的估计量相同。
1.2.2 MA 模型
最小二乘参数估计的思路：
1.2.3 ARMA 模型
最小二乘参数估计的思路：
1.3 极大似然估计与无条件最小二乘估计
2 R 中如何实现时间序列参数估计
2.1 对于AR 模型
ar(x, aic = TRUE, order.max = NULL,
method=c("yule-walker", "burg", "ols", "mle", "yw"), na.action, series, ...)
> ar(ar1.s,order.max=1,AIC=F,method='yw')#即矩估计 Call:
ar(x = ar1.s, order.max = 1, method = "yw", AIC = F) Coefficients: 1
0.8314
Order selected 1 sigma^2 estimated as 1.382
> ar(ar1.s,order.max=1,AIC=F,method='ols')#最小二乘估计 Call:
ar(x = ar1.s, order.max = 1, method = "ols", AIC = F) Coefficients: 1
0.857
Intercept: 0.02499 (0.1308)
Order selected 1 sigma^2 estimated as 1.008
> ar(ar1.s,order.max=1,AIC=F,method='mle')#极大似然估计 Call:
ar(x = ar1.s, order.max = 1, method = "mle", AIC = F) Coefficients: 1
0.8924
Order selected 1 sigma^2 estimated as 1.041
采用自编函数总结三个不同的估计值
> Myar(ar2.s,order.max=3)
最小二乘估计 矩估计 极大似然估计
1 1.5137146 1.4694476 1.5061369
2 -0.8049905 -0.7646034 -0.7964453
2.2 对于ARMA 模型
arima(x, order = c(0, 0, 0), seasonal = list(order = c(0, 0, 0), period = NA),
xreg = NULL, include.mean = TRUE, transform.pars = TRUE, fixed = NULL, init = NULL, method = c("CSS-ML", "ML", "CSS"), n.cond, optim.control = list(), kappa = 1e+06, io = NULL, xtransf, transfer = NULL)
order 的三个参数分别代表AR ，差分 MA的阶数
> arima(arma11.s,order=c(1,0,1),method='CSS') Call: arima(x = arma11.s, order = c(1, 0, 1), method = "CSS") Coefficients: ar1 ma1 intercept
0.5586 0.3669 0.3928 s.e. 0.1219 0.1564 0.3380
sigma^2 estimated as 1.199: part log likelihood = -150.98
> arima(arma11.s,order=c(1,0,1),method='ML') Call: arima(x = arma11.s, order = c(1, 0, 1), method = "ML") Coefficients:
ar1 ma1 intercept
0.5647 0.3557 0.3216
s.e. 0.1205 0.1585 0.3358
sigma^2 estimated as 1.197: log likelihood = -151.33, aic = 308.65
采用自编函数总结三个不同的估计值
> Myarima(arma11.s,order=c(1,0,1)) $coef
条件SS 估计 极大似然估计 条件似然估计
ar1 0.5585828 0.5647477 0.5647498
ma1 0.3668814 0.3556965 0.3556973
intercept 0.3927654 0.3216166 0.3216152
$log
条件SS 估计 极大似然估计 条件似然估计
[1,] -150.984 -151.3268 -151.3268
$sigma2
条件SS 估计 极大似然估计 条件似然估计
[1,] 1.199378 1.196984 1.196984
$aic
条件SS 估计 极大似然估计 条件似然估计
[1,] NA 308.6537 308.6537
2.3 采用自助法arima.boot()
此函数估计的是参数的取值置信区间，而不是指具体的某个值，与arima 是不同的。
> res=arima(sqrt(hare),order=c(3,0,0),include.mean=T) > set.seed(12345)
> # Method I以最初三个观测为条件，并假设误差服从正态分布，得到95%的置信区间quantile 用于计算置信区间值，signif 类似于四舍五入函数，保留有效数值。
>
coefm.cond.norm =arima.boot(res,cond.boot=T,is.normal=T,B=1000,init=sqrt(hare))
> signif(apply(coefm.cond.norm,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3) ar1 ar2 ar3 intercept noise var 2.5% 0.593 -0.667 -0.6740 5.12 0.548
97.5% 1.280 0.244 -0.0135 6.38 1.540 >
> # Method II 假设误差并不服从正态分布，而是需要从样本抽样中得到coefm.cond.replace =arima.boot(res,cond.boot=T,is.normal=F,B=1000,init=sqrt(hare))
>
signif(apply(coefm.cond.replace,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3)
ar1 ar2 ar3 intercept noise var
2.5% 0.611 -0.700 -0.6720 4.98 0.516
97.5% 1.300 0.241 -0.0417 6.32 1.500
> # Method III基于平稳自助法的置信区间，且误差服从正态分布 >
coefm.norm =arima.boot(res,cond.boot=F,is.normal=T,ntrans=100,B=1000,init=sqrt(hare))
>
signif(apply(coefm.norm,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3)
ar1 ar2 ar3 intercept noise var
2.5% 0.687 -0.747 -0.6600 4.99 0.508 97.5% 1.380 0.192 -0.0168 6.33 1.500 >
> # Method IV基于平稳自助法的置信区间，且误差不服从正态分布
coefm.replace =arima.boot(res,cond.boot=F,is.normal=F,ntrans=100,B=1000,init=sqrt(hare))
>
signif(apply(coefm.replace,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3)
ar1 ar2 ar3 intercept noise var
2.5% 0.70 -0.715 -0.6620 4.98 0.47
97.5% 1.36 0.183 -0.0187 6.30 1.50
3 附自编函数
3.1 Myar
#用于自回归模型的参数估计，整合矩估计，最小二乘估计，以及极大似然估计 #该函数用于对时间序列中心化数据（因此截距项一定为0）估计AR 模型的参数，AIC 为真时，滞后项根据AIC 准则确定，为假时则根据设置的order.max 设定
Myar=function(tsdata, order.max = 1,AIC = F){
library(TSA)
ols
yw
mle
olscoef
ywcoef
mlecoef
result=data.frame(olscoef,ywcoef,mlecoef)
colnames(result)=c('最小二乘估计',' 矩估计',' 极大似然估计') return(result)
}
3.2 Myarima
#用于自回归模型的参数估计，整合矩估计，最小二乘估计，以及极大似然估计 #该函数用于对时间序列中心化数据（因此截距项一定为0）估计AR 模型的参数，AIC 为真时，滞后项根据AIC 准则确定，为假时则根据设置的order.max 设定
Myarima=function(tsdata, order=c(0,0,0)){
library(TSA)
result=NULL
css
ml
cssml
result$coef=cbind(css$coef,ml$coef,cssml$coef)
result$log=cbind(css$log,ml$log,cssml$log)
result$sigma2=cbind(css$sigma2,ml$sigma2,cssml$sigma2)
result$aic=cbind(NA,ml$aic,cssml$aic)
colnames(result$coef)=c('条件SS 估计',' 极大似然估计',' 条件似然估计') colnames(result$log)=c('条件SS 估计',' 极大似然估计',' 条件似然估计') colnames(result$aic)=c('条件SS 估计',' 极大似然估计',' 条件似然估计') colnames(result$sigma2)=c('条件SS 估计',' 极大似然估计',' 条件似然估计')
return(result) }
时间序列模型参数估计
1 理论基础
1.1 矩估计
1.1.1 AR 模型
矩估计法参数估计的思路：
即从样本中依次求中r k 然后求其对应的参数Φk 值
方差：
1.1.2 MA 模型
对于MA 模型采用矩估计是比较不精确的，所以这里不予讨论
1.1.3 ARMA （1，1）
矩估计法参数估计的思路：
方差：
1.2 最小二乘估计
1.2.1 AR 模型
最小二乘参数估计的思路：
对于AR （P ）而言也可以得到类似矩估计得到的方程，即最小二乘与矩估计得到的估计量相同。
1.2.2 MA 模型
最小二乘参数估计的思路：
1.2.3 ARMA 模型
最小二乘参数估计的思路：
1.3 极大似然估计与无条件最小二乘估计
2 R 中如何实现时间序列参数估计
2.1 对于AR 模型
ar(x, aic = TRUE, order.max = NULL,
method=c("yule-walker", "burg", "ols", "mle", "yw"), na.action, series, ...)
> ar(ar1.s,order.max=1,AIC=F,method='yw')#即矩估计 Call:
ar(x = ar1.s, order.max = 1, method = "yw", AIC = F) Coefficients: 1
0.8314
Order selected 1 sigma^2 estimated as 1.382
> ar(ar1.s,order.max=1,AIC=F,method='ols')#最小二乘估计 Call:
ar(x = ar1.s, order.max = 1, method = "ols", AIC = F) Coefficients: 1
0.857
Intercept: 0.02499 (0.1308)
Order selected 1 sigma^2 estimated as 1.008
> ar(ar1.s,order.max=1,AIC=F,method='mle')#极大似然估计 Call:
ar(x = ar1.s, order.max = 1, method = "mle", AIC = F) Coefficients: 1
0.8924
Order selected 1 sigma^2 estimated as 1.041
采用自编函数总结三个不同的估计值
> Myar(ar2.s,order.max=3)
最小二乘估计 矩估计 极大似然估计
1 1.5137146 1.4694476 1.5061369
2 -0.8049905 -0.7646034 -0.7964453
2.2 对于ARMA 模型
arima(x, order = c(0, 0, 0), seasonal = list(order = c(0, 0, 0), period = NA),
xreg = NULL, include.mean = TRUE, transform.pars = TRUE, fixed = NULL, init = NULL, method = c("CSS-ML", "ML", "CSS"), n.cond, optim.control = list(), kappa = 1e+06, io = NULL, xtransf, transfer = NULL)
order 的三个参数分别代表AR ，差分 MA的阶数
> arima(arma11.s,order=c(1,0,1),method='CSS') Call: arima(x = arma11.s, order = c(1, 0, 1), method = "CSS") Coefficients: ar1 ma1 intercept
0.5586 0.3669 0.3928 s.e. 0.1219 0.1564 0.3380
sigma^2 estimated as 1.199: part log likelihood = -150.98
> arima(arma11.s,order=c(1,0,1),method='ML') Call: arima(x = arma11.s, order = c(1, 0, 1), method = "ML") Coefficients:
ar1 ma1 intercept
0.5647 0.3557 0.3216
s.e. 0.1205 0.1585 0.3358
sigma^2 estimated as 1.197: log likelihood = -151.33, aic = 308.65
采用自编函数总结三个不同的估计值
> Myarima(arma11.s,order=c(1,0,1)) $coef
条件SS 估计 极大似然估计 条件似然估计
ar1 0.5585828 0.5647477 0.5647498
ma1 0.3668814 0.3556965 0.3556973
intercept 0.3927654 0.3216166 0.3216152
$log
条件SS 估计 极大似然估计 条件似然估计
[1,] -150.984 -151.3268 -151.3268
$sigma2
条件SS 估计 极大似然估计 条件似然估计
[1,] 1.199378 1.196984 1.196984
$aic
条件SS 估计 极大似然估计 条件似然估计
[1,] NA 308.6537 308.6537
2.3 采用自助法arima.boot()
此函数估计的是参数的取值置信区间，而不是指具体的某个值，与arima 是不同的。
> res=arima(sqrt(hare),order=c(3,0,0),include.mean=T) > set.seed(12345)
> # Method I以最初三个观测为条件，并假设误差服从正态分布，得到95%的置信区间quantile 用于计算置信区间值，signif 类似于四舍五入函数，保留有效数值。
>
coefm.cond.norm =arima.boot(res,cond.boot=T,is.normal=T,B=1000,init=sqrt(hare))
> signif(apply(coefm.cond.norm,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3) ar1 ar2 ar3 intercept noise var 2.5% 0.593 -0.667 -0.6740 5.12 0.548
97.5% 1.280 0.244 -0.0135 6.38 1.540 >
> # Method II 假设误差并不服从正态分布，而是需要从样本抽样中得到coefm.cond.replace =arima.boot(res,cond.boot=T,is.normal=F,B=1000,init=sqrt(hare))
>
signif(apply(coefm.cond.replace,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3)
ar1 ar2 ar3 intercept noise var
2.5% 0.611 -0.700 -0.6720 4.98 0.516
97.5% 1.300 0.241 -0.0417 6.32 1.500
> # Method III基于平稳自助法的置信区间，且误差服从正态分布 >
coefm.norm =arima.boot(res,cond.boot=F,is.normal=T,ntrans=100,B=1000,init=sqrt(hare))
>
signif(apply(coefm.norm,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3)
ar1 ar2 ar3 intercept noise var
2.5% 0.687 -0.747 -0.6600 4.99 0.508 97.5% 1.380 0.192 -0.0168 6.33 1.500 >
> # Method IV基于平稳自助法的置信区间，且误差不服从正态分布
coefm.replace =arima.boot(res,cond.boot=F,is.normal=F,ntrans=100,B=1000,init=sqrt(hare))
>
signif(apply(coefm.replace,2,function(x){quantile(x,c(.025,.975),na.rm=T)}),3)
ar1 ar2 ar3 intercept noise var
2.5% 0.70 -0.715 -0.6620 4.98 0.47
97.5% 1.36 0.183 -0.0187 6.30 1.50
3 附自编函数
3.1 Myar
#用于自回归模型的参数估计，整合矩估计，最小二乘估计，以及极大似然估计 #该函数用于对时间序列中心化数据（因此截距项一定为0）估计AR 模型的参数，AIC 为真时，滞后项根据AIC 准则确定，为假时则根据设置的order.max 设定
Myar=function(tsdata, order.max = 1,AIC = F){
library(TSA)
ols
yw
mle
olscoef
ywcoef
mlecoef
result=data.frame(olscoef,ywcoef,mlecoef)
colnames(result)=c('最小二乘估计',' 矩估计',' 极大似然估计') return(result)
}
3.2 Myarima
#用于自回归模型的参数估计，整合矩估计，最小二乘估计，以及极大似然估计 #该函数用于对时间序列中心化数据（因此截距项一定为0）估计AR 模型的参数，AIC 为真时，滞后项根据AIC 准则确定，为假时则根据设置的order.max 设定
Myarima=function(tsdata, order=c(0,0,0)){
library(TSA)
result=NULL
css
ml
cssml
result$coef=cbind(css$coef,ml$coef,cssml$coef)
result$log=cbind(css$log,ml$log,cssml$log)
result$sigma2=cbind(css$sigma2,ml$sigma2,cssml$sigma2)
result$aic=cbind(NA,ml$aic,cssml$aic)
colnames(result$coef)=c('条件SS 估计',' 极大似然估计',' 条件似然估计') colnames(result$log)=c('条件SS 估计',' 极大似然估计',' 条件似然估计') colnames(result$aic)=c('条件SS 估计',' 极大似然估计',' 条件似然估计') colnames(result$sigma2)=c('条件SS 估计',' 极大似然估计',' 条件似然估计')
return(result) }