    PROGRAM g13fefe

!      G13FEF Example Program Text

!      Mark 23 Release. NAG Copyright 2011.

!      .. Use Statements ..
       USE nag_library, ONLY : g13fef, g13fff, nag_wp
!      .. Implicit None Statement ..
       IMPLICIT NONE
!      .. Parameters ..
       INTEGER, PARAMETER              :: nin = 5, nout = 6
!      .. Local Scalars ..
       REAL (KIND=nag_wp)              :: gamma, hp, lgf, tol
       INTEGER                         :: i, ifail, ip, iq, l, ldcovr, ldx,    &
                                          lwork, maxit, mn, npar, nreg, nt,    &
                                          num, pgamma, tdx
       LOGICAL                         :: tdist
       CHARACTER (1)                   :: dist
!      .. Local Arrays ..
       REAL (KIND=nag_wp), ALLOCATABLE :: covr(:,:), et(:), fht(:), ht(:),     &
                                          sc(:), se(:), theta(:), work(:),     &
                                          x(:,:), yt(:)
       LOGICAL                         :: copts(2)
!      .. Intrinsic Functions ..
       INTRINSIC                          max
!      .. Executable Statements ..
       WRITE (nout,*) 'G13FEF Example Program Results'
       WRITE (nout,*)

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

!      Read in the problem size
       READ (nin,*) num, mn, nreg

       ldx = num
       tdx = max(nreg+mn,1)
       ALLOCATE (yt(num),x(ldx,tdx))

!      Read in the series
       READ (nin,*) yt(1:num)

!      Read in the exogenous variables
       IF (nreg>0) THEN
          READ (nin,*) (x(i,1:nreg),i=1,num)
       END IF

!      Read in details of the model to fit
       READ (nin,*) dist, ip, iq

!      Read in control parameters
       READ (nin,*) copts(1:2), maxit, tol

!      Calculate NPAR
       npar = 2 + iq + ip + mn + nreg
       IF (dist=='T' .OR. dist=='t') THEN
          npar = npar + 1
          tdist = .TRUE.
       ELSE
          tdist = .FALSE.
       END IF

       ldcovr = npar
       lwork = (nreg+3)*num + npar + 403
       ALLOCATE (theta(npar),se(npar),sc(npar),covr(ldcovr,npar),et(num), &
          ht(num),work(lwork))

!      Read in initial values
!      alpha_0
       READ (nin,*) theta(1)
       l = 2
!      alpha_i       
       IF (iq>0) THEN
          READ (nin,*) theta(l:(l+iq-1))
          l = l + iq
       END IF
!      beta_i
       IF (ip>0) THEN
          READ (nin,*) theta(l:(l+ip-1))
          l = l + ip
       END IF
!      gamma
       READ (nin,*) theta(l)
       pgamma = l
       l = l + 1
!      degrees of freedom
       IF (tdist) THEN
          READ (nin,*) theta(l)
          l = l + 1
       END IF
!      mean
       IF (mn==1) THEN
          READ (nin,*) theta(l)
          l = l + 1
       END IF
!      Regression parameters and pre-observed conditional variance
       IF ( .NOT. copts(2)) THEN
          READ (nin,*) theta(l:(l+nreg-1))
          READ (nin,*) hp
       END IF

!      Fit the GARCH model
       ifail = 0
       CALL g13fef(dist,yt,x,ldx,num,ip,iq,nreg,mn,npar,theta,se,sc,covr, &
          ldcovr,hp,et,ht,lgf,copts,maxit,tol,work,lwork,ifail)

!      Read in forecast horizon
       READ (nin,*) nt

       ALLOCATE (fht(nt))

!      Extract the estimate of the asymmetry parameter from theta
       gamma = theta(pgamma)

!      Calculate the volatility forecast
       ifail = 0
       CALL g13fff(num,nt,ip,iq,theta,gamma,fht,ht,et,ifail)

!      Output the results
       WRITE (nout,*) '               Parameter        Standard'
       WRITE (nout,*) '               estimates         errors '
!      Output the coefficient alpha_0
       WRITE (nout,99999) 'Alpha', 0, theta(1), se(1)
       l = 2
!      Output the coefficients alpha_i
       IF (iq>0) THEN
          WRITE (nout,99999) ('Alpha',i-1,theta(i),se(i),i=l,l+iq-1)
          l = l + iq
       END IF
       WRITE (nout,*)
!      Output the coefficients beta_j
       IF (ip>0) THEN
          WRITE (nout,99999) (' Beta',i-l+1,theta(i),se(i),i=l,l+ip-1)
          l = l + ip
          WRITE (nout,*)
       END IF
!      Output the estimated asymmetry parameter, gamma
       WRITE (nout,99998) ' Gamma', theta(l), se(l)
       WRITE (nout,*)
       l = l + 1
!      Output the estimated degrees of freedom, df
       IF (dist=='T') THEN
          WRITE (nout,99998) '    DF', theta(l), se(l)
          WRITE (nout,*)
          l = l + 1
       END IF
!      Output the estimated mean term, b_0
       IF (mn==1) THEN
          WRITE (nout,99999) '    B', 0, theta(l), se(l)
          l = l + 1
       END IF
!      Output the estimated linear regression coefficients, b_i
       IF (nreg>0) THEN
          WRITE (nout,99999) ('    B',i-l+1,theta(i),se(i),i=l,l+nreg-1)
       END IF
       WRITE (nout,*)

!      Display the volatility forecast
       WRITE (nout,*)
       WRITE (nout,99997) 'Volatility forecast = ', fht(nt)
       WRITE (nout,*)

99999  FORMAT (1X,A,I0,1X,2F16.2)
99998  FORMAT (1X,A,1X,2F16.2)
99997  FORMAT (1X,A,F12.2)
    END PROGRAM g13fefe
