Program g02kbfe

!     G02KBF Example Program Text

!     Mark 25 Release. NAG Copyright 2014.

!     .. Use Statements ..
      Use nag_library, Only: g02kbf, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Integer                          :: i, ifail, ip, ldb, ldpe, ldvf, ldx,  &
                                          lh, lpec, m, n, pl, tdb, tdpe, tdvf, &
                                          wantb, wantvf
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: b(:,:), h(:), nep(:), pe(:,:),       &
                                          vf(:,:), x(:,:), y(:)
      Integer, Allocatable             :: isx(:)
      Character (1), Allocatable       :: pec(:)
!     .. Intrinsic Procedures ..
      Intrinsic                        :: count, min
!     .. Executable Statements ..
      Write (nout,*) 'G02KBF Example Program Results'
      Write (nout,*)

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

!     Read in the problem size 
      Read (nin,*) n, m, lh, lpec, wantb, wantvf

      ldx = n
      Allocate (x(ldx,m),isx(m),y(n),h(lh),pec(lpec))

!     Read in data
      If (lpec>0) Then
        Read (nin,*) pec(1:lpec)
      End If
      Read (nin,*)(x(i,1:m),y(i),i=1,n)

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

!     Read in the ridge coefficients
      Read (nin,*) h(1:lh)

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

      If (wantb/=0) Then
        ldb = ip + 1
        tdb = lh
      Else
        ldb = 0
        tdb = 0
      End If
      If (wantvf/=0) Then
        ldvf = ip
        tdvf = lh
      Else
        ldvf = 0
        tdvf = 0
      End If
      If (lpec>0) Then
        ldpe = lpec
        tdpe = lh
      Else
        ldpe = 0
        tdpe = 0
      End If
      Allocate (nep(lh),b(ldb,tdb),vf(ldvf,tdvf),pe(ldpe,tdpe))

!     Fit ridge regression
      ifail = 0
      Call g02kbf(n,m,x,ldx,isx,ip,y,lh,h,nep,wantb,b,ldb,wantvf,vf,ldvf,lpec, &
        pec,pe,ldpe,ifail)

!     Display results
      Write (nout,99994) 'Number of parameters used = ', ip + 1
      Write (nout,*) 'Effective number of parameters (NEP):'
      Write (nout,*) '   Ridge   '
      Write (nout,*) '   Coeff.  ', 'NEP'
      Write (nout,99993)(h(i),nep(i),i=1,lh)

!     Parameter estimates
      If (wantb/=0) Then
        Write (nout,*)
        If (wantb==1) Then
          Write (nout,*) 'Parameter Estimates (Original scalings)'
        Else
          Write (nout,*) 'Parameter Estimates (Standarised)'
        End If
        pl = min(ip,4)
        Write (nout,*) '  Ridge  '
        Write (nout,99997) '   Coeff.  ', ' Intercept ', (i,i=1,pl)
        If (pl<ip-1) Then
          Write (nout,99996)(i,i=pl+1,ip-1)
        End If
        pl = min(ip+1,5)
        Do i = 1, lh
          Write (nout,99999) h(i), b(1:pl,i)
          If (pl<ip) Then
            Write (nout,99998) b((pl+1):ip,i)
          End If
        End Do
      End If

!     Variance inflation factors
      If (wantvf/=0) Then
        Write (nout,*)
        Write (nout,*) 'Variance Inflation Factors'
        pl = min(ip,5)
        Write (nout,*) '  Ridge  '
        Write (nout,99995) '  Coeff.  ', (i,i=1,pl)
        If (pl<ip) Then
          Write (nout,99996)(i,i=pl+1,ip)
        End If
        Do i = 1, lh
          Write (nout,99999) h(i), vf(1:pl,i)
          If (pl<ip) Then
            Write (nout,99998) vf((pl+1):ip,i)
          End If
        End Do
      End If

!     Prediction error criterion
      If (lpec>0) Then
        Write (nout,*)
        Write (nout,*) 'Prediction error criterion'
        pl = min(lpec,5)
        Write (nout,*) '  Ridge  '
        Write (nout,99995) '  Coeff.  ', (i,i=1,pl)
        If (pl<lpec) Then
          Write (nout,99996)(i,i=pl+1,lpec)
        End If
        Do i = 1, lh
          Write (nout,99999) h(i), pe(1:pl,i)
          If (pl<ip) Then
            Write (nout,99998) pe((pl+1):ip,i)
          End If
        End Do
        Write (nout,*)
        Write (nout,*) 'Key:'
        Do i = 1, lpec
          Select Case (pec(i))
          Case ('L')
            Write (nout,99992) i, 'Leave one out cross-validation'
          Case ('G')
            Write (nout,99992) i, 'Generalised cross-validation'
          Case ('U')
            Write (nout,99992) i, 'Unbiased estimate of variance'
          Case ('F')
            Write (nout,99992) i, 'Final prediction error'
          Case ('B')
            Write (nout,99992) i, 'Bayesian information criterion'
          End Select
        End Do
      End If

99999 Format (1X,F10.4,5F10.4)
99998 Format (1X,10X,5F10.4)
99997 Format (1X,A,A,4I10)
99996 Format (10X,5I10)
99995 Format (1X,A,5I10)
99994 Format (1X,A,I10)
99993 Format (1X,F10.4,F10.4)
99992 Format (1X,1X,I5,1X,A)
    End Program g02kbfe