Program g11cafe ! G11CAF Example Program Text ! Mark 24 Release. NAG Copyright 2012. ! .. Use Statements .. Use nag_library, Only: g11caf, nag_wp ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Integer, Parameter :: nin = 5, nout = 6 ! .. Local Scalars .. Real (Kind=nag_wp) :: dev, tol Integer :: cm, i, ifail, ip, iprint, ldz, lwk, & m, maxit, n, n0, ns ! .. Local Arrays .. Real (Kind=nag_wp), Allocatable :: b(:), cov(:), sc(:), se(:), wk(:), & z(:,:) Integer, Allocatable :: cnt(:), ic(:), isi(:), isz(:), & nca(:), nct(:) ! .. Intrinsic Procedures .. Intrinsic :: count, maxval, sum ! .. Executable Statements .. Write (nout,*) 'G11CAF Example Program Results' Write (nout,*) ! Skip heading in data file Read (nin,*) ! Read in problem size and control parameters Read (nin,*) n, m, ns, maxit, iprint, tol ldz = n Allocate (z(ldz,m),isz(m),ic(n),isi(n),nca(ns),nct(ns),cnt(ns)) ! Read in data Read (nin,*)(isi(i),ic(i),z(i,1:m),i=1,n) ! Read in variable inclusion flags Read (nin,*) isz(1:m) ! Calculate IP ip = count(isz(1:m)>0) ! Calculate number of observations in the model and maximum number of ! observations in any stratum cnt(1:ns) = 0 Do i = 1, n If (isi(i)>0 .And. isi(i)<=ns) Then cnt(isi(i)) = cnt(isi(i)) + 1 End If End Do cm = maxval(cnt(1:ns)) n0 = sum(cnt(1:ns)) lwk = ip*n0 + (cm+1)*(ip+1)*(ip+2)/2 + cm Allocate (b(ip),se(ip),sc(ip),cov(ip*(ip+1)/2),wk(lwk)) ! Read in initial estimate for B Read (nin,*) b(1:ip) ! Calculate parameter estimates ifail = 0 Call g11caf(n,m,ns,z,ldz,isz,ip,ic,isi,dev,b,se,sc,cov,nca,nct,tol, & maxit,iprint,wk,lwk,ifail) ! Display results Write (nout,99999) ' Deviance = ', dev Write (nout,*) Write (nout,*) ' Strata No. Cases No. Controls' Write (nout,*) Write (nout,99998)(i,nca(i),nct(i),i=1,ns) Write (nout,*) Write (nout,*) ' Parameter Estimate', ' Standard Error' Write (nout,*) Write (nout,99997)(i,b(i),se(i),i=1,ip) 99999 Format (A,E13.4) 99998 Format (3X,I3,10X,I2,10X,I2) 99997 Format (I6,10X,F8.4,10X,F8.4) End Program g11cafe