PROGRAM g08chfe ! Mark 23 Release. NAG Copyright 2011. ! .. Use Statements .. USE nag_library, ONLY : g05kff, g05sff, g08chf, nag_wp ! .. Implicit None Statement .. IMPLICIT NONE ! .. Parameters .. INTEGER, PARAMETER :: genid = 1, lseed = 1, nin = 5, & nout = 6 ! .. Local Scalars .. REAL (KIND=nag_wp) :: a2, aa2, beta, nupper, p, sa2, sbeta INTEGER :: i, ifail, j, k, lstate = 17, n, & nsim, n_pseudo, subid = -1 LOGICAL :: issort ! .. Local Arrays .. REAL (KIND=nag_wp), ALLOCATABLE :: x(:), xsim(:), y(:) INTEGER :: seed(lseed), state(17) ! .. Intrinsic Functions .. INTRINSIC exp, real, sum ! .. Executable Statements .. WRITE (nout,*) 'G08CHF Example Program Results' WRITE (nout,*) ! Skip heading in data file READ (nin,*) ! Read number of observations READ (nin,*) n ! Memory allocation ALLOCATE (x(n),y(n)) ! Read observations READ (nin,*) (x(i),i=1,n) ! Maximum likelihood estimate of mean beta = sum(x(1:n))/real(n,kind=nag_wp) ! PIT, using exponential CDF with mean beta DO i = 1, n y(i) = 1.0E0_nag_wp - exp(-x(i)/beta) END DO ! Let g08chf sort the (approximately) uniform variates issort = .FALSE. ! Calculate A-squared ifail = 0 a2 = g08chf(n,issort,y,ifail) aa2 = (1.0E0_nag_wp+0.6E0_nag_wp/real(n,kind=nag_wp))*a2 ! Number of simulations nsim = 888 ! Generate exponential variates using a repeatable seed ALLOCATE (xsim(n*nsim)) seed(1) = 206033 ifail = 0 CALL g05kff(genid,subid,seed,lseed,state,lstate,ifail) n_pseudo = n*nsim ifail = 0 CALL g05sff(n_pseudo,beta,state,xsim,ifail) ! Simulations loop nupper = 0.0E0_nag_wp DO j = 1, nsim k = (j-1)*n ! Maximum likelihood estimate of mean sbeta = sum(xsim(k+1:k+n))/real(n,kind=nag_wp) ! PIT DO i = 1, n y(i) = 1.0E0_nag_wp - exp(-xsim(k+i)/sbeta) END DO ! Calculate A-squared ifail = 0 sa2 = g08chf(n,issort,y,ifail) IF (sa2>aa2) THEN nupper = nupper + 1.0E0_nag_wp END IF END DO ! Simulated upper tail probability value p = nupper/real(nsim+1,kind=nag_wp) ! Results WRITE (nout,'(1X,A,E11.4)') & 'H0: data from exponential distribution with mean', beta WRITE (nout,'(1X,A,1X,F8.4)') 'Test statistic, A-squared: ', a2 WRITE (nout,'(1X,A,1X,F8.4)') 'Upper tail probability: ', p END PROGRAM g08chfe