hfnplt(3f) - [M_datapac:LINE_PLOT] generate a half-normal probability
plot
SUBROUTINE HFNPLT(X,N)
REAL(kind=wp),intent(in) :: X(:)
INTEGER,intent(in) :: N
HFNPLT(3f) generates a halfnormal probability plot.
The prototype halfnormal distribution used herein has mean = sqrt(2/pi)
= 0.79788456 and standard deviation = 1.
This distribution is defined for all non-negative X and has the
probability density function
f(X) = (2/sqrt(2*pi)) * exp(-X*X/2)
The prototype halfnormal distribution used herein is the distribution
of the variate X = abs(Z) where the variate Z is normally distributed
with mean = 0 and standard deviation = 1.
As used herein, a probability plot for a distribution is a plot
of the ordered observations versus the order statistic medians for
that distribution.
The halfnormal probability plot is useful in graphically testing
the composite (that is, location and scale parameters need not be
specified) hypothesis that the underlying distribution from which
the data have been randomly drawn is the halfnormal distribution.
If the hypothesis is true, the probability plot should be near-linear.
A measure of such linearity is given by the calculated probability
plot correlation coefficient.
X The vector of (unsorted or sorted) observations.
N The integer number of observations in the vector X.
The maximum allowable value of N for this subroutine is 7500.
A one-page halfnormal probability plot.
Sample program:
program demo_hfnplt
use M_datapac, only : hfnplt
implicit none
! call hfnplt(x,y)
end program demo_hfnplt
Results:
The original DATAPAC library was written by James Filliben of the
Statistical Engineering Division, National Institute of Standards
and Technology.
John Urban, 2022.05.31
CC0-1.0
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| real(kind=wp), | intent(in) | :: | X(:) | |||
| integer, | intent(in) | :: | N |
| Type | Attributes | Name | Initial | |||
|---|---|---|---|---|---|---|
| real | :: | WS(15000) |
SUBROUTINE HFNPLT(X,N) REAL(kind=wp),intent(in) :: X(:) INTEGER,intent(in) :: N REAL(kind=wp) :: W(7500), Y(7500) REAL(kind=wp) :: an , cc , hold , q , sum1 , sum2 , sum3 , tau , wbar , WS , ybar , yint , yslope INTEGER :: i , iupper COMMON /BLOCK2_real32/ WS(15000) EQUIVALENCE (Y(1),WS(1)) EQUIVALENCE (W(1),WS(7501)) ! DATA tau/1.41223913_wp/ ! iupper = 7500 ! ! CHECK THE INPUT ARGUMENTS FOR ERRORS ! IF ( N<1 .OR. N>iupper ) THEN WRITE (G_IO,99001) iupper 99001 FORMAT (' ***** FATAL ERROR--THE SECOND INPUT ARGUMENT TO HFNPLT(3f) IS OUTSIDE THE ALLOWABLE (1,', & & I0,') INTERVAL *****') WRITE (G_IO,99002) N 99002 FORMAT (' ','***** THE VALUE OF THE ARGUMENT IS ',I0,' *****') RETURN ELSEIF ( N==1 ) THEN WRITE (G_IO,99003) 99003 FORMAT (' ***** NON-FATAL DIAGNOSTIC--THE SECOND INPUT ARGUMENT TO HFNPLT(3f) HAS THE VALUE 1 *****') RETURN ELSE hold = X(1) DO i = 2 , N IF ( X(i)/=hold ) GOTO 50 ENDDO WRITE (G_IO,99004) hold 99004 FORMAT (' ***** NON-FATAL DIAGNOSTIC--THE FIRST INPUT ARGUMENT (A VECTOR) TO HFNPLT(3f) HAS ALL ELEMENTS = ',& & E15.8,' *****') ! !-----START POINT----------------------------------------------------- ! 50 continue an = N ! ! SORT THE DATA ! CALL SORT(X,N,Y) ! ! GENERATE UNIFORM ORDER STATISTIC MEDIANS ! CALL UNIMED(N,W) ! ! COMPUTE HALFNORMAL ORDER STATISTIC MEDIANS ! DO i = 1 , N q = W(i) q = (q+1.0_wp)/2.0_wp CALL NORPPF(q,W(i)) ENDDO ! ! PLOT THE ORDERED OBSERVATIONS VERSUS ORDER STATISTICS MEDIANS. ! WRITE OUT THE TAIL LENGTH MEASURE OF THE DISTRIBUTION ! AND THE SAMPLE SIZE. ! CALL PLOT(Y,W,N) WRITE (G_IO,99005) tau , N 99005 FORMAT (' ','HALFNORMAL PROBABILITY PLOT (TAU = ',E15.8,')',52X,'THE SAMPLE SIZE N = ',I0) ! ! COMPUTE THE PROBABILITY PLOT CORRELATION COEFFICIENT. COMPUTE LOCATION AND SCALE ESTIMATES ! FROM THE INTERCEPT AND SLOPE OF THE PROBABILITY PLOT. THEN WRITE THEM OUT. ! sum1 = 0.0_wp sum2 = 0.0_wp DO i = 1 , N sum1 = sum1 + Y(i) sum2 = sum2 + W(i) ENDDO ybar = sum1/an wbar = sum2/an sum1 = 0.0_wp sum2 = 0.0_wp sum3 = 0.0_wp DO i = 1 , N sum1 = sum1 + (Y(i)-ybar)*(Y(i)-ybar) sum2 = sum2 + (Y(i)-ybar)*(W(i)-wbar) sum3 = sum3 + (W(i)-wbar)*(W(i)-wbar) ENDDO cc = sum2/SQRT(sum3*sum1) yslope = sum2/sum3 yint = ybar - yslope*wbar WRITE (G_IO,99006) cc , yint , yslope 99006 FORMAT (' ','PROBABILITY PLOT CORRELATION COEFFICIENT = ',F8.5,& & 5X,'ESTIMATED INTERCEPT = ',E15.8,3X,'ESTIMATED SLOPE = ',E15.8) ENDIF ! END SUBROUTINE HFNPLT