Keyword: VARMA, 多変量, 時系列, 実現値
概要
本サンプルはVARMAモデルの多変量時系列の実現値の生成を行うFortranによるサンプルプログラムです。 本サンプルは以下の式で示されるVARMAモデルに対して１つの実現値が48個の観測値をもつ２変量時系列の２つの実現値を生成し出力します。
※本サンプルはnAG Fortranライブラリに含まれるルーチン g05pjf() のExampleコードです。本サンプル及びルーチンの詳細情報は g05pjf のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで
入力データ
（本ルーチンの詳細はg05pjf のマニュアルページを参照）| このデータをダウンロード |
G05PJF Example Program Data 1 1 1762543 :: GENID,SUBID,SEED(1) 48 2 :: N,NREAL 2 1 0 :: K, IP, IQ 0.80 0.07 0.00 0.58 :: End of PHI 5.00 9.00 :: XMEAN 2.97 0.64 5.38 :: End of VAR (lower triangle)
- １行目はタイトル行で読み飛ばされます。
- ２行目に使用する生成器（genid=1:nAG基本生成器）、生成器に関する情報（subid=1:GENIDが1の場合この値は参照されません）、生成器の初期値（seed=1762543）を指定しています。
- ３行目に生成される観測値の数（n=48）と実現値の数（nreal=2）を指定しています。
- ３行目に多変量時系列の次数（k=2）、自己回帰パラメータ行列の数（ip=1）、移動平均パラメータ行列の数（iq=0）を指定しています。
- ５〜６行目にモデルの自己回帰パラメータ行列(phi)を指定しています。
- ７行目に多変量時系列の平均ベクトル(xmean)を指定しています。
- ８〜９行目に分散共分散行列の下三角部分の要素(var)を指定しています。
出力結果
（本ルーチンの詳細はg05pjf のマニュアルページを参照）| この出力例をダウンロード |
G05PJF Example Program Results
Realization Number 1
Series number 1
-------------
4.833 2.813 3.224 3.825 1.023 1.415 2.184 3.005
5.547 4.832 4.705 5.484 9.407 10.335 8.495 7.478
6.373 6.692 6.698 6.976 6.200 4.458 2.520 3.517
3.054 5.439 5.699 7.136 5.750 8.497 9.563 11.604
9.020 10.063 7.976 5.927 4.992 4.222 3.982 7.107
3.554 7.045 7.025 4.106 5.106 5.954 8.026 7.212
Series number 2
-------------
8.458 9.140 10.866 10.975 9.245 5.054 5.023 12.486
10.534 10.590 11.376 8.793 14.445 13.237 11.030 8.405
7.187 8.291 5.920 9.390 10.055 6.222 7.751 10.604
12.441 10.664 10.960 8.022 10.073 12.870 12.665 14.064
11.867 12.894 10.546 12.754 8.594 9.042 12.029 12.557
9.746 5.487 5.500 8.629 9.723 8.632 6.383 12.484
Realization Number 2
Series number 1
-------------
5.396 4.811 2.685 5.824 2.449 3.563 5.663 6.209
3.130 4.308 4.333 4.903 1.770 1.278 1.340 -0.527
1.745 3.211 4.478 5.170 5.365 4.852 6.080 6.464
2.765 2.148 6.641 7.224 10.316 7.102 5.604 3.934
4.839 3.698 5.210 5.384 7.652 7.315 7.332 7.561
7.537 7.788 6.868 7.575 6.108 6.188 8.132 10.310
Series number 2
-------------
11.345 10.070 13.654 12.409 11.329 13.054 12.465 9.867
10.263 13.394 10.553 10.331 7.814 8.747 10.025 11.167
10.626 9.366 9.607 9.662 10.492 10.766 11.512 10.813
10.799 8.780 9.221 14.245 11.575 10.620 8.282 5.447
9.935 9.386 11.627 10.066 11.394 7.951 7.907 12.616
15.246 9.962 13.216 11.350 11.227 6.021 6.968 12.428
- ８〜１３行目に生成された１つ目の実現値の時系列１が出力されています。
- １８〜２３行目に生成された１つ目の実現値の時系列２が出力されています。
- ３０〜３５行目に生成された２つ目の実現値の時系列１が出力されています。
- ４０〜４５行目に生成された２つ目の実現値の時系列２が出力されています。
ソースコード
（本ルーチンの詳細はg05pjf のマニュアルページを参照）
※本サンプルソースコードは科学技術・統計計算ライブラリである「nAG Fortranライブラリ」のルーチンを呼び出します。
サンプルのコンパイル及び実行方法
| このソースコードをダウンロード |
PROGRAM g05pjfe
! G05PJF Example Program Text
! Mark 23 Release. nAG Copyright 2011.
! .. Use Statements ..
USE nag_library, ONLY : g05kff, g05pjf, nag_wp
! .. Implicit None Statement ..
IMPLICIT NONE
! .. Parameters ..
INTEGER, PARAMETER :: lseed = 1, nin = 5, nout = 6
! .. Local Scalars ..
INTEGER :: genid, i, ifail, ii, ip, iq, j, k, &
k2, l, ldvar, ldx, lphi, lr, lstate, &
ltheta, mode, n, nreal, rn, subid
! .. Local Arrays ..
REAL (KIND=nag_wp), ALLOCATABLE :: phi(:), r(:), theta(:), var(:,:), &
x(:,:), xmean(:)
INTEGER :: seed(lseed)
INTEGER, ALLOCATABLE :: state(:)
! .. Intrinsic Functions ..
INTRINSIC max
! .. Executable Statements ..
WRITE (nout,*) 'G05PJF Example Program Results'
WRITE (nout,*)
! Skip heading in data file
READ (nin,*)
! Read in the base generator information and seed
READ (nin,*) genid, subid, seed(1)
! Initial call to initialiser to get size of STATE array
lstate = 0
ALLOCATE (state(lstate))
ifail = 0
CALL g05kff(genid,subid,seed,lseed,state,lstate,ifail)
! Reallocate STATE
DEALLOCATE (state)
ALLOCATE (state(lstate))
! Initialize the generator to a repeatable sequence
ifail = 0
CALL g05kff(genid,subid,seed,lseed,state,lstate,ifail)
! Read in the sample size and number of realizations
READ (nin,*) n, nreal
! Read in the number of coefficients
READ (nin,*) k, ip, iq
k2 = k**2
rn = max(ip,iq)
l = k*(k+1)/2
IF (ip>0) THEN
l = l + (ip-1)*k2
END IF
IF (k>=6) THEN
lr = (5*rn**2+1)*k2 + (4*rn+3) + 4
ELSE
lr = ((ip+iq)**2+1)*k2 + (4*(ip+iq)+3)*k + &
max(k*rn*(k*rn+2),k2*(ip+iq)**2+l*(l+3)+k2*(iq+1)) + 4
END IF
lphi = ip*k*k
ltheta = iq*k*k
ldvar = k
ldx = k
ALLOCATE (phi(lphi),theta(ltheta),var(ldvar,k),r(lr),x(ldx,n),xmean(k))
! Read in the AR parameters
DO l = 1, ip
DO i = 1, k
ii = (l-1)*k*k + i
READ (nin,*) (phi(ii+k*(j-1)),j=1,k)
END DO
END DO
! Read in the MA parameters
DO l = 1, iq
DO i = 1, k
ii = (l-1)*k*k + i
READ (nin,*) (theta(ii+k*(j-1)),j=1,k)
END DO
END DO
! Read in the means
READ (nin,*) xmean(1:k)
! Read in the variance / covariance matrix
READ (nin,*) (var(i,1:i),i=1,k)
! For the first realization we need to set up the reference vector
! as well as generate the series
mode = 2
! Generate NREAL realizations
D_LP: DO rn = 1, nreal
ifail = 0
CALL g05pjf(mode,n,k,xmean,ip,phi,iq,theta,var,ldvar,r,lr,state,x, &
ldx,ifail)
! Display the results
WRITE (nout,99999) ' Realization Number ', rn
DO i = 1, k
WRITE (nout,*)
WRITE (nout,99999) ' Series number ', i
WRITE (nout,*) ' -------------'
WRITE (nout,*)
WRITE (nout,99998) x(i,1:n)
END DO
WRITE (nout,*)
! For subsequent realizations we use previous reference vector
mode = 3
END DO D_LP
99999 FORMAT (1X,A,I0)
99998 FORMAT (8(2X,F8.3))
END PROGRAM g05pjfe