CID BMDX63 DOUBLE PRECISION PP,P1,P2,NA DOUBLE PRECISION XX,X DIMENSION NA(9) COMMON XX(5000),X(100),FMT(324),NI,MM,ND,NH,MT,NO, 1KA,KB,NT,NMP,QA,QC,QD DATA P1,P2/6HFINISH,6HPROBLM/ CALL USAGEB('BMDX63') 20 READ (5,1) PP,NI,ND,NH,NO,MT,NF,NA 1 FORMAT(A6,3I2,I4,2I2,8A6,A4) IF(PP.EQ.P1)GO TO 22 23 IF(PP.EQ.P2)GO TO 21 24 WRITE (6,25) CALL EXIT 25 FORMAT(30H0ILLEGAL PROBLM OR FINISH CARD ) 21 NT=NI+ND WRITE (6,82) NA,NI,ND,NO,NH 82 FORMAT('1BMDX63 - MULTIVARIATE GENERAL LINEAR HYPOTHESIS - REVISED 1 ','APRIL 14, 1969'/ X41H HEALTH SCIENCES COMPUTING FACILITY, UCLA//1H0,8X,9A6// 129H NO. OF INDEPENDENT VARIABLES I5/ 229H NO. OF DEPENDENT VARIABLES I5/ 329H NO. OF OBSERVATIONS I5/ 429H NO. OF HYPOTHESES I5///) IF(NT*NT+ND*ND- 5000)3,2,2 4 FORMAT(1H1,9A6//22H TOO MANY VARIABLES ) 2 WRITE (6,4) NA 22 CALL EXIT 3 LLLL=0 IF(MT)40,40,41 40 MT=5 GO TO 42 41 IF(MT-5)43,42,43 43 LLLL=1 REWIND MT 42 NF=18*MAX0(NF,1) WRITE(6,50) 50 FORMAT(///38H0D-TYPE VARIABLE FORMAT CARDS FOR DATA) READ (5,8) (FMT(I),I=1,NF) 8 FORMAT(18A4) WRITE(6,9)(FMT(I),I=1,NF) 9 FORMAT(1X,18A4) L=NI+1 M=NI*NT+1 N=M+NI K=NT*NT+1 CALL DOIT1(XX,XX(M),XX(L),XX(N),XX(K),NT,ND) IF(LLLL)20,20,46 46 REWIND MT GO TO 20 END SUBROUTINE DOIT1(XX,XY,YX,YY,P,NT1,NT2) DOUBLE PRECISION XX(NT1,NT1),XY(NT1,NT1) DOUBLE PRECISION YX(NT1,NT1),YY(NT1,NT1),P(NT2,NT2) DOUBLE PRECISION XXX,X 1,D,TT DOUBLE PRECISION H,UIQ,HYP COMMON XXX(5000),X(100),FMT(324),NI,MM,ND,NH,MT,NO, 1KA,KB,NT,NMP,QA,QC,QD DATA HYP/6HHYPOTH/ DO 85 I=1,NT NNN=I DO 85 J=1,NNN 85 XX(I,J)=0.0 DO 6 N=1,NO READ (MT,FMT) (X(I),I=1,NT) DO 6 I=1,NT DO 6 J=1,I 6 XX(I,J)=XX(I,J)+X(I)*X(J) DO 7 I=1,NT DO 7 J=1,I 7 XX(J,I)=XX(I,J) WRITE (6,5000) 5000 FORMAT(34H-CROSS PRODUCT MATRIX (X,Y)'(X,Y)//) DO 400 I=1,NT 400 WRITE (6,99) I,(XX(I,J),J=1,NT) CALL MATINV(XX,NI,D,NT,TT) WRITE (6,5005) TT 5005 FORMAT(1H-,10X,2H-1/6X,5H(X'X)5X,10H TOLERANCE F12.8//) DO 405 I=1,NI 405 WRITE (6,99) I,(XX(I,J),J=1,NI) 50 CALL MPY(XX,XY,XY,NI,NI,ND,NT,NT,NT,1) CALL MPY(YX,XY,P,ND,NI,ND,NT,NT,ND,1) DO 10 I=1,ND DO 10 J=1,ND YY(I,J)=YY(I,J)-P(I,J) 10 P(I,J)=YY(I,J) WRITE (6,5002) 5002 FORMAT(1H-32X,2H-1/38H REGRESSION COEFFICIENTS B=(X'X) X'Y//) DO 402 I=1,NI 402 WRITE (6,99) I,(XY(I,J),J=1,ND) 99 FORMAT(I4,10F12.5/(4X,10F12.5)) WRITE (6,5003) 5003 FORMAT(43H-RESIDUAL SUM OF PRODUCT MATRIX E=Y'Y-Y'XB//) DO 71 I=1,ND 71 WRITE (6,99) I,(YY(J,I),J=1,ND) NMP=NO-NI READ (5,66) (FMT(I),I=1,18) 66 FORMAT(18A4) WRITE(6,6666) 6666 FORMAT(39H0D-TYPE VARIABLE FORMAT CARD FOR HYPOTH) WRITE(6,77)(FMT(I),I=1,18) 77 FORMAT(1X,18A4) IF(NH)80,80,82 80 CALL EXIT 82 DO 2 IQ=1,NH READ (5,1) H,UIQ,MM,MNM,NMN 1 FORMAT(2A6,3I2) QA=MM QC=MNM QD=NMN IF(MM.NE.0) GO TO 300 MM=NI 300 IF(MNM.NE.0) GO TO 301 MNM=ND 301 IF(H.NE.HYP)GO TO 1111 302 IF(MM.GT.NI)GO TO 2222 303 IF(MNM.GT.ND)GO TO 3333 KB=MAX0(ND,NI) KA=NI L=NT*NT+1 M=KA*MAX0(ND,MM)+L N=KB*MM+M NNNN=N+ND*MAX0(MNM,MM) IF(ND*MNM+NNNN-5000)4,3,3 1111 WRITE(6,1119)UIQ 1119 FORMAT(//40H0HYPOTH IS PUNCHED WRONG ON HYPOTH CARD ,A6) STOP 2222 WRITE(6,2229)UIQ 2229 FORMAT(//46H COL.13,14 ARE SPECIFIED WRONG ON HYPOTH CARD ,A6) STOP 3333 WRITE(6,3339)UIQ 3339 FORMAT(//46H COL.15,16 ARE SPECIFIED WRONG ON HYPOTH CARD ,A6) STOP 8 FORMAT(1X,I4,2X,A6,5X,30HTHIS HYPOTHESIS IS TOO LARGE ) 3 WRITE (6,8) IQ,UIQ DO 200 I=1,MM 200 READ (5,FMT) (X(J),J=1,NI) DO 201 I=1,MNM 201 READ (5,FMT) (X(J),J=1,ND) DO 202 I=1,MM 202 READ (5,FMT) (X(J),J=1,MNM) GO TO 2 4 WRITE (6,6000) IQ,UIQ 6000 FORMAT(1H080(1H*)/29H-OUTPUT FOR HYPOTHESIS NUMBERI4,4X,A6) CALL DOIT2(XX,XY,YX,YY,XXX(N),XXX(L),XXX(M),XXX(NNNN),NT,ND,KA,KB, 1MNM) 2 CONTINUE RETURN END SUBROUTINE DOIT2(XX,XY,YX,YY,C,A,B,CEC,NT1,ND1,KA1,KB1,MMM) DOUBLE PRECISION XXX,X,DE,XX,XY,YX,YY 1,C,A,B,DEH,TT,CEC COMMON XXX(5000),X(100),FMT(324),NI,MM,ND,NH,MT,NO, 1KA,KB,NT,NMP,QA,QC,QD DIMENSION XX(NT1,NT1),XY(NT1,NT1), YX(NT1,NT1),YY(NT1,NT1),C(ND1,N 1D1),A(KA1,KA1),B(KB1,KB1),CEC(ND1,ND1) DATA Z1,Z2,Z3/1HA,1HC,1HD/ WRITE (6,201) Z1 201 FORMAT(1H-A2,6HMATRIX//) DO 3 I=1,MM IF(QA)500,500,501 500 DO 502 J=1,NI 502 A(J,I)=0. A(I,I)=1. GO TO 3 501 READ (5,FMT) (A(J,I),J=1,NI) 3 WRITE (6,200) I,(A(J,I),J=1,NI) 200 FORMAT(I4,10F12.5/(4X,10F12.5)) CALL MPY(XY,A,B,ND,NI,MM,NT,KA,KB,2) WRITE (6,201) Z2 DO 20 I=1,MMM IF(QC)508,508,509 508 DO 510 J=1,ND 510 C(J,I)=0. C(I,I)=1. GO TO 20 509 READ (5,FMT) (C(J,I),J=1,ND) 20 WRITE (6,200) I,(C(J,I),J=1,ND) CALL MPY(YY,C,CEC,ND,ND,MMM,NT,ND,ND,2) CALL MPY(C,CEC,CEC,MMM,ND,MMM,ND,ND,ND,2) CALL MPY(C,B,B,MMM,ND,MM,ND,KB,KB,2) WRITE (6,201) Z3 DO 21 J=1,MM IF(QD)564,564,566 564 DO 562 I=1,MMM 562 C(I,J)=0. GO TO 21 566 READ (5,FMT) (C(I,J),I=1,MMM) 21 WRITE (6,200) J,(C(I,J),I=1,MMM) NDD=ND ND=MMM DO 1 I=1,ND DO 1 J=1,MM 1 C(I,J)=B(I,J)-C(I,J) WRITE (6,6001) 6001 FORMAT(17H-CONTRAST ABC'-D//) DO 801 I=1,MM 801 WRITE (6,200) I,(C(J,I),J=1,ND) CALL MPY(XX,A,B,NI,NI,MM,NT,KA,KB,1) CALL MPY(A,B,B,MM,NI,MM,KA,KB,KB,2) WRITE (6,6000) 6000 FORMAT(1H-10X,2H-1/5X,10HA(X'X) A'//) DO 800 I=1,MM 800 WRITE (6,200) I,(B(I,J),J=1,MM) CALL MATINV(B,MM,D,KB,TT) 80 CALL MPY(B,C,A,MM,MM,ND,KB,NDD,KA,3) DO 10 I=1,ND DO 11 J=1,MM 11 X(J)=C(I,J) DO 12 J=1,ND 12 C(I,J)=0. DO 10 K=1,ND DO 10 J=1,MM 10 C(I,K)=C(I,K)+X(J)*A(J,K) WRITE (6,6009) 6009 FORMAT(1H-4X,4HCEC'//) DO 301 I=1,ND 301 WRITE (6,200) I,(CEC(I,J),J=1,ND) WRITE (6,6007) 6007 FORMAT(1H-52X,7H-1 -1/68H HYPOTHESIS SUM OF PRODUCT MATRIX H=(A 1BC'-D)'(A(X'X) A') (ABC'-D)//) DO 810 I=1,ND 810 WRITE (6,200) I,(C(I,J),J=1,ND) DO 2 I=1,ND DO 2 J=1,ND 2 C(I,J)=C(I,J)+CEC(I,J) CALL MATINV(C,ND,DEH,NDD,TT) CALL MATINV(CEC,ND,DE,NDD,TT) V=DE/DEH WRITE (6,310) DE,DEH,V,ND,MM,NMP 310 FORMAT(20H-DETERMINANT(CEC') F20.6/20H0DETERMINANT(CEC'+H) 1F20.6/20H0U-STATISTIC F20.6,3X,18HDEGREES OF FREEDOM3I4) ZS=ND ZR=MM ZN=NMP ZD=ZR*ZR+ZS*ZS-5. ZT=1. IF(ZD.EQ.0.) GO TO 50 ZT=SQRT((ZR*ZR*ZS*ZS-4.)/ZD) 50 ZH=(ZN-.5*(ZS-ZR+1.))*ZT-ZR*ZS*.5+1. IDF1=ZR*ZS ZY=V**(1./ZT) ZF=(1.-ZY)/ZY*ZH/(ZR*ZS) WRITE (6,55) ZF,IDF1,ZH 55 FORMAT(20H0F-STATISTIC F20.6,3X,18HDEGREES OF FREEDOM XI4,F7.0) ND=NDD RETURN END SUBROUTINE MATINV(A,N,D,ND,T) DOUBLE PRECISION A,U,C,T,D,X,Y,G DIMENSION A(ND,ND),IND(100),U(100),C(100) DO 1 I=1,N C(I)=A(I,I) 1 IND(I)=0 K=1 T=1. D=1. Y=1. DO 7 L=1,N IF(Y.LT.T)T=Y IF(K)9,9,10 10 IND(K)=1 DO 2 I=1,K U(I)=A(K,I) 2 A(K,I)=0. X=U(K) DO 3 I=K,N U(I)=A(I,K) 3 A(I,K)=0. U(K)=1. D=D*X Y=0. DO 7 I=1,N DO 4 J=1,I 4 A(I,J)=A(I,J)-U(I)*U(J)/X IF(IND(I))5,5,7 5 G=A(I,I)/C(I) IF(Y-G)6,7,7 6 Y=G K=I 7 CONTINUE DO 11 I=1,N DO 11 J=1,I A(I,J)=-A(I,J) 11 A(J,I)=A(I,J) RETURN 9 D=0. RETURN END SUBROUTINE MPY(A,B,C,II,JJ,KK,L,M,N,NNN) DOUBLE PRECISION A,B,C,U DIMENSION A(L,L),B(M,M),C(N,N),U(100) DO 1 K=1,KK GO TO (2,2,7),NNN 2 DO 4 J=1,JJ 4 U(J)=B(J,K) 3 GO TO (5,6,7),NNN 5 DO 8 I=1,II C(I,K)=0. DO 8 J=1,JJ 8 C(I,K)=C(I,K)+A(I,J)*U(J) GO TO 1 6 DO 9 I=1,II C(I,K)=0. DO 9 J=1,JJ 9 C(I,K)=C(I,K)+A(J,I)*U(J) GO TO 1 7 DO 10 I=1,II C(I,K)=0. DO 10 J=1,JJ 10 C(I,K)=C(I,K)+A(I,J)*B(K,J) 1 CONTINUE RETURN END