#########################################
#  第10章 状態空間モデルによる時系列の解析
#########################################
#
# ARMA モデルの推定
#
data(Sunspot)         # Sun spot number data
y <- log10(Sunspot)   # y = log( sunspot data )

#
# ARMA(4,3) model
#
armafit(y, ar.order=4, ma.order=2)

#
# ARMA(3,1)
#
z1 <- armafit(y, ar.order=3, ma.order=1)
z1
#
# ARMA(3,2) model
# Use the parameters of ARMA(3,1) as the initial estimates
# of ARMA(3,2)
#
z1 <- armafit(y, ar.order=3, ma.order=2, ar=c(1.42758,-0.69843,0.00998), ma=c(0.3974,0.0))
z1

#####################################
#
# Estiamte ARMA(2,4) and compute Impuse response, autocovariance,
# PARCOR, power spectrum and characteristic roots.
#
z1 <- armafit(y, ar.order=2, ma.order=4)
z1
len <- length(Sunspot)
z <- armaimp(arcoef=z1$arcoef, macoef=z1$macoef, v=z1$sigma2, n=len, lag=20)


####################################
# 第11章　トレンドの推定
####################################
#
# Polinomial regression model 
#   highest order = 7
#
  data(TemperData)   # Highest Temperature Data of Tokyo
  polreg(TemperData, 7)
#
# plot AIC's of polinomial regression models
#
plot(aic,type="b",lwd=2,ylim=c(2450,2520))

#
#  Trigonometric regression models
# 
 lsqr(TemperData)
#

####################################
#  Trend Model
####################################
#
#  Second order trend model
#
trend(TemperData, 2)

#
#  First order trend model
#     TAU2 = 0.223, Sig2 = 0.001
#
 trend(TemperData, 1, 0.223, 0.001)