C	WESTERN MICHIGAN UNIVERSITY
C	COLCHI.F4 (FILENAME ON LIBRARY DECTAPE)
C	COLCHI, 1.1.3 (CALLING NAME, SUBLST NO.)
C	CHI SQUARE, GAMMA, TAU-A,B,C, AND SOMER'S D WITH
C	 COLLAPSING OF CONTINGENCY TABLE CAPABILITIES.
C	COLCHI WAS PROGRAMMED BY SAM ANEMA AND LATER MODIFIED BY
C	 R. R. BARR.
C	LIBRARY DECTAPE PROGS. USED:  USAGE.MAC
C	APLIB.F4 PROGS. USED:  CHIPRB, NORMCV
C	INTERNAL SUBR. USED:  CLAP, CHIS, SOMER, THETA
C	ABOVE COMMENTS AND RIGHT ADJUSTED COMMENTS WERE PUT IN BY WG
C
      DIMENSION NFRE(40,40),NFQ(40,40),ID(14)
      DIMENSION NSR(40),NSC(40)
      IOUT=-1
      INP=-4
C	CALL USAGE('COLCHI')
	WRITE(IOUT,99)
99	FORMAT(//,' WMU - COLLAPSING CHI SQUARE',//)
100   WRITE(IOUT,101)
101   FORMAT(' ','HOW MANY ROWS?'/)
      READ(INP,102)NR
102   FORMAT(I)
      WRITE(IOUT,103)
103   FORMAT(' ','HOW MANY COLUMNS?'/)
      READ(INP,102)NC
107   WRITE(IOUT,109)
109   FORMAT(' ','ENTER IDENTIFICATION.'/)
      READ(INP,110)(ID(I),I=1,14)
110   FORMAT(14A5)
106   WRITE(IOUT,213)
213   FORMAT(' ','ENTER FREQUENCIES.'/)
      DO 290 J=1,NR
290   READ(INP,113)(NFRE(I,J),I=1,NC)
113   FORMAT(40I)
      WRITE(IOUT,300)(ID(I),I=1,14)
300   FORMAT(' ',14A5)
      DO 301 I=1,NR
      DO 301 J=1,NC
301   NFQ(J,I)=NFRE(J,I)
      NNR=NR
      NNC=NC
302   CALL CHIS(NFQ,NNR,NNC,NTOT,NSR)
304   WRITE(IOUT,303)
303   FORMAT(' ','TYPE:'/3X,'1 TO TERMINATE'/3X,'2 TO ENTER',
     1' MORE DATA'/3X,'3 TO COLLAPSE'/3X,'4 FOR GAMMA STATISTICS'/
     23X,'5 FOR THETA'/)
      READ(INP,102)NTO
      GO TO (1000,100,305,306,320),NTO
305   CALL CLAP(NFRE,NR,NC,NFQ,NNR,NNC)
      GO TO 302
306   CALL SOMER(NFQ,NNR,NNC,NTOT)
      GO TO 304
320   CALL THETA(NFQ,NNR,NNC,NSR)
      GO TO 304
1000  CALL EXIT
      END
C---------------NFRE, NC ARE INPUT.  NFG, NNR, NNC ARE RETURNED.
C--------------- NR APPARENTLY NOT USED.
      SUBROUTINE CLAP(NFRE,NR,NC,NFQ,NNR,NNC)
      DIMENSION NFRE(40,40),IR(40),NFQ(40,40)
      IOUT=-1
      INP=-4
      WRITE(IOUT,10)
10    FORMAT(' ','WHAT IS THE NEW NUMBER OF ROW CATEGORIES?'/)
      READ(INP,20)NNR
20    FORMAT(40I)
      WRITE(IOUT,30) 
30    FORMAT(' ','ENTER NEW ROW CATEGORIZATION'/)
      DO 50 I=1,NNR
      READ(INP,20)(IR(K),K=1,40)
      KK=0
      DO 40 K=1,40
      IF(IR(K))70,70,60
60    KK=KK+1
40    CONTINUE
70    DO 80 K=1,NC
      LM=0
      DO 100 KR=1,KK
      IRT=IR(KR)
100   LM=LM+NFRE(K,IRT)
80    NFQ(K,I)=LM
50    CONTINUE
      DO 400 I=1,40
400   IR(I)=0
      WRITE(IOUT,200)
200   FORMAT(' ','WHAT IS THE NEW NUMBER OF COLUMN CATEGORIES?'/)
      READ(INP,20)NNC
      WRITE(IOUT,210)
210   FORMAT(' ','ENTER NEW COLUMN CATEGORIZATION.'/)
      DO 220 I=1,NNC
      READ(INP,20)(IR(K),K=1,40)
      KK=0
      DO 230 K=1,40
      IF(IR(K))240,240,250
250   KK=KK+1
230   CONTINUE
240   DO 260 K=1,NNR
      LM=0
      DO 270 KR=1,KK
      IRT=IR(KR)
270   LM=LM+NFQ(IRT,K)
260   NFQ(I,K)=LM
220   CONTINUE
      RETURN
      END
C---------------NTOT, NSR RETURNED.  OTHER ARGS. ARE INPUT.
      SUBROUTINE CHIS(NFRE,NR,NC,NTOT,NSR)
      DIMENSION NFRE(40,40),NSR(40),NSC(40)
      INP=-4
      IOUT=-1
      DO 200 J=1,NR
200   NSR(J)=0
      DO 210 I=1,NC
210   NSC(I)=0
      NTOT=0
      DO 220 J=1,NR
      DO 220 I=1,NC
      NSR(J)=NSR(J)+NFRE(I,J)
      NSC(I)=NSC(I)+NFRE(I,J)
220   CONTINUE
      DO 230 J=1,NR
230   NTOT=NTOT+NSR(J)
      CHI=0.
	IDF=(NR-1)*(NC-1)
      DO 300 I=1,NC
      DO 300 J=1,NR
      IF(NR.EQ.1.OR.NC.EQ.1)GO TO 377
      GO TO 378
377   E=FLOAT(NTOT)/FLOAT(NR*NC)
      GO TO 379
378   E=FLOAT(NSR(J)*NSC(I))/FLOAT(NTOT)
379   CHI=CHI+((FLOAT(NFRE(I,J))-E)**2)/E
300   CONTINUE
	CALL CHIPRB(CHI,IDF,CPRB)
      WRITE(IOUT,310)
310   FORMAT(' ',20X,'CONTINGENCY TABLE'//)
      WRITE(IOUT,320)(I,I=1,NC)
320   FORMAT(' ','VAR',I3,40I6)
      DO 330 J=1,NR
330   WRITE(IOUT,340)J,(NFRE(I,J),I=1,NC),NSR(J)
340   FORMAT(/' ',I2,I4,41I6)
      WRITE(IOUT,350)(NSC(I),I=1,NC),NTOT
350   FORMAT(/' ',41I6)
      WRITE(IOUT,360)CHI,CPRB
360   FORMAT(//' ','CHI-SQUARE = ',F12.5,4X,'PROB =',F8.5)
      IF(NR.LT.2.OR.NC.LT.2)RETURN
      CONCOF=SQRT(CHI/(FLOAT(NTOT)+CHI))
      PS=CHI/FLOAT(NTOT)
      L=MIN0(NC,NR)
      PP=SQRT(PS/(FLOAT(L)-1.))
      IF(NC.NE.2.OR.NR.NE.2)GO TO 371
      A=NFRE(1,1)
      B=NFRE(1,2)
      C=NFRE(2,1)
      D=NFRE(2,2)
      TOT=NTOT
      CHI=TOT*(ABS(A*D-B*C)-TOT/2.0)**2/((A+B)*(C+D)*(A+C)*(B+D))
	CALL CHIPRB(CHI,IDF,CPRB)
      WRITE(IOUT,372)CHI,CPRB
372   FORMAT(' 2X2 CORRECTED CHI-SQUARE = 'F12.5,4X,'PROB =',F8.5)
371   WRITE(IOUT,370)CONCOF,PS,PP,IDF
370   FORMAT(' ','CONTINGENCY COEFFICIENT = ',F12.5/' PHI-SQUARE = ',
     1F12.5/' PHI-PRIME = ',F12.5/ ' DEGREES OF FREEDOM = ',I4)
      RETURN
      END
      SUBROUTINE SOMER(NFRE,NR,NC,NTOT)
      DIMENSION NFRE(40,40)
      IOUT=-1
      INP=-4
      NC1=NC-1
      NR1=NR-1
      NP=0
      DO 400 K2=1,NC1
      DO 400 K=1,NR1
      NT=0
      LL=K2+1
      KK=K+1
      DO 410 I=LL,NC
      DO 410 J=KK,NR
      NT=NT+NFRE(I,J)
410   CONTINUE
      NP=NP+NT*NFRE(K2,K)
400   CONTINUE
      NQ=0
      DO 420 K2=2,NC
      DO 420 K=1,NR1
      NT=0
      LL=K2-1
      KK=K+1
      DO 430 I=1,LL
      DO 430 J=KK,NR
      NT=NT+NFRE(I,J)
430   CONTINUE
      NQ=NQ+NT*NFRE(K2,K)
420   CONTINUE
      NX=0
      DO 440 K=1,NC
      DO 440 L=1,NR1
      L1=L+1
       DO 440 J=L1,NR
      NX=NX+NFRE(K,L)*NFRE(K,J)
440   CONTINUE
      NY=0
      DO 450 K=1,NR
      DO 450 L=1,NC1
      L1=L+1
      DO 450 J=L1,NC
      NY=NY+NFRE(L,K)*NFRE(J,K)
450   CONTINUE
500   FORMAT(' ',8I7)
      P=NP
      Q=NQ
      TOT=NTOT
      X=NX
      Y=NY
      GAMMA=(P-Q)/(P+Q)
      TAUA=(2.*(P-Q))/(TOT*(TOT-1.))
      TAUB=(P-Q)/SQRT((P+Q+X)*(P+Q+Y))
      EM=AMIN0(NR,NC)
      TAUC=((P-Q)*2.*EM)/(TOT**2*(EM-1.))
      DYX=(P-Q)/(P+Q+Y)
      DXY=(P-Q)/(P+Q+X)
      WRITE(IOUT,460)GAMMA,TAUA,TAUB,TAUC,DYX,DXY
460   FORMAT(/' ','GAMMA = ',F12.5/' TAU-A = ',F12.5/' TAU-B = ',
     1F12.5/' TAU-C = ',F12.5/' DYX = ',F12.5/' DXY = ',F12.5)
      RETURN
      END
C---------------ALL ARGS. ARE INPUT
      SUBROUTINE THETA(NFQ,NR,NC,NSR)
      DIMENSION NFQ(40,40),NSR(40)
      IOUT=-1
      INP=-4
      NTT=0
      NC1=NC-1
      NR1=NR-1
      DO 500 I=1,NR1
      II=I+1
      DO 500 J=II,NR
      NTT=NTT+NSR(I)*NSR(J)
500   CONTINUE
      ND=0
      DO 510 J=1,NR1
      JJ=J+1
      DO 510 K=JJ,NR
      NB=0
      DO 530 I=1,NC1
      II=I+1
      DO 530 L=II,NC
      NB=NB+NFQ(I,J)*NFQ(L,K)
530   CONTINUE
      NA=0
      DO 540 I=2,NC
      II=I-1
      DO 540 L=1,II
      NA=NA+NFQ(I,J)*NFQ(L,K)
540   CONTINUE
      ND=ND+IABS(NB-NA)
510   CONTINUE
      THET=FLOAT(ND)/FLOAT(NTT)
      WRITE(IOUT,520)ND,NTT,THET
520   FORMAT(' ','D = ',I5/' T2 = ',I5/' THETA = ',F12.5)
      RETURN
      END
