Ridgeパラメータを最適化するRidge回帰

ホーム > 製品 > nAG数値計算ライブラリ > nAG Fortranライブラリ > サンプルソースコード集 > Ridgeパラメータを最適化するRidge回帰

Keyword: Ridge, リッジ, パラメータ, 回帰

概要

本サンプルはRidgeパラメータを最適化するRidge回帰の計算を行うFortranによるサンプルプログラムです。本サンプルは以下に示されるデータについてRidge回帰の計算を行います。

Ridge回帰のデータ　

※本サンプルはnAG Fortranライブラリに含まれるルーチン g02kaf() のExampleコードです。本サンプル及びルーチンの詳細情報は g02kaf のマニュアルページをご参照ください。
ご相談やお問い合わせはこちらまで

入力データ

このデータをダウンロード

G02KAF Example Program Data
20 3 0.5 1 1.0e-4 25 2 2 0.0 : N, M, H, OPT, TOL, NITER, ORIG, OPTLOO, TAU
19.5 43.1 29.1 11.9
24.7 49.8 28.2 22.8
30.7 51.9 37.0 18.7
29.8 54.3 31.1 20.1
19.1 42.2 30.9 12.9
25.6 53.9 23.7 21.7
31.4 58.5 27.6 27.1
27.9 52.1 30.6 25.4
22.1 49.9 23.2 21.3
25.5 53.5 24.8 19.3
31.1 56.6 30.0 25.4
30.4 56.7 28.3 27.2
18.7 46.5 23.0 11.7
19.7 44.2 28.6 17.8
14.6 42.7 21.3 12.8
29.5 54.4 30.1 23.9
27.7 55.3 25.7 22.6
30.2 58.6 24.6 25.4
22.7 48.2 27.1 14.8
25.2 51.0 27.5 21.1 : End of data
1 1 1 : ISX

１行目はタイトル行で読み飛ばされます。
２行目は観測値の数（n=20）、データ行列で有効な独立変数の数（m=3）、Ridge回帰パラメータの初期値（h=0.5）、Ridge回帰パラメータの最適化に使用される予測誤差の測定法（opt=1:GCV）、回帰パラメータの反復が中止となる値（tol=1.0e-4）、Ridg回帰パラメータの最適化が可能な反復の最大値（niter=25）、パラメータ推定値が元データに対して計算されるか標準化されたデータに対して計算されるか（orig=2:標準化されたデータに対して計算）、予測誤差の交差検定（leave-one-out cross-validation）の推定値が計算されるかどうか（optloo=2:交差検定の推定値が計算される）、データ行列の特異値分解でこれより小さい特異値をゼロにするその基準値（tau=0.0）を指定しています。
３～２２行目は独立変数の観測値(x)と従属変数の観測値(y)を指定しています。
２３行目はどの独立変数がモデル(isx)に含まれるかを指定しています。

出力結果

（本ルーチンの詳細はg02kaf のマニュアルページを参照）

この出力例をダウンロード

 G02KAF Example Program Results

 Value of ridge parameter:     0.0712

 Sum of squares of residuals: 0.1092E+03
 Degrees of freedom:     16
 Number of effective parameters:     2.9059

 Parameter estimates
    1     20.1950
    2      9.7934
    3      9.9576
    4     -2.0125

 Number of iterations:               6

 Ridge parameter minimises GCV

 Estimated prediction errors:
 GCV    =     7.4718
 UEV    =     6.3862
 FPE    =     7.3141
 BIC    =     8.2380
 LOO CV =     7.5495

 Residuals
    1     -1.9894
    2      3.5469
    3     -3.0392
    4     -3.0309
    5     -0.1899
    6     -0.3146
    7      0.9775
    8      4.0157
    9      2.5332
   10     -2.3560
   11      0.5446
   12      2.3989
   13     -4.0876
   14      3.2778
   15      0.2894
   16      0.7330
   17     -0.7116
   18     -0.6092
   19     -2.9995
   20      1.0110

 Variance inflation factors
    1      0.2928
    2      0.4162
    3      0.8089

３行目にRidgeパラメータの値が出力されています。
５行目に残差平方和が出力されています。
６行目に自由度が出力されています。
７行目に効果パラメータの数が出力されています。
１０～１３行目にパラメータ推定値が出力されています。
１５行目に反復数が出力されています。
１７行目にRidge回帰パラメータの最小化の手法GCVが出力されています。
２０～２４行目に推定予測誤差が出力されています。"GCV" は "Generalised cross-validation" を表します。"UEV" は "Unbiased estimate of variance" を表します。"FPE" は "Future prediction error" を表します。"BIC" は "Bayesian information criterion" を表します。"LOO CV" は "Leave-one-out cross-validation" を表します。
２７～４６行目に残差が出力されています。
４９～５１行目に分散拡大要因が出力されています。

ソースコード

（本ルーチンの詳細はg02kaf のマニュアルページを参照）

※本サンプルソースコードは科学技術・統計計算ライブラリである「nAG Fortranライブラリ」のルーチンを呼び出します。
サンプルのコンパイル及び実行方法

このソースコードをダウンロード

    PROGRAM g02kafe

!      G02KAF Example Program Text

!      Mark 23 Release. nAG Copyright 2011.

!      .. Use Statements ..
       USE nag_library, ONLY : g02kaf, nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       INTEGER, PARAMETER              :: nin = 5, nout = 6
!      .. Local Scalars ..
       REAL (KIND=nag_wp)              :: h, nep, rss, tau, tol
       INTEGER                         :: df, i, ifail, ip, ldx, m, n, niter,  &
                                          opt, optloo, orig
!      .. Local Arrays ..
       REAL (KIND=nag_wp), ALLOCATABLE :: b(:), res(:), vif(:), x(:,:), y(:)
       REAL (KIND=nag_wp)              :: perr(5)
       INTEGER, ALLOCATABLE            :: isx(:)
!      .. Intrinsic Functions ..
       INTRINSIC                          count
!      .. Executable Statements ..
       WRITE (nout,*) 'G02KAF Example Program Results'
       WRITE (nout,*)

!      Skip heading in data file
       READ (nin,*)

!      Read in the problem size
       READ (nin,*) n, m, h, opt, tol, niter, orig, optloo, tau

       ldx = n
       ALLOCATE (x(ldx,m),y(n),isx(m))

!      Read in data
       READ (nin,*) (x(i,1:m),y(i),i=1,n)

!      Read in variable inclusion flags
       READ (nin,*) isx(1:m)

!      Calculate IP
       ip = count(isx(1:m)==1)

       ALLOCATE (b(ip+1),vif(ip),res(n))

!      Fit ridge regression model
       ifail = -1
       CALL g02kaf(n,m,x,ldx,isx,ip,tau,y,h,opt,niter,tol,nep,orig,b,vif,res, &
          rss,df,optloo,perr,ifail)
       IF (ifail/=0) THEN
          IF (ifail/=-1) THEN
             GO TO 20
          END IF
       END IF

!      Display results
       WRITE (nout,99999) 'Value of ridge parameter:', h
       WRITE (nout,*)
       WRITE (nout,99998) 'Sum of squares of residuals:', rss
       WRITE (nout,99997) 'Degrees of freedom: ', df
       WRITE (nout,99999) 'Number of effective parameters:', nep
       WRITE (nout,*)
       WRITE (nout,*) 'Parameter estimates'
       WRITE (nout,99995) (i,b(i),i=1,ip+1)
       WRITE (nout,*)
       WRITE (nout,99996) 'Number of iterations:', niter
       WRITE (nout,*)
       IF (opt==1) THEN
          WRITE (nout,*) 'Ridge parameter minimises GCV'
       ELSE IF (opt==2) THEN
          WRITE (nout,*) 'Ridge parameter minimises UEV'
       ELSE IF (opt==3) THEN
          WRITE (nout,*) 'Ridge parameter minimises FPE'
       ELSE IF (opt==4) THEN
          WRITE (nout,*) 'Ridge parameter minimises BIC'
       END IF
       WRITE (nout,*)
       WRITE (nout,*) 'Estimated prediction errors:'
       WRITE (nout,99999) 'GCV    =', perr(1)
       WRITE (nout,99999) 'UEV    =', perr(2)
       WRITE (nout,99999) 'FPE    =', perr(3)
       WRITE (nout,99999) 'BIC    =', perr(4)
       IF (optloo==2) THEN
          WRITE (nout,99999) 'LOO CV =', perr(5)
       END IF
       WRITE (nout,*)
       WRITE (nout,*) 'Residuals'
       WRITE (nout,99995) (i,res(i),i=1,n)
       WRITE (nout,*)
       WRITE (nout,*) 'Variance inflation factors'
       WRITE (nout,99995) (i,vif(i),i=1,ip)

20     CONTINUE

99999  FORMAT (1X,A,1X,F10.4)
99998  FORMAT (1X,A,E11.4)
99997  FORMAT (1X,A,1X,I5)
99996  FORMAT (1X,A,I16)
99995  FORMAT (1X,I4,1X,F11.4)
    END PROGRAM g02kafe

Fortranによるサンプルソースコード : 使用ルーチン名：g02kaf

概要

入力データ

出力結果

ソースコード