*     G02GPF Example Program Text
*     Mark 22 Revised. NAG Copyright 2008.
*     .. Parameters ..
      INTEGER          NIN, NOUT
      PARAMETER        (NIN=5,NOUT=6)
      INTEGER          NMAX, MMAX, LDX, LDV, PNMAX
      PARAMETER        (NMAX=5,MMAX=2,LDX=NMAX,LDV=NMAX,PNMAX=5)
*     .. Local Scalars ..
      DOUBLE PRECISION A, EPS, RSS, S, TOL
      INTEGER          I, IDF, IFAIL, IP, IPRINT, IRANK, J, M, MAXIT, N
      LOGICAL          VFOBS
      CHARACTER        ERROR, LINK, MEAN, OFFSET, WEIGHT
*     .. Local Arrays ..
      DOUBLE PRECISION B(MMAX), COV((MMAX*MMAX+MMAX)/2), ETA(PNMAX),
     +                 FWT(NMAX), OFF(PNMAX), PRED(PNMAX), SE(MMAX),
     +                 SEETA(PNMAX), SEPRED(PNMAX), T(1), V(LDV,7+MMAX),
     +                 WK((MMAX*MMAX+3*MMAX+22)/2), WT(PNMAX),
     +                 X(LDX,MMAX), Y(NMAX)
      INTEGER          ISX(MMAX)
*     .. External Subroutines ..
      EXTERNAL         G02GAF, G02GPF
*     .. Executable Statements ..
      CONTINUE
      WRITE (NOUT,*) 'G02GPF Example Program Results'
*
*     Skip headings in data file
      READ (NIN,*)
      READ (NIN,*)
*     Read in training data for model that will be used for prediction
      READ (NIN,*) LINK, MEAN, OFFSET, WEIGHT, N, M, S, IPRINT
      IF (N.GT.NMAX .OR. M.GE.MMAX) THEN
         WRITE (NOUT,99994)
         GO TO 200
      END IF
      IF (WEIGHT.EQ.'W' .OR. WEIGHT.EQ.'w') THEN
         DO 20 I = 1, N
            READ (NIN,*) (X(I,J),J=1,M), Y(I), FWT(I)
   20    CONTINUE
      ELSE
         DO 40 I = 1, N
            READ (NIN,*) (X(I,J),J=1,M), Y(I)
   40    CONTINUE
      END IF
      READ (NIN,*) (ISX(J),J=1,M)
*     Calculate IP
      IP = 0
      DO 60 J = 1, M
         IF (ISX(J).GT.0) THEN
            IP = IP + 1
         END IF
   60 CONTINUE
      IF (MEAN.EQ.'M' .OR. MEAN.EQ.'m') THEN
         IP = IP + 1
      END IF
      IF (LINK.EQ.'E' .OR. LINK.EQ.'e') THEN
         READ (NIN,*) A
      END IF
*     Set control parameters
      EPS = 0.000001D0
      TOL = 0.00005D0
      MAXIT = 10
      IFAIL = 1
*     Call routine to fit model to training data
      CALL G02GAF(LINK,MEAN,OFFSET,WEIGHT,N,X,LDX,M,ISX,IP,Y,FWT,S,A,
     +            RSS,IDF,B,IRANK,SE,COV,V,LDV,TOL,MAXIT,IPRINT,EPS,WK,
     +            IFAIL)
      IF (IFAIL.NE.0) THEN
         WRITE (NOUT,99996) IFAIL
         IF (IFAIL.LT.6) THEN
            GO TO 200
         END IF
      END IF
*     Display parameter estimates for training data
      WRITE (NOUT,*)
      WRITE (NOUT,99999) RSS, IDF
      WRITE (NOUT,*)
      WRITE (NOUT,*) '      Estimate     Standard error'
      WRITE (NOUT,*)
      DO 80 I = 1, IP
         WRITE (NOUT,99998) B(I), SE(I)
   80 CONTINUE
*
*     Skip second lot of headings in data file
      READ (NIN,*)
*     Read in data to predict from and check array sizes
      READ (NIN,*) N, VFOBS, OFFSET, WEIGHT
      IF (N.GT.PNMAX) THEN
         WRITE (NOUT,*) 99994
         GO TO 200
      END IF
      IF (OFFSET.EQ.'Y' .AND. WEIGHT.EQ.'W') THEN
         DO 100 I = 1, N
            READ (NIN,*) (X(I,J),J=1,M), OFF(I), WT(I)
  100    CONTINUE
      ELSE IF (OFFSET.EQ.'Y') THEN
         DO 120 I = 1, N
            READ (NIN,*) (X(I,J),J=1,M), OFF(I)
  120    CONTINUE
      ELSE IF (WEIGHT.EQ.'W') THEN
         DO 140 I = 1, N
            READ (NIN,*) (X(I,J),J=1,M), WT(I)
  140    CONTINUE
      ELSE
         DO 160 I = 1, N
            READ (NIN,*) (X(I,J),J=1,M)
  160    CONTINUE
      END IF
*     Using G02GAF to fit training model, so error structure is normal
      ERROR = 'N'
*     Call prediction routine
      IFAIL = 1
      CALL G02GPF(ERROR,LINK,MEAN,OFFSET,WEIGHT,N,X,LDX,M,ISX,IP,T,OFF,
     +            WT,S,A,B,COV,VFOBS,ETA,SEETA,PRED,SEPRED,IFAIL)
      IF (IFAIL.NE.0) THEN
         WRITE (NOUT,99995) IFAIL
         IF (IFAIL.NE.22) THEN
            GO TO 200
         END IF
      END IF
*     Display predicted values
      WRITE (NOUT,*)
      WRITE (NOUT,*) '  I      ETA         SE(ETA)      Predicted  ',
     +  ' SE(Predicted)'
      WRITE (NOUT,*)
      DO 180 I = 1, N
         WRITE (NOUT,99997) I, ETA(I), SEETA(I), PRED(I), SEPRED(I)
  180 CONTINUE
*
  200 CONTINUE
*
99999 FORMAT (1X,'Residual sum of squares = ',E12.4,' Degrees of freed',
     +       'om = ',I2)
99998 FORMAT (1X,2F14.4)
99997 FORMAT (1X,I3,')',4(F10.5,4X))
99996 FORMAT (1X,/1X,' ** G02GAF returned with IFAIL = ',I5)
99995 FORMAT (1X,/1X,' ** G02GPF returned with IFAIL = ',I5)
99994 FORMAT (1X,' ** Problem size too large, increase array limits')
      END
