Program g02qgfe ! G02QGF Example Program Text ! Mark 24 Release. NAG Copyright 2012. ! .. 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 Procedures .. 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 100 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,*) 100 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