    PROGRAM g05hnfe

!      G05HNF Example Program Text

!      Mark 23 Release. NAG Copyright 2011.

!      .. Use Statements ..
       USE nag_library, ONLY : g05hnf, g05kbf, g13fgf, g13fhf, nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       REAL (KIND=nag_wp), PARAMETER   :: zero = 0.0E0_nag_wp
       INTEGER, PARAMETER              :: nout = 6, nparmx = 10, nregmx = 10,  &
                                          num = 1500, num1 = 3000
       INTEGER, PARAMETER              :: ldcovr = nparmx
!      .. Local Scalars ..
       REAL (KIND=nag_wp)              :: df, fac1, hp, lgf, mean, tol, xterm
       INTEGER                         :: i, ifail, igen, ip, iq, k, ldx,      &
                                          lwork, maxit, mn, npar, npar2, nreg, &
                                          nt
       LOGICAL                         :: copt, fcall
       CHARACTER (1)                   :: dist
!      .. Local Arrays ..
       REAL (KIND=nag_wp)              :: bx(10), covr(ldcovr,nparmx),         &
                                          cvar(100), etm(num1), ht(num1+10),   &
                                          htm(num1), param(nparmx), rvec(40),  &
                                          rwsav(9), sc(nparmx), se(nparmx),    &
                                          theta(nparmx),                       &
                                      work(num1*3+nparmx+nregmx*num1+20*20+1), &
                                          x(num1,10), yt(num1+10)
       INTEGER                         :: iseed(4)
!      .. Intrinsic Functions ..
       INTRINSIC                          real, sin
!      .. Executable Statements ..
       WRITE (nout,*) 'G05HNF Example Program Results'
       WRITE (nout,*)

!      Initialize the generator
       iseed(1) = 111
       igen = 0
       CALL g05kbf(igen,iseed)

!      Example using Normal errors
       ldx = num1
       bx(1) = 1.5E0_nag_wp
       bx(2) = 2.5E0_nag_wp
       mean = 3.0E0_nag_wp

       DO i = 1, num
          fac1 = real(i,kind=nag_wp)*0.01E0_nag_wp
          x(i,1) = 0.01E0_nag_wp + 0.7E0_nag_wp*sin(fac1)
          x(i,2) = 0.5E0_nag_wp + fac1*0.1E0_nag_wp
          x(i,3) = 1.0E0_nag_wp
       END DO

       nreg = 2
       mn = 1
       ip = 1
       iq = 1
       npar = ip + 2*iq + 1

       param(1) = 0.1E0_nag_wp
       param(2) = -0.3E0_nag_wp
       param(3) = 0.1E0_nag_wp
       param(4) = 0.9E0_nag_wp

       df = 5.0E0_nag_wp
       dist = 'N'

!      Call the routine for the first time
       fcall = .TRUE.
       ifail = 0
       CALL g05hnf(dist,800,ip,iq,param,df,ht,yt,fcall,rvec,igen,iseed,rwsav, &
          ifail)
!      Discard the first 800 values, and call the routine again
       fcall = .FALSE.
       ifail = 0
       CALL g05hnf(dist,num,ip,iq,param,df,ht,yt,fcall,rvec,igen,iseed,rwsav, &
          ifail)

       DO i = 1, num
          xterm = zero
          DO k = 1, nreg
             xterm = xterm + x(i,k)*bx(k)
          END DO
          IF (mn==1) THEN
             yt(i) = mean + xterm + yt(i)
          ELSE
             yt(i) = xterm + yt(i)
          END IF
       END DO

       copt = .TRUE.

       maxit = 50
       tol = 1.0E-5_nag_wp

       DO i = 1, npar
          theta(i) = param(i)*0.5E0_nag_wp
       END DO

       IF (mn==1) THEN
          theta(npar+mn) = mean*0.5E0_nag_wp
       END IF

       DO i = 1, nreg
          theta(npar+mn+i) = bx(i)*0.5E0_nag_wp
       END DO

!      Fit an EGARCH model to the generated data
       lwork = nreg*num + 3*num + 3
       npar2 = 1 + iq*2 + ip + mn + nreg
       ifail = 0
       CALL g13fgf(dist,yt,x,ldx,num,ip,iq,nreg,mn,npar2,theta,se,sc,covr, &
          ldcovr,hp,etm,htm,lgf,copt,maxit,tol,work,lwork,ifail)
!      Display the output
       WRITE (nout,*)
       WRITE (nout,*) 'Normal distribution'
       WRITE (nout,*)
       WRITE (nout,*) '             Correct    Parameter   Standard'
       WRITE (nout,*) '             values     estimates    errors'
       DO i = 1, npar
          WRITE (nout,99998) 'PARAM(', i, ')', param(i), theta(i), se(i)
       END DO
       IF (mn==1) THEN
          WRITE (nout,99999) 'Mean    ', mean, theta(npar+1), se(npar+1)
       END IF
       DO i = 1, nreg
          WRITE (nout,99998) 'BX(', i, ')   ', bx(i), theta(npar+mn+i), &
             se(npar+mn+i)
       END DO

       nt = 4
       ifail = 0
       CALL g13fhf(num,nt,ip,iq,theta,cvar,htm,etm,ifail)
       WRITE (nout,*)
       WRITE (nout,99997) 'Volatility forecast = ', cvar(nt)
       WRITE (nout,*)

!      Example using Students T errors
       nreg = 2
       mn = 1
       ip = 1
       iq = 2
       npar = ip + 2*iq + 1

       param(1) = 0.1E0_nag_wp
       param(2) = -0.3E0_nag_wp
       param(3) = -0.1E0_nag_wp
       param(4) = 0.1E0_nag_wp
       param(5) = 0.3E0_nag_wp
       param(6) = 0.7E0_nag_wp

       dist = 'T'

!      Call the routine the first time
       fcall = .TRUE.
       ifail = 0
       CALL g05hnf(dist,num,ip,iq,param,df,ht,yt,fcall,rvec,igen,iseed,rwsav, &
          ifail)
!      Call the routine the a second time
       fcall = .FALSE.
       ifail = 0
       CALL g05hnf(dist,num,ip,iq,param,df,ht,yt,fcall,rvec,igen,iseed,rwsav, &
          ifail)
!      Discard all previously generated values, and call the routine again
       fcall = .FALSE.
       ifail = 0
       CALL g05hnf(dist,num,ip,iq,param,df,ht,yt,fcall,rvec,igen,iseed,rwsav, &
          ifail)

       DO i = 1, num
          xterm = zero
          DO k = 1, nreg
             xterm = xterm + x(i,k)*bx(k)
          END DO
          IF (mn==1) THEN
             yt(i) = mean + xterm + yt(i)
          ELSE
             yt(i) = xterm + yt(i)
          END IF
       END DO

       copt = .TRUE.
       maxit = 50
       tol = 1.0E-5_nag_wp

       DO i = 1, npar
          theta(i) = param(i)*0.5E0_nag_wp
       END DO

       theta(npar+1) = df*0.65E0_nag_wp

       IF (mn==1) THEN
          theta(npar+1+mn) = mean*0.5E0_nag_wp
       END IF

       DO i = 1, nreg
          theta(npar+1+mn+i) = bx(i)*0.5E0_nag_wp
       END DO

!      Fit an EGARCH model to the generated data
       lwork = nreg*num + 3*num + 3
       npar2 = 2 + iq*2 + ip + mn + nreg
       ifail = 0
       CALL g13fgf(dist,yt,x,ldx,num,ip,iq,nreg,mn,npar2,theta,se,sc,covr, &
          ldcovr,hp,etm,htm,lgf,copt,maxit,tol,work,lwork,ifail)

       WRITE (nout,*)
       WRITE (nout,*) 'Students t distribution'
       WRITE (nout,*)
       WRITE (nout,*) '             Correct    Parameter   Standard'
       WRITE (nout,*) '             values     estimates    errors'
       DO i = 1, npar
          WRITE (nout,99998) 'PARAM(', i, ')', param(i), theta(i), se(i)
       END DO
       WRITE (nout,99999) 'DF      ', df, theta(npar+1), se(npar+1)
       IF (mn==1) THEN
          WRITE (nout,99999) 'Mean    ', mean, theta(npar+1+mn), se(npar+1+mn)
       END IF
       DO i = 1, nreg
          WRITE (nout,99998) 'BX(', i, ')   ', bx(i), theta(npar+1+mn+i), &
             se(npar+1+mn+i)
       END DO

       nt = 4
       ifail = 0
       CALL g13fhf(num,nt,ip,iq,theta,cvar,htm,etm,ifail)
       WRITE (nout,*)
       WRITE (nout,99997) 'Volatility forecast = ', cvar(nt)
       WRITE (nout,*)

99999  FORMAT (1X,A,3(F10.2,2X))
99998  FORMAT (1X,A,I1,A,3(F10.2,2X))
99997  FORMAT (1X,A,F12.2)
    END PROGRAM g05hnfe
