* 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