PROGRAM g02qgfe ! G02QGF Example Program Text ! Mark 23 Release. NAG Copyright 2011. ! .. Use Statements .. USE nag_library, ONLY : g02qgf, g02zkf, g02zlf, g05kff, nag_wp ! .. Implicit None Statement .. IMPLICIT NONE ! .. Parameters .. INTEGER, PARAMETER :: lseed = 1, nin = 5, nout = 6 ! .. Local Scalars .. REAL (KIND=nag_wp) :: df, rvalue INTEGER :: genid, i, ifail, ip, ivalue, j, l, & ldbl, lddat, ldres, liopts, lopts, & lstate, lwt, m, n, ntau, optype, & sorder, subid, tdch CHARACTER (1) :: c1, weight CHARACTER (30) :: cvalue, semeth CHARACTER (100) :: optstr ! .. Local Arrays .. REAL (KIND=nag_wp), ALLOCATABLE :: b(:,:), bl(:,:), bu(:,:), ch(:,:,:), & dat(:,:), opts(:), res(:,:), tau(:), & wt(:), y(:) INTEGER, ALLOCATABLE :: info(:), iopts(:), isx(:), state(:) INTEGER :: seed(lseed) ! .. Intrinsic Functions .. INTRINSIC count, len_trim, min ! .. Executable Statements .. WRITE (nout,*) 'G02QGF Example Program Results' WRITE (nout,*) FLUSH (nout) ! Skip heading in data file READ (nin,*) ! Read in the problem size READ (nin,*) sorder, c1, weight, n, m, ntau ! Read in the data IF (weight=='W' .OR. weight=='w') THEN lwt = n ELSE lwt = 0 END IF ALLOCATE (wt(lwt),isx(m),y(n),tau(ntau)) IF (sorder==1) THEN ! DAT(N,M) lddat = n ALLOCATE (dat(lddat,m)) IF (lwt==0) THEN READ (nin,*) (dat(i,1:m),y(i),i=1,n) ELSE READ (nin,*) (dat(i,1:m),y(i),wt(i),i=1,n) END IF ELSE ! DAT(M,N) lddat = m ALLOCATE (dat(lddat,n)) IF (lwt==0) THEN READ (nin,*) (dat(1:m,i),y(i),i=1,n) ELSE READ (nin,*) (dat(1:m,i),y(i),wt(i),i=1,n) END IF END IF ! Read in variable inclusion flags READ (nin,*) isx(1:m) ! Calculate IP ip = count(isx(1:m)==1) IF (c1=='Y' .OR. c1=='y') THEN ip = ip + 1 END IF ! Read in the quantiles required READ (nin,*) tau(1:ntau) liopts = 100 lopts = 100 ALLOCATE (iopts(liopts),opts(lopts)) ! Initialize the optional argument array ifail = 0 CALL g02zkf('INITIALIZE = G02QGF',iopts,liopts,opts,lopts,ifail) C_LP: DO ! Read in any optional arguments. Reads in to the end of ! the input data, or until a blank line is reached ifail = 1 READ (nin,99994,IOSTAT=ifail) optstr IF (ifail/=0) THEN EXIT C_LP ELSE IF (len_trim(optstr)==0) THEN EXIT C_LP END IF ! Set the supplied option ifail = 0 CALL g02zkf(optstr,iopts,liopts,opts,lopts,ifail) END DO C_LP ! Assume that no intervals or output matrices are required ! unless the optional argument state differently ldbl = 0 tdch = 0 ldres = 0 lstate = 0 ! Query the optional arguments to see what output is required ifail = 0 CALL g02zlf('INTERVAL METHOD',ivalue,rvalue,cvalue,optype,iopts,opts, & ifail) semeth = cvalue IF (semeth/='NONE') THEN ! Require the intervals to be output ldbl = ip IF (semeth=='BOOTSTRAP XY') THEN ! Need to find the length of the state array for the random ! number generator ! Read in the generator ID and a seed READ (nin,*) genid, subid, seed(1) ! Query the length of the state array ALLOCATE (state(lstate)) ifail = 0 CALL g05kff(genid,subid,seed,lseed,state,lstate,ifail) ! Deallocate STATE so that it can reallocated later DEALLOCATE (state) END IF ifail = 0 CALL g02zlf('MATRIX RETURNED',ivalue,rvalue,cvalue,optype,iopts, & opts,ifail) IF (cvalue=='COVARIANCE') THEN tdch = ntau ELSE IF (cvalue=='H INVERSE') THEN IF (semeth=='BOOTSTRAP XY' .OR. semeth=='IID') THEN ! NB: If we are using bootstrap or IID errors then any request for ! H INVERSE is ignored tdch = 0 ELSE tdch = ntau + 1 END IF END IF ifail = 0 CALL g02zlf('RETURN RESIDUALS',ivalue,rvalue,cvalue,optype,iopts, & opts,ifail) IF (cvalue=='YES') THEN ldres = n END IF END IF ! Allocate memory for output arrays ALLOCATE (b(ip,ntau),info(ntau),bl(ldbl,ntau),bu(ldbl,ntau), & ch(ldbl,ldbl,tdch),state(lstate),res(ldres,ntau)) IF (lstate>0) THEN ! Doing bootstrap, so initialise the RNG ifail = 0 CALL g05kff(genid,subid,seed,lseed,state,lstate,ifail) END IF ! Call the model fitting routine ifail = -1 CALL g02qgf(sorder,c1,weight,n,m,dat,lddat,isx,ip,y,wt,ntau,tau,df,b, & bl,bu,ch,res,iopts,opts,state,info,ifail) IF (ifail/=0) THEN IF (ifail==231) THEN WRITE (nout,*) 'Additional error information (INFO): ', & info(1:ntau) ELSE GO TO 20 END IF END IF ! Display the parameter estimates DO l = 1, ntau WRITE (nout,99999) 'Quantile: ', tau(l) WRITE (nout,*) IF (ldbl>0) THEN WRITE (nout,*) ' Lower Parameter Upper' WRITE (nout,*) ' Limit Estimate Limit' ELSE WRITE (nout,*) ' Parameter' WRITE (nout,*) ' Estimate' END IF DO j = 1, ip IF (ldbl>0) THEN WRITE (nout,99998) j, bl(j,l), b(j,l), bu(j,l) ELSE WRITE (nout,99998) j, b(j,l) END IF END DO WRITE (nout,*) IF (tdch==ntau) THEN WRITE (nout,*) 'Covariance matrix' DO i = 1, ip WRITE (nout,99997) ch(1:i,i,l) END DO WRITE (nout,*) ELSE IF (tdch==ntau+1) THEN WRITE (nout,*) 'J' DO i = 1, ip WRITE (nout,99997) ch(1:i,i,1) END DO WRITE (nout,*) WRITE (nout,*) 'H inverse' DO i = 1, ip WRITE (nout,99997) ch(1:i,i,l+1) END DO WRITE (nout,*) END IF WRITE (nout,*) END DO IF (ldres>0) THEN WRITE (nout,*) 'First 10 Residuals' WRITE (nout,*) ' Quantile' WRITE (nout,99995) 'Obs.', tau(1:ntau) DO i = 1, min(n,10) WRITE (nout,99996) i, res(i,1:ntau) END DO ELSE WRITE (nout,*) 'Residuals not returned' END IF WRITE (nout,*) 20 CONTINUE 99999 FORMAT (1X,A,F6.3) 99998 FORMAT (1X,I3,3(3X,F7.3)) 99997 FORMAT (1X,10(E10.3,3X)) 99996 FORMAT (2X,I3,10(1X,F10.5)) 99995 FORMAT (1X,A,10(3X,F6.3,2X)) 99994 FORMAT (A100) END PROGRAM g02qgfe