C ASYMETRICAL CORRELATION WITH MISSING DATA OCTOBER 21, 1965 DOUBLE PRECISION P,FIN,SUBPRO,PROBLM,VARSEM DIMENSION F(162),DATE(4),A(12000),Q(400) DATA MXST/12000/ C DATA PER,ONO,FIN,YES,SUBPRO,PROBLM,VARSEM/1H.,2HNO,6HFINISH,3HYES, X6HSUBPRO,6HPROBLM,6HVARSEL/ DIMENSION PC(2) LOGICAL NBL NPR=0 CALL USAGEB('BMDX84') NERR=0 GO TO 100 203 NERR=NERR+1 202 NC=NC1 NV=NV1 NVA=NVA1 IT=IT1 NF=NF1 ON=ON1 GO TO 105 100 READ (5,1) P,PC,NC,NV,NVA,IT,NF,ON 1 FORMAT(A6,2A3,5I6,4X,A2) 105 IF(P.EQ.FIN) STOP IF(P.NE.PROBLM)GO TO 399 C C C C NERR=0 405 NPR=NPR+1 IF(IT.EQ.1.OR.IT.EQ.6)WRITE(6, 502) IF(IT.EQ.0) IT=5 IF(IT.EQ.5) ON=ONO IF(ON.NE.ONO) ON=YES IF(ON.NE.ONO) REWIND IT IF(NF.GE.1.AND.NF.LE.9)GO TO 500 NF=1 WRITE(6, 501) 500 NF=18*NF IF(IT.GE.0) READ (5,10) (F(I),I=1,NF) 10 FORMAT(18A4) C C WRITE(6, HEADING WRITE (6,2) PC,NC,NV,NVA,IT,ON 2 FORMAT ('1BMDX84 - ASYMMETRICAL CORRELATION WITH MISSING DATA', * ' - REVISED MAY 10, 1968' / X41H0HEALTH SCIENCES COMPUTING FACILITY, UCLA// X31H0PROBLEM CODE 2A3/ X31H0NUMBER OF CASES I6/ X31H0NUMBER OF VARIABLES READ IN I6/ X31H0NUMBER OF VARIABLES ADDED I6/ X31H0INPUT TAPE NUMBER I6/ X35H0REWIND INPUT TAPE A4) WRITE (6,22) (F(I),I=1,NF) 22 FORMAT(31H0INPUT FORMAT 18A4/(31X,18A4)) C C SET ARGUMENTS AND CALL PASS1 NV1=NV+NVA L1=1+NV1 L11=L1+MAX0(NV1,NV) L2=L11+NV1 L21=L2+NV1 CALL PASS1(A,A(L11),A(L1),A(L2),A(L21),NV,NV1,NC,IT,F,NPR) NV=NV1 IPP=0 C C READ SUBPROBLEM CARD 205 READ (5,1) P,PC,NC1,NV1,NVA1,IT1,NF1,ON1 IF(P.NE.SUBPRO)GO TO 203 IPP=IPP+1 N1=NC1 N2=NV1 SSS=PC(2) IF(SSS.NE.YES) SSS=ONO C WRITE(6, SUBPROBLEM INFORMATION WRITE (6,49) IPP,SSS 49 FORMAT(31H1SUBPROBLEM NUMBER I6/ X35H0BLANKS TREATED AS MISSING A4) C C SET DIAGONAL LIMITS KH=1 IF(NBL(N1).OR.NBL(N2)) GO TO 52 N2=NV NRZ=0 52 MM1=66 M1=0 99 M0=M1+1 M1=M1+MM1 C C READ IN VARIABLE SELECTION CARDS READ (5,76) P,(Q(M),M=M0,M1) 76 FORMAT(A6,66A1) IF(P.NE.VARSEM)GO TO 402 NRZ=NRZ+1 DO 98 M=M0,M1 98 IF(Q(M).EQ.PER) GO TO (50,89),KH GO TO 99 50 CALL VARSEL(A(L2),NXX,NV,Q,NRZ) L3=L2+NXX KH=2 K2=L3 J1=1 WRITE (6,60) (Q(I),I=1,M1) 60 FORMAT(31H0VARIABLE SELECTION 100A1/(31X,100A1)) GO TO 52 89 WRITE (6,61) (Q(I),I=1,M1) 61 FORMAT(/(31X,100A1)) WRITE (6,685) N1,N2 685 FORMAT(31H0LOWER DIAGONAL LIMIT I6/ X31H0UPPER DIAGONAL LIMIT I6) 88 NP=0 CALL VARSEL(A(L3),NYY,NV,Q,NRZ) II0=MIN0(MAX0(1,J1+N1),NXX) DO 13 J=J1,NYY II1=MIN0(MAX0(1,J+N1),NXX) II2=MIN0(MAX0(1,J+N2),NXX) NP1=NP+II2-II1+1 NX=II2-II0+1 NY=J-J1+1 NST=3*NV+NXX+2*NX+3*NY+NP1 NST1=NST+5*NP1+NX+NY IF(SSS.NE.YES) GO TO 41 IF(NST1.GT.MXST) GO TO 30 GO TO 13 41 IF(NST.GT.MXST) GO TO 30 13 NP=NP1 J=NYY+1 30 NY=J-J1 II2=MIN0(MAX0(1,J-1+N2),NXX) NX=II2-II0+1 L22=L2+II0-1 M1=N1+J1-II0 M2=N2+J1-II0 J1=J K1=L3 DO 18 I=1,NY A(K1)=A(K2) K1=K1+1 18 K2=K2+1 L4=L3+NY L5=L4+NX L6=L5+NY L7=L6+NX L8=L7+NY L9=L8+NX L10=L9+NY L13=L10+NP L14=L13+NP L15=L14+NP L16=L15+NP L17=L16+NP REWIND 1 IF(SSS.NE.YES) GO TO 44 CALL CORR(A,A(L1),A(L22),A(L3),A(L4),A(L5),A(L6),A(L7),A(L8),A(L9) X,A(L10),A(L13),A(L14),A(L15),A(L16),A(L17),NX,NY,M1,M2,NV,NC,NP,SS XS) GO TO 46 44 CALL COR(A(L11),A(L1),A(L22),A(L3),A(L4),A(L5),A(L6),A(L7),A(L8),N XX,NY,M1,M2,NV,NC,NP,SSS) 46 REWIND 1 IF(J1.LE.NYY) GO TO 88 GO TO 205 399 IF(NERR)400,400,403 400 WRITE(6, 4000)P 4000 FORMAT(' PROGRAM EXPECTED PROBLM CARD INSTEAD READ THE FOLLOWING'/ 11X,A6) 501 FORMAT(' NUMBER OF VARIABLE FORMAT CARDS INCORRECTLY SPECIFIED,ASS 1UMED TO BE 1') 502 FORMAT(' INPUT TAPE CANNOT EQUAL 1 OR 6') STOP 402 WRITE(6, 4002)NRZ,P 4002 FORMAT(' VARSEL CARD',1X,I4,' READS AS FOLLOWS',1X,A6) STOP 403 WRITE(6, 4003)P 4003 FORMAT(' PROGRAM EXPECTED SUBPRO OR PROBLM CARD INSTEAD READ THE 1FOLLOWING '/1X,A6) STOP END C SUBROUTINE CON FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE CON(N,B) DIMENSION A(10),B(6) DATA A,BL,ALP/1H0,1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9,1H ,1H(/ M=N DO 1 L=1,6 1 B(L)=BL L=6 DO 2 I=1,5 J=MOD(M,10) B(L)=A(J+1) M=M/10 L=L-1 2 IF(M.EQ.0) GO TO 7 7 B(L)=ALP RETURN END C SUBROUTINE COR FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE COR(EAN,H,IX,IY,X,Y,SSX,SSY,SXY,NX,NY,N1,N2,NV,NC,NP,SS XS) DIMENSION EAN(2),H(2),IX(2),IY(2),X(2),Y(2),SSX(2),SSY(2),SXY(2) DATA BL/1H / DO 1 I=1,NX 1 SSX(I)=0. DO 2 J=1,NY 2 SSY(J)=0. DO 3 K=1,NP 3 SXY(K)=0. DO 6 L=1,NC READ (1) (H(I),I=1,NV) 400 FORMAT ( 20A4 ) DO 5 I=1,NX I1=IX(I) IF(H(I1).EQ.BL) H(I1)=0. X(I)=H(I1)-EAN(I1) 5 SSX(I)=SSX(I)+X(I)*X(I) K=1 DO 6 J=1,NY II1=MIN0(MAX0(1,J+N1),NX) II2=MIN0(MAX0(1,J+N2),NX) J1=IY(J) IF(H(J1).EQ.BL) H(J1)=0. Y(J)=H(J1)-EAN(J1) SSY(J)=SSY(J)+Y(J)*Y(J) DO 6 I=II1,II2 SXY(K)=SXY(K)+X(I)*Y(J) 6 K=K+1 K=1 DO 7 J=1,NY II1=MIN0(MAX0(1,J+N1),NX) II2=MIN0(MAX0(1,J+N2),NX) DO 7 I=II1,II2 A=-0. IF(SSX(I)*SSY(J).LE.0.0) GO TO 8 A=SXY(K)/SQRT(SSX(I)*SSY(J)) 8 SXY(K)=A 7 K=K+1 CALL OUTPUT(SXY,MXY,IX,IY,NX,NY,N1,N2,NP,NC,SSS) RETURN END C SUBROUTINE CORR FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE CORR(EAN,H,IX,IY,X,Y,SSX,SSY,MX,MY,SXMY,SYMX,SSXMY,SSYM XX,SXY,MXY,NX,NY,N1,N2,NV,NC,NP,SSS) DIMENSION IX(2),IY(2),X(2),Y(2),SSX(2),SSY(2),MX(2),MY(2),SXMY(2), XSYMX(2),SSXMY(2),SSYMX(2),SXY(2),MXY(2),EAN(2),H(2) DATA AMISS/1H / DO 1 I=1,NX SSX(I)=0. 1 MX(I)=0 DO 2 I=1,NY SSY(I)=0. 2 MY(I)=0 DO 3 I=1,NP SXY(I)=0. MXY(I)=0 SXMY(I)=0. SYMX(I)=0. SSXMY(I)=0. 3 SSYMX(I)=0. DO 8 L=1,NC READ (1) (H(I),I=1,NV) 400 FORMAT ( 20A4 ) DO 6 I=1,NX I1=IX(I) X(I)=H(I1) IF(X(I).EQ.AMISS) GO TO 7 X(I)=X(I)-EAN(I1) SSX(I)=SSX(I)+X(I)*X(I) GO TO 6 7 MX(I)=MX(I)+1 6 CONTINUE K=1 DO 8 J=1,NY II1=MIN0(MAX0(1,J+N1),NX) II2=MIN0(MAX0(1,J+N2),NX) J1=IY(J) Y(J)=H(J1) IF(Y(J).EQ.AMISS) GO TO 9 Y(J)=Y(J)-EAN(J1) SSY(J)=SSY(J)+Y(J)*Y(J) DO 10 I=II1,II2 IF(X(I).EQ.AMISS) GO TO 11 SXY(K)=SXY(K)+X(I)*Y(J) GO TO 10 11 SYMX(K)=SYMX(K)+Y(J) SSYMX(K)=SSYMX(K)+Y(J)*Y(J) 10 K=K+1 GO TO 8 9 MY(J)=MY(J)+1 DO 12 I=II1,II2 IF(X(I).EQ.AMISS) GO TO 13 SXMY(K)=SXMY(K)+X(I) SSXMY(K)=SSXMY(K)+X(I)*X(I) GO TO 12 13 MXY(K)=MXY(K)+1 12 K=K+1 8 CONTINUE K=1 DO 16 J=1,NY II1=MIN0(MAX0(1,J+N1),NX) II2=MIN0(MAX0(1,J+N2),NX) DO 15 I=II1,II2 MXY(K)=NC-MX(I)-MY(J)+MXY(K) XY=MXY(K) A=SSX(I)-SSXMY(K) B=SSY(J)-SSYMX(K) C=SXY(K) IF(XY.LE.0.0)GO TO 14 A=A-SXMY(K)*SXMY(K)/XY B=B-SYMX(K)*SYMX(K)/XY C=C-SXMY(K)*SYMX(K)/XY 14 SXY(K)=-0.0 IF(A*B.LE.0.) GO TO 15 SXY(K)=C/SQRT(A*B) 15 K=K+1 16 CONTINUE CALL OUTPUT(SXY,MXY,IX,IY,NX,NY,N1,N2,NP,NC,SSS) RETURN END C FUNCTION NBL FOR BMDX84 OCTOBER 21, 1965 LOGICAL FUNCTION NBL(X) EXTERNAL SIGN NBL=.TRUE. IF(SIGN(1.0,X).NE.1.0.AND.X.EQ.0.0)NBL=.FALSE. RETURN END C SUBROUTINE OUTPUT FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE OUTPUT(R,MC,IX,IY,NX,NY,N1,N2,NP,NC,SSS) DATA YES/3HYES/ DIMENSION R(2),MC(2),IX(2),IY(2),II(2),JJ(2) DIMENSION RR(10),I1(10),J1(10),A1(6,10) JJ(1)=0 II(1)=0 II2=0 K1=0 L1=0 6 L1=MIN0(L1+200,NP) WRITE (6,3) 3 FORMAT(1H1,4X,10(12H I CORR)/5X,10(12H J COUNT)) 5 K0=K1+1 K1=MIN0(K1+10,L1) I=0 DO 8 K=K0,K1 I=I+1 II(1)=II(1)+1 IF(II(1).LE.II2) GO TO 9 JJ(1)=JJ(1)+1 JJJ=JJ(1) II(1)=MIN0(MAX0(1,JJJ+N1),NX) II2=MIN0(MAX0(1,JJJ+N2),NX) 9 I1(I)=IX(II(1)) J1(I)=IY(JJ(1)) RR(I)=R(K) MCK=NC IF(SSS.EQ.YES) MCK=MC(K) 8 CALL CON(MCK,A1(1,I)) WRITE (6,4) (I1(K),RR(K),K=1,I) 4 FORMAT(1H0,4X,10(I5,F7.3)) WRITE (6,7) (J1(K),(A1(J,K),J=1,6),K=1,I) 7 FORMAT(5X,10(I5,6A1,1H))) IF(K1.LT.L1) GO TO 5 IF(L1.LT.NP) GO TO 6 RETURN END C SUBROUTINE PASS1 FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE PASS1(XM,XM1,X,SD,SD1,NV,NV1,NC,IT,F1,NPR) C C READ IN DATA MATRIX, INITIATE TRANSGENERATIONS AND SELECTIONS C WRITE RESULTS ON INTERMEDIATE TAPE C COMPUTE AND WRITE(6, MEANS AND STANDARD DEVIATIONS C DIMENSION XM(2),XM1(2),X(2),F1(162),SD(2),SD1(2) LOGICAL NBL DATA BIN,BL/3HBIN,1H / BBCD=1. IF(IT.LT.0) BBCD=BIN IT=IABS(IT) REWIND 1 DO 11 I=1,NV1 SD(I)=0. XM1(I)=0. 11 XM(I)=0. MY=0 DO 12 J=1,NC SELECT=1. IF(BBCD.EQ.BIN) READ (IT) (X(I),I=1,NV) IF(BBCD.NE.BIN) READ (IT,F1) (X(I),I=1,NV) CALL TRANS(X,NV1,J,NPR,SELECT) IF(SELECT.EQ.0.) GO TO 12 MY=MY+1 DO 14 I=1,NV1 IF(NBL(X(I))) GO TO 16 XM1(I)=XM1(I)+1. X(I)=BL GO TO 14 16 XM(I)=XM(I)+X(I) SD(I)=SD(I)+X(I)*X(I) 14 CONTINUE WRITE (1) (X(I),I=1,NV1) 410 FORMAT ( 20A4 ) 12 CONTINUE NC=MY H=NC DO 15 I=1,NV1 S=XM(I) X(I)=H-XM1(I) XM(I)=0.0 IF(X(I).GT.0.0)XM(I)=S/X(I) XM1(I)=S/H SD1(I)=SQRT((SD(I)-S*XM1(I))/(H-1.)) 15 SD(I)=SQRT((SD(I)-S*XM(I))/(X(I)-1.)) WRITE (6,1) (I,XM1(I),SD1(I),H,XM(I),SD(I),X(I),I=1,NV1) 1 FORMAT(76H1 BLANKS TREATED AS ZEROS BLANK XS CONSIDERED MISSING/ X77H0VARIABLE MEAN STD. DEV. COUNT MEAN ST XD. DEV. COUNT/(I6,3X,2F13.5,F8.0,2F13.5,F8.0)) ENDFILE 1 REWIND 1 RETURN END C SUBROUTINE VARSEL FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE VARSEL(LL,NX,NV,Q,NRZ) C C SELECT VARIABLES SPECIFIED FOR SUBPROBLEM C DIMENSION Q(2), LL(2), R(18) DATA R/1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9,1H0,1H,,1H-,1H.,1H ,1H( X,1H),1H(,1H)/ L=0 I=0 1 NN=0 19 M=0 M1=1 6 L=L+1 DO 2 J=1,18 2 IF(Q(L).EQ.R(J)) GO TO 3 L=L-(NRZ-1)*66+6 WRITE(6, 101)L 101 FORMAT(' CHARACTER IN COLUMN',I2,' IS ILLEGAL') 30 WRITE (6,4) NRZ 4 FORMAT(21H1ERROR ON VARSEL CARDI4) 60 STOP 3 IF(J.GT.10) GO TO 5 J=MOD(J,10) M=M*10+J GO TO 6 9 M1=M M=0 GO TO 6 5 J=J-10 GO TO (11,11,11,6,11,9,11,9), J 40 WRITE(6, 100) 100 FORMAT(' INDEX OF VARIABLE SELECTED EXCEEDS TOTAL NUMBER OF VARIAB 1LES') GO TO 30 50 WRITE(6,102),NRZ 102 FORMAT(' INDICES ARE REVERSED ON VARSEL CARD',I4) STOP 11 IF(M.GT.NV) GO TO 40 I=I+1 LL(I)=M IF(NN.EQ.0) GO TO 10 I=I-1 L1=LL(I)+M1 IF(M.LT.L1) GO TO 50 DO 13 K=L1,M,M1 I=I+1 13 LL(I)=K 10 NN=1 GO TO (1,19,20,6,19,9,19,9), J 20 NX=I RETURN END C SUBROUTINE TRANS FOR BMDX84 OCTOBER 21, 1965 SUBROUTINE TRANS(X,NV,NC,NPR,SELECT) DIMENSION X(NV) RETURN END