Program g13fafe

!     G13FAF Example Program Text

!     Mark 28.6 Release. NAG Copyright 2022.

!     .. Use Statements ..
      Use nag_library, Only: g13faf, g13fbf, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=nag_wp)               :: gamma, lgf, pht, tol
      Integer                          :: i, ifail, ip, iq, isym, l, ldcovr,   &
                                          ldx, lwork, maxit, mn, npar, nreg,   &
                                          nt, num, pgamma
      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)
!     .. Executable Statements ..
      Write (nout,*) 'G13FAF Example Program Results'
      Write (nout,*)

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

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

      ldx = num
      Allocate (yt(num),x(ldx,nreg))

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

!     Read in the exogenous variables
      If (nreg>0) Then
        Read (nin,*,Iostat=ifail)(x(i,1:nreg),i=1,num)
      End If

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

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

!     Calculate NPAR
      npar = 1 + iq + ip + isym + 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
      If (isym==1) Then
        Read (nin,*) theta(l)
        pgamma = l
        l = l + 1
      End If
!     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,*) pht
      End If

!     Fit the GARCH model
      ifail = -1
      Call g13faf(dist,yt,x,ldx,num,ip,iq,nreg,mn,isym,npar,theta,se,sc,covr,  &
        ldcovr,pht,et,ht,lgf,copts,maxit,tol,work,lwork,ifail)
      If (ifail/=0) Then
        If (ifail/=5 .And. ifail/=6) Then
          Go To 100
        End If
      End If

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

      Allocate (fht(nt))

!     Extract the estimate of the asymmetry parameter from theta
      If (isym==1) Then
        gamma = theta(pgamma)
      Else
        gamma = 0.0E0_nag_wp
      End If

!     Calculate the volatility forecast
      ifail = 0
      Call g13fbf(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
      If (isym==1) Then
        Write (nout,99998) ' Gamma', theta(l), se(l)
        Write (nout,*)
        l = l + 1
      End If
!     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,99997) 'Volatility forecast = ', fht(nt)
      Write (nout,*)

100   Continue

99999 Format (1X,A,I0,1X,2F16.2)
99998 Format (1X,A,1X,2F16.2)
99997 Format (1X,A,F12.2)
    End Program g13fafe