Program g08chfe ! Mark 24 Release. NAG Copyright 2012. ! .. 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 Procedures .. 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