Keyword: 多変量時系列, 多入力モデル, 推定
概要
本サンプルは時系列モデルの推定の計算を行うFortranによるサンプルプログラムです。 本サンプルは以下に示される時系列データを分析し、時系列モデルの推定値を出力します。
※本サンプルはnAG Fortranライブラリに含まれるルーチン g13bef() のExampleコードです。本サンプル及びルーチンの詳細情報は g13bef のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで
入力データ
（本ルーチンの詳細はg13bef のマニュアルページを参照）| このデータをダウンロード |
G13BEF Example Program Data 1 1 40 2 3 20 0 0 :: KZEF,KFC,NXXY,NSER,KEFF,NIT,KZSP,KPRIV 1 0 0 0 0 1 4 :: MR 1 0 1 3 :: Transfer fun. for series 1, MT(:,1) 0.0 0.0 2.0 0.5 0.0 :: PARA 8.075 105.0 7.819 119.0 7.366 119.0 8.113 109.0 7.380 117.0 7.134 135.0 7.222 126.0 7.768 112.0 7.386 116.0 6.965 122.0 6.478 115.0 8.105 115.0 8.060 122.0 7.684 138.0 7.580 135.0 7.093 125.0 6.129 115.0 6.026 108.0 6.679 100.0 7.414 96.0 7.112 107.0 7.762 115.0 7.645 123.0 8.639 122.0 7.667 128.0 8.080 136.0 6.678 140.0 6.739 122.0 5.569 102.0 5.049 103.0 5.642 89.0 6.808 77.0 6.636 89.0 8.241 94.0 7.968 104.0 8.044 108.0 7.791 119.0 7.024 126.0 6.102 119.0 6.053 103.0 :: End of XXY
- １行目はタイトル行で読み飛ばされます。
- ２行目に入力時系列の値と出力ノイズ時系列の値がXXYの値を上書きしないかどうか（kzef=1：上書きする）、定数Cが初期値で固定されるかどうか（kfc=1：推定される）、入力・出力時系列の長さ（nxxy=40）、入力・出力時系列の合計数（nser=2）、尤度オプション（keff=3：周辺尤度）、最大反復数（nit=20）、最小化処理でZSPの入力値を使用するかどうか（kzsp=0：デフォルト値を使用）、最適化処理の経過を出力するかどうか（kpriv=0：出力しない）を指定しています。
- ３行目はARIMAモデルの次数ベクトルmr（自己回帰の数、非季節階差の次数、移動平均、季節自己回帰、季節階差の次数、季節移動平均と季節期間）を指定しています。
- ４行目は伝達関数モデルの次数(mt)を指定しています。
- ５行目は多入力モデルの引数(para)を指定しています。
- ６〜４５行目に入力時系列と出力時系列の値(xxy)を指定しています。
出力結果
（本ルーチンの詳細はg13bef のマニュアルページを参照）| この出力例をダウンロード |
G13BEF Example Program Results
The number of iterations carried out is 11
Final values of the parameters and their standard deviations
I PARA(I) SD
1 0.380924 0.166379
2 -0.257786 0.178178
3 8.956084 0.948061
4 0.659641 0.060239
5 -75.435521 33.505341
The correlation matrix is
1 2 3 4 5
1 1.0000 -0.1839 -0.1775 -0.0340 0.1394
2 -0.1839 1.0000 0.0518 0.2547 -0.2860
3 -0.1775 0.0518 1.0000 -0.3070 -0.2926
4 -0.0340 0.2547 -0.3070 1.0000 -0.8185
5 0.1394 -0.2860 -0.2926 -0.8185 1.0000
The residuals and the z and n values are
I RES(I) z(t) n(t)
1 0.397 180.567 -75.567
2 3.086 191.430 -72.430
3 -2.818 196.302 -77.302
4 -9.941 195.460 -86.460
5 -5.061 201.594 -84.594
6 14.053 199.076 -64.076
7 2.624 195.211 -69.211
8 -5.823 193.450 -81.450
9 -2.147 197.179 -81.179
10 -0.216 196.217 -74.217
11 -2.517 191.812 -76.812
12 7.916 184.544 -69.544
13 1.423 194.322 -72.322
14 11.936 200.369 -62.369
15 5.117 200.990 -65.990
16 -5.672 200.468 -75.468
17 -5.681 195.763 -80.763
18 -1.637 184.025 -76.025
19 -1.019 175.360 -75.360
20 -2.623 175.492 -79.492
21 3.283 182.162 -75.162
22 6.896 183.857 -68.857
23 5.395 190.797 -67.797
24 0.875 194.327 -72.327
25 -4.153 205.558 -77.558
26 6.206 204.261 -68.261
27 4.208 207.104 -67.104
28 -2.387 196.423 -74.423
29 -11.803 189.924 -87.924
30 6.435 175.158 -72.158
31 1.342 160.761 -71.761
32 -4.924 156.575 -79.575
33 4.799 164.256 -75.256
34 -0.074 167.783 -73.783
35 -6.023 184.483 -80.483
36 -6.427 193.055 -85.055
37 -2.527 199.390 -80.390
38 2.039 201.302 -75.302
39 0.243 195.695 -76.695
40 -3.166 183.738 -80.738
The state set consists of 6 values
6.0530 183.7384 -5.7855 -0.1645 0.1800 -3.0977
The number of degrees of freedom is 34
- ３行目には実行した反復数が出力されています。
- ９〜１３行目にはパラメータの最終値と標準偏差が出力されています。
- １６〜２１行目には相関行列が出力されています。
- ２６〜６６行目には残差、Z値、ノイズが出力されています。
- ７０行目には状態セットが出力されています。
- ７２行目には自由度が出力されています。
ソースコード
（本ルーチンの詳細はg13bef のマニュアルページを参照）
※本サンプルソースコードは科学技術・統計計算ライブラリである「nAG Fortranライブラリ」のルーチンを呼び出します。
サンプルのコンパイル及び実行方法
| このソースコードをダウンロード |
PROGRAM g13befe
! G13BEF Example Program Text
! Mark 23 Release. nAG Copyright 2011.
! .. Use Statements ..
USE nag_library, ONLY : g13bef, nag_wp, x04abf, x04caf
! .. Implicit None Statement ..
IMPLICIT NONE
! .. Parameters ..
INTEGER, PARAMETER :: iset = 1, nin = 5, nout = 6
! .. Local Scalars ..
REAL (KIND=nag_wp) :: d, s
INTEGER :: dp, i, ifail, imwa, inc, isttf, itc, &
iwa, kef, kfc, kpriv, kzef, kzsp, &
ldcm, ldxxy, mx, nadv, ncd, nce, &
ncf, ncg, ndf, ndv, nis, nit, npara, &
nser, nsttf, nxxy, qp, qx, smx
! .. Local Arrays ..
REAL (KIND=nag_wp), ALLOCATABLE :: cm(:,:), para(:), res(:), sd(:), &
sttf(:), wa(:), xxy(:,:)
REAL (KIND=nag_wp) :: zsp(4)
INTEGER :: mr(7)
INTEGER, ALLOCATABLE :: mt(:,:), mwa(:)
! .. Intrinsic Functions ..
INTRINSIC max, sum
! .. Executable Statements ..
WRITE (nout,*) 'G13BEF Example Program Results'
WRITE (nout,*)
! Skip heading in data file
READ (nin,*)
! Read in problem size
READ (nin,*) kzef, kfc, nxxy, nser, kef, nit, kzsp, kpriv
IF (kzsp/=0) THEN
READ (nin,*) zsp
END IF
! Number of input series
nis = nser - 1
! Set the advisory channel to NOUT for monitoring information
IF (kpriv/=0) THEN
nadv = nout
CALL x04abf(iset,nadv)
END IF
ALLOCATE (mt(4,nser))
! Read in orders
READ (nin,*) mr(1:7)
! Read in transfer function
DO i = 1, nis
READ (nin,*) mt(1:4,i)
END DO
! Calculate NPARA and various other quantities required
! for calculate array sizes
npara = 0
ncg = 0
qx = 0
smx = 0
ncf = nser
inc = 1
DO i = 1, nis
npara = npara + mt(2,i) + mt(3,i)
IF (mt(4,i)>1) THEN
ncg = ncg + sum(mt(1:3,i))
IF (mt(4,i)==3) THEN
mx = max(mt(1,i)+mt(2,i),mt(3,i))
qx = max(qx,mx)
smx = smx + mx
END IF
ELSE IF (mt(4,i)==1 .AND. kef==3) THEN
IF (mt(3,i)>0) THEN
ncf = ncf + 1
END IF
inc = inc + 1
END IF
END DO
npara = npara + mr(1) + mr(3) + mr(4) + mr(6) + nser
! Calculate size of arrays
isttf = mr(4)*mr(7) + mr(2) + mr(5)*mr(7) + mr(3) + &
max(mr(1),mr(6)*mr(7)) + ncg
ldxxy = nxxy
ldcm = npara
qp = mr(3) + mr(6)*mr(7)
dp = mr(2) + mr(5)*mr(7)
IF (mr(3)>0 .AND. kef>1) THEN
inc = inc + 1
END IF
IF (kfc>0 .AND. kef==3) THEN
inc = inc + 1
END IF
qx = qp
ncd = npara + kfc + smx
IF (mr(1)>0) THEN
ncf = ncf + inc
END IF
IF (mr(3)>0) THEN
ncf = ncf + inc
END IF
IF (mr(4)>0) THEN
ncf = ncf + inc
END IF
IF (mr(6)>0) THEN
ncf = ncf + inc
END IF
IF (qx>0) THEN
ncf = ncf + 1
END IF
IF (kfc>0) THEN
ncf = ncf + 1
END IF
ncd = ncd + qx
nce = nxxy + dp + 6*qx
iwa = 2*ncd**2 + nce*(ncf+4)
iwa = 2*iwa
imwa = 16*nser + 7*ncd + 3*npara + 3*kfc + 27
ALLOCATE (xxy(ldxxy,nser),para(npara),sd(npara),cm(ldcm,npara), &
res(nxxy),sttf(isttf),wa(iwa),mwa(imwa))
! Read in rest of data
READ (nin,*) para(1:npara)
READ (nin,*) (xxy(i,1:nser),i=1,nxxy)
ifail = -1
CALL g13bef(mr,nser,mt,para,npara,kfc,nxxy,xxy,ldxxy,kef,nit,kzsp,zsp, &
itc,sd,cm,ldcm,s,d,ndf,kzef,res,sttf,isttf,nsttf,wa,iwa,mwa,imwa, &
kpriv,ifail)
IF (ifail/=0) THEN
IF (ifail/=8 .AND. ifail/=9 .AND. ifail/=11) THEN
GO TO 20
END IF
END IF
! Display results
WRITE (nout,99999) 'The number of iterations carried out is', itc
WRITE (nout,*)
WRITE (nout,*) &
'Final values of the parameters and their standard deviations'
WRITE (nout,*)
WRITE (nout,*) ' I PARA(I) SD'
WRITE (nout,*)
WRITE (nout,99998) (i,para(i),sd(i),i=1,npara)
WRITE (nout,*)
FLUSH (nout)
ifail = 0
CALL x04caf('General',' ',npara,npara,cm,ldcm, &
'The correlation matrix is',ifail)
WRITE (nout,*)
WRITE (nout,*) 'The residuals and the z and n values are'
WRITE (nout,*)
WRITE (nout,*) ' I RES(I) z(t) n(t)'
WRITE (nout,*)
ndv = nxxy - mr(2) - mr(5)*mr(7)
DO i = 1, nxxy
IF (i<=ndv) THEN
WRITE (nout,99997) i, res(i), xxy(i,1:nser)
ELSE
WRITE (nout,99996) i, xxy(i,1:nser)
END IF
END DO
IF (mr(2)/=0 .OR. mr(5)/=0) THEN
WRITE (nout,*)
WRITE (nout,*) &
'** Note that the residuals relate to differenced values **'
END IF
WRITE (nout,*)
WRITE (nout,99995) 'The state set consists of', nsttf, ' values'
WRITE (nout,*)
WRITE (nout,99994) sttf(1:nsttf)
WRITE (nout,*)
WRITE (nout,99999) 'The number of degrees of freedom is', ndf
20 CONTINUE
99999 FORMAT (1X,A,I4)
99998 FORMAT (1X,I4,2F20.6)
99997 FORMAT (1X,I4,3F15.3)
99996 FORMAT (1X,I4,F30.3,F15.3)
99995 FORMAT (1X,A,I4,A)
99994 FORMAT (1X,6F10.4)
END PROGRAM g13befe