C	WESTERN MICHIGAN UNIVERSITY
C	KENTAU.F4 (FILE NAME ON LIBRARY DECTAPE)
C	KENTAU, 1.10.1 (CALLING NAME, SUBLST. NO.)
C	KENDALL'S TAU, RANK ORDER CORRELATION
C	PROGRAMMED BY SAM ANEMA
C	LIBRARY DECTAPE PROGS. USED:  USAGE.MAC
C	FORWMU PROGS. USED:  TTYPTY, DEVCHG, EXISTS, PRINTS
C	APLIB PROGS. USED:  IO, GETFOR
C	INTERNAL SUBR. USED:  RATAU
C	ABOVE COMMENTS AND RIGHT ADJUSTED COMMENTS PUT IN BY WG
C
      DIMENSION Z(20,20)
C---------------T CONTAINS CORRELATION MATRIX. V CONTAINS RAW DATA.
C--------------- X, Y ARE WORKING STORAGES FOR CALCULATION OF
C--------------- CORRELATION MATRIX.  Z CONTAINS Z-VALUES FOR A
C--------------- KENDALL'S TAU.  ID CONTAINS IDENT. FOR OUTPUT.
      DIMENSION V(20,200),X(200),Y(200),T(20,20)
      DIMENSION IFMT(96),ID(16)
      IOUT=3
      INP=2
      NOUTD=-1
      IND=-4
      WRITE(NOUTD,1)
1     FORMAT(20X,'WMU'/5X,'KENDALL TAU CORRELATION PROGRAM'/)
C	CALL USAGE('KENTAU')
C---------------ICODE RETURNED, =0 MEANS TERMINAL JOB, =1 MEANS
C--------------- PSEUDO-TERMINAL
      CALL TTYPTY(ICODE)
C---------------NDV1, NDV2 ARE RETURNED.  OTHER ARGS. ARE INPUT.
C--------------- 1 CAUSES OUTPUT? TO PRINT.  0 CAUSES INPUT? TO PRINT.
C--------------- ICODE COMES FROM CALL TTYPTY.
      CALL IO(IOUT,NDV1,NOUTD,IND,1,ICODE)
9     CALL IO(INP,NDV2,NOUTD,IND,0,ICODE)
C---------------IFMT, ISTD ARE RETURNED. OTHER ARGS. ARE INPUT.
C---------------96=NO. OF FMT. WORDS FOR OBJ. TIME FORMAT (6 LINES)
C--------------- 2 MEANS F-TYPE ONLY
      CALL GETFOR(NOUTD,IND,IFMT,ISTD,96,2)
      IF(ISTD.EQ.1)IFMT(1)='(10F)'
104   WRITE(NOUTD,10)
10    FORMAT(' ','ENTER NO. OF VARIABLES.'/)
      READ(IND,20)NVAR
      IF(NVAR.GT.20.OR.NVAR.LT.2)GO TO 101
      GO TO 202
101   WRITE(NOUTD,102)
102   FORMAT(' ILLEGAL NUMBER OF VARIABLES.'/)
      IF(ICODE.EQ.-1)CALL EXIT
      GO TO 104
20    FORMAT(I)
202   WRITE(NOUTD,25)
25    FORMAT(' ','ENTER MISSING DATA SYMBOL.'/)
      READ(IND,27)GONE
27    FORMAT(F)
1003  WRITE(NOUTD,1004)
1004  FORMAT(' ','ENTER IDENTIFICATION'/)
      READ(IND,1005)(ID(I),I=1,16)
1005  FORMAT(16A5)
      IF(NDV2.EQ.'TTY')GO TO 2
      WRITE(NOUTD,31)
31    FORMAT(' YOUR DATA IS BEING READ.'/)
      GO TO 32
2     IF(ISTD.EQ.1)WRITE(NOUTD,33)
33    FORMAT(' ENTER DATA (AT MOST 10 PER LINE)'/)
      IF(ISTD.NE.1)WRITE(NOUTD,30)
30    FORMAT(' ','ENTER DATA.'/)
32    K=0
60    K=K+1
39    READ(INP,IFMT,END=120)(V(J,K),J=1,NVAR)
      IF(K.GT.200)GO TO 41
       GO TO 60
41    WRITE(NOUTD,42)
42    FORMAT(' TOO MANY OBSERVATIONS'/)
      IF(ICODE.EQ.-1)CALL EXIT
      GO TO 9
120   K=K-1
      WRITE(IOUT,130)
130   FORMAT(15X,'KENDALL TAU CORRELATION MATRIX'/)
      WRITE(IOUT,2000)(ID(I),I=1,16)
2000  FORMAT(' ',16A5)
      DO 220 I=1,NVAR
      DO 220 J=1,I
      DO 140 L=1,K
      X(L)=V(I,L)
140   Y(L)=V(J,L)
      LS=0
      DO 160 L=1,K
      IF(X(L).EQ.GONE)GO TO 160
150   IF(Y(L).EQ.GONE)GO TO 160
170   LS=LS+1
      X(LS)=X(L)
      Y(LS)=Y(L)
160   CONTINUE
      NGO=LS
      IF(NGO-1)5001,5001,5002
5001  WRITE(30,5003)J,I
5003  FORMAT(' ','FEWER THAN 2 DATA PAIRS AT ',2I3)
      Z(I,J)=999.99
      T(I,J)=0.
      T(J,I)=0.
      GO TO 220
5000  FORMAT(' ',10 F 7.0)
5002  CALL RATAU(X,NGO,TX)
400   FORMAT(' ',10F7.2)
      CALL RATAU(Y,NGO,TY)
      S=0.
      DO 210 IM=1,NGO-1
      DO 210 IJ=IM+1,NGO
      SS=(X(IM)-X(IJ))*(Y(IM)-Y(IJ))
      IF(SS)200,210,211
200   S=S-1.
      GO TO 210
211   S=S+1.
210   CONTINUE
      XNGO=NGO
      XNGO=.5*XNGO*(XNGO-1.)
      IF(XNGO-TX)705,705,704
704   IF(XNGO-TY)705,705,706
705   WRITE(IOUT,710)J,I
710   FORMAT(' ','S.D. = 0 AT',2I3)
      T(I,J)=0.0
      GO TO 1124
706   T(I,J)=S/(SQRT(XNGO-TX)*SQRT(XNGO-TY))
1124  T(J,I)=T(I,J)
      IF(NGO-10)503,502,502
502   Z(I,J)=T(I,J)/SQRT(FLOAT(4*NGO+10)/FLOAT(9*NGO*(NGO-1)))
      GO TO 220
503   Z(I,J)=999.99
220   CONTINUE
      DO 230 I=1,NVAR
230   WRITE(IOUT,240)(T(I,J),J=1,I)
240   FORMAT(' ',10F7.4)
      WRITE(NOUTD,51)
51    FORMAT(' WOULD YOU LIKE Z-VALUES?'/)
      READ(IND,52)IH
52    FORMAT(A3)
      IF(IH.NE.'YES')GO TO 53
500   WRITE(IOUT,501)
501   FORMAT(///' ',20X,'Z-VALUES'//)
      DO 504 I=2,NVAR
504   WRITE(IOUT,510)(Z(I,J),J=1,I-1)
510   FORMAT(' ',10F7.2)
      WRITE(IOUT,723)
723   FORMAT(///)
53    WRITE(NOUTD,54)
54    FORMAT(' WOULD YOU LIKE PARTIAL CORRELATIONS?'/)
      READ(IND,52)IH
      IF(IH.NE.'YES')GO TO 9
260   WRITE(NOUTD,270)
270   FORMAT(' ','ENTER VARIABLES TO BE PARTIALLY CORRELATED'/)
      READ(IND,280)IX,IY
      IF(IX.EQ.0)GO TO 9
      IF(IX.GT.NVAR.OR.IX.LT.1.OR.IY.LT.1.OR.IY.GT.NVAR.OR.IX.EQ.IY)
     1GO TO 105
      GO TO 106
105   WRITE(NOUTD,107)
107   FORMAT(' ERROR ABOVE!'/)
      IF(ICODE.EQ.-1)CALL EXIT
      GO TO 260
280   FORMAT(2I)
106   WRITE(NOUTD,290)
290   FORMAT(' ','ENTER VARIABLE TO BE PARTIALLED OUT'/)
      READ(IND,20)IZ
      IF(IZ.EQ.IX.OR.IZ.EQ.IY.OR.IZ.GT.NVAR.OR.IZ.LT.1)GO TO 110
      GO TO 111
110   WRITE(NOUTD,107)
      IF(ICODE.EQ.-1)CALL EXIT
      GO TO 106
111   TXYZ=(T(IX,IY)-T(IZ,IY)*T(IX,IZ))/SQRT((1.-T(IZ,IY)**2)
     1*(1.-T(IZ,IX)**2))
      WRITE(IOUT,295)IX,IY,IZ,TXYZ
295   FORMAT(' ','PARTIAL OF VAR',I3,' WITH VAR',I3,' WITH VAR',I3,
     1' PARTIALED OUT IS ',F8.5/)
      GO TO 260
300   CALL EXIT
      END
      SUBROUTINE RATAU(X,NGO,TX)
      DIMENSION X(1),Y(200)
      TX=0.
      DO 1 I=1,NGO
      DO 2 J=1,NGO
      IF(I-J)3,2,3
3     IF(X(I)-X(J))2,4,2
4     TX=TX+1.
2      CONTINUE
1      CONTINUE
      TX=TX/2.
      RETURN
      END

