Program g02defe

!     G02DEF Example Program Text

!     Mark 27.3 Release. NAG Copyright 2021.

!     .. Use Statements ..
      Use nag_library, Only: g02ddf, g02def, nag_wp
!     .. Implicit None Statement ..
      Implicit None
!     .. Parameters ..
      Integer, Parameter               :: nin = 5, nout = 6
!     .. Local Scalars ..
      Real (Kind=nag_wp)               :: rss, rsst, tol
      Integer                          :: i, idf, ifail, ip, irank, ldq, lwt,  &
                                          m, n
      Logical                          :: svd
      Character (1)                    :: weight
!     .. Local Arrays ..
      Real (Kind=nag_wp), Allocatable  :: b(:), cov(:), p(:), q(:,:), se(:),   &
                                          wk(:), wt(:), x(:)
!     .. Executable Statements ..
      Write (nout,*) 'G02DEF Example Program Results'
      Write (nout,*)

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

!     Read in the problem size
      Read (nin,*) n, m, weight

      If (weight=='W' .Or. weight=='w') Then
        lwt = n
      Else
        lwt = 0
      End If
      ldq = n
      Allocate (b(m),cov(m*(m+1)/2),p(m*(m+2)),q(ldq,m+1),se(m),wk(m*m+5*m),wt &
        (n),x(n))

!     Read in the dependent variable, Y, and store in first column of Q
      Read (nin,*) q(1:n,1)

!     Read in weights
      If (lwt>0) Then
        Read (nin,*) wt(1:n)
      End If

!     Use suggested value for tolerance
      tol = 0.000001E0_nag_wp

!     Loop over each of the supplied variables
      ip = 0
u_lp: Do
        Read (nin,*,Iostat=ifail) x(1:n)
        If (ifail/=0) Then
          Exit u_lp
        End If

!       Add the new variable to the model
        ifail = -1
        Call g02def(weight,n,ip,q,ldq,p,wt,x,rss,tol,ifail)
        If (ifail/=0) Then
          If (ifail==3) Then
            Write (nout,99999) ' * Variable ', ip,                             &
              ' is linear combination of previous columns'
            Write (nout,*) '   so it has not been added'
            Write (nout,*)
            Cycle u_lp
          Else
            Go To 100
          End If
        End If

        ip = ip + 1
        Write (nout,99999) 'Variable ', ip, ' added'

!       Get G02DDF to recalculate RSS
        rsst = 0.0E0_nag_wp

!       Calculate the parameter estimates
        ifail = 0
        Call g02ddf(n,ip,q,ldq,rsst,idf,b,se,cov,svd,irank,p,tol,wk,ifail)

        If (svd) Then
          Write (nout,*) 'Model not of full rank'
          Write (nout,*)
        End If
        Write (nout,99998) 'Residual sum of squares = ', rsst
        Write (nout,99999) 'Degrees of freedom = ', idf
        Write (nout,*)
        Write (nout,*) 'Variable   Parameter estimate   Standard error'
        Write (nout,*)
        Write (nout,99997)(i,b(i),se(i),i=1,ip)
        Write (nout,*)
      End Do u_lp

100   Continue

99999 Format (1X,A,I0,A)
99998 Format (1X,A,E13.4)
99997 Format (1X,I6,2E20.4)
    End Program g02defe