R言語と時系列学習ノート(2)


ARMAモデルのパラメータ推定方法
ARMAパラメータ推定は,前述した点推定と類似しており,モーメント推定と最小二乗推定の2つの方法も紹介した.
前回のポイント推定と同じように、今回私が共有した内容は主に:モーメント推定、最小二乗推定、1つの応用例題
モーメント推定と最小二乗推定の基本思想については、前の点推定の紹介を参照する.
A RMAモデル(Auto-regressive and Moving Average Model)は時間系列を研究する重要な方法であり、自己回帰モデル(ARモデルと略称する)とスライド平均モデル(MAモデルと略称する)を基礎とする「混合」構成.市場研究では長期追跡資料の研究によく用いられ、例えばパネル研究では消費行為モデルの変遷研究に用いられる.小売研究では,季節変動の特徴を持つ販売量,市場規模の予測などに用いられる.
ARMAモデルの基本原理:予測指標が時間とともに形成されるデータ系列をランダム系列と見なし,このランダム変数のグループが持つ依存関係は元のデータの時間的継続性を体現している.一方、影響要因の影響、一方、自身の変動法則があり、影響要因をx 1,x 2,...,xkと仮定し、回帰分析から、
armaモデルのモーメント推定と最小二乗推定にはyule‐walker方程式を理解する必要がある.以下を参照してください.http://blog.sina.com.cn/s/blog_4b700c4c0102e728.htmlyule-walker方程式についての内容を理解します.
ここでは、データからar,arma,モデルを直接推定する2つの関数:ar(),arima().ここでar()、arima()は拡張パッケージTSAの関数であり、呼び出す場合はこの拡張パッケージをインストールしてロードする必要があります.
まずar()の使い方を見てみましょう.
ar(x, aic = TRUE, order.max = NULL,
  method=c("yule-walker", "burg", "ols","mle", "yw"),
  na.action, series, ...)
ここでarモデルをモーメント推定する場合,methodは「yw」を選択し,極大尤度推定は「mle」を用いる.
例を見てみましょう
やはり、前のブログで述べたarモデルデータを生成する方法でデータを生成し、モーメント推定を行います.
>w<-rnorm(550)
>x<-filter(w,filter=c(1,-0.9),"recursive")
>ar(x,order=2,method="yw")
 
Call:
ar(x = x,order.max = 2, method = "yw")
 
Coefficients:
      1       2 
 0.9870 -0.8905 
 
Orderselected 2  sigma^2 estimated as  1.085
コントラスト係数は,この推定が非常に近いことが分かった.これはarモデルのモーメント推定に利用するyule‐waler方程式が線形方程式と見なすことができ,推定効果が非常に良好であるからである.最小二乗推定値を比較できます.
> ar(x,order=2,method="ols")
 
Call:
ar(x = x, order.max = 2, method ="ols")
 
Coefficients:
     1        2 
 0.9953 -0.8988 
 
Intercept: -0.003315 (0.04332)
 
Order selected 2  sigma^2 estimated as  1.028
推定効果の差は多くないことが分かるが,同様に,極大尤度推定を選択することができ,ここでは後述しない.
 
maモデルでは,最小二乗,極大尤度推定,モーメント推定にかかわらず,効果は望ましくないことを指摘しなければならない.Rではarma(0,q)にmaモデルを表すことができるほか(arma,arimaが提供する方法は比較的細かく、モーメント推定のような粗い方法は採用されていない)、私もmaモデルのモーメント推定関数を見つけていません(効果が悪いためかもしれません.自分で書くしかありません).
armaモデルでは、関数arima()、arma()を用いることができる.パラメータ推定をします.ここでは、2つの関数の使い方を以下に抜粋します.
arma(x, order = c(1, 1), lag = NULL, coef = NULL,
     include.intercept = TRUE, series = NULL, qr.tol = 1e-07, ...)
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.method = "BFGS",
      optim.control = list(), kappa = 1e6)

最後にarmaのパラメータ推定を例に挙げて説明する.
やはり太陽黒子のデータで問題を説明します.
>data(sunspots)
>arma(sunspots,order=c(2,1))
 
Call:
arma(x =sunspots, order = c(2, 1))
 
Coefficient(s):
      ar1       ar2        ma1  intercept 
   1.1984   -0.2116    -0.6214     0.6673