CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC C C C P R I M A L S (VERSION 3.0 - DECEMBER 1984) C C (VERSION 2.0 - DECEMBER 1981) C C (VERSION 1.0 - MARCH 1981) C C C C O N E D I M E N S I O N A L H O M O G E N E I T Y A N A L Y S I S C C C C C C C C INPUT DESCRIPTION AND OPTIONS C C C C C C -CARD 1: TITLE (20A4) C C C C -CARD 2: NOBJ,NOBJA,NVAR,NVARA,MAXCAT,ISUMCA (6I5) C C C C NOBJ = NUMBER OF OBJECTS C C NOBJA = NUMBER OF OBJECTS ACTIVE IN THE ANALYSIS C C NVAR = NUMBER OF VARIABLES C C NVARA = NUMBER OF VARIABLES ACTIVE IN THE ANALYSIS C C MAXCAT = MAXIMUM NUMBER OF CATEGORIES C C ISUMCA = TOTAL NUMBER OF CATEGORIES C C C C -CARD 3: NCAT(J),J=1,NVAR (I5) C C C C NCAT(J)= NUMBER OF CATEGORIES OF VARIABLE J C C C C -CARD 4: NORM,IEIG,NITER,CRIT,MISSOL,MISEST (6I5) C C C C NORM = 1: CATEGORIES ARE IN THE CENTROID C C = 2: OBJECTS ARE IN THE CENTROID C C IEIG = 0: MAXIMIZE LARGEST EIGENVALUE C C = 1: MINIMIZE SMALLEST EIGENVALUE C C NITER = MAXIMUM NUMBER OF ITERATIONS C C CRIT = CONVERGENCE CRITERION C C MISSOL = 0: MISSING DATA PASSIVE CATEGORIES C C = 1: MISSING DATA MULTIPLE CATEGORIES C C MISEST = 0: SUBSTITUTE OBJECT SCORES FOR MISSING DATA C C = 1: SUBSTITUTE 0.0 FOR MISSING DATA C C = 2: AVERAGE OF OBJECT SCORES WITH MISSING DATA C C C C -CARD 5: LUN(K),K=1,5 (5I5) C C C C LUN(1) = INPUT UNIT DATA C C LUN(2) = OUTPUT UNIT OBJECT SCORES/CATEGORY QUANTIFIC. C C LUN(3) = OUTPUT UNIT TRANSFORMED DATA MATRIX C C LUN(4) = OUTPUT UNIT CORRELATION MATRICES C C LUN(5) = OUTPUT UNIT PRINCIPAL COMPONENT SOLUTIONS C C C C -CARD 6: IPRI(K),K=1,10 (10I5) C C C C IPRI(1)= PRINT DATA MATRIX C C IPRI(2)= PRINT OBJECT SCORES C C IPRI(3)= PRINT CATEGORY QUANTIFICATIONS C C IPRI(4)= PRINT TRANSFORMED DATA MATRIX C C IPRI(5)= PRINT CORRELATION MATRIX BEFORE TRANSFORMATION C C IPRI(6)= PRINT CORRELATION MATRIX AFTER TRANSFORMATION C C IPRI(7)= COMPONENT SCORES BEFORE TRANSFORMATION C C IPRI(8)= COMPONENT SCORES AFTER TRANSFORMATION C C IPRI(9)= COMPONENT LOADINGS BEFORE TRANSFORMATION C C IPRI(10)=COMPONENT LOADINGS AFTER TRANSFORMATION C C C C -CARD 7: IPLO(K),K=1,11 (11I5) C C C C IPLO(K)= PLOT OPTION K, K IDENTICAL TO IPRI(K) C C IPLO(11)=PLOT COMPONENT LOADINGS B&A TRANSFORMATION C C C C -CARD 8: IWRI(K),K=1,10 (10I5) C C C C IWRI(K)= WRITE OPTION K, K IDENTICAL TO IPRI(K) C C C -CARD 9: FORMAT OF THE DATA MATRIX (20A4) C C C C -FOLLOWING CARDS: THE DATA MATRIX C C C C CARDS 1 - 9 ARE READ FROM INPUT MEDIUM C C INPARA (IN THIS VERSION FIXED TO BE 5). C C C C TO ESTIMATE THE TOTAL NUMBER OF WORDS NEEDED FOR THE ARRAY AREA, C C USE THE FORMULA NWORDS= C C 3*(NOBJ+ISUMCA) + 5*NVAR + NVAR*(NOBJ+NVAR) + NW1 + NIND C C C C C C IN THE STATIC ALLOCATION VERSION, NWORDS IS FIXED AT 22482. C C C C AUTHOR : JACQUELINE MEULMAN C C DEPARTMENT OF DATA THEORY / RUL - FSW C C MIDDELSTEGRACHT 4 - 2312 TW LEIDEN - THE NETHERLANDS C C C C REFERENCES: C C VAN DE GEER, J.P. & MEULMAN, J. PRIMALS USER'S GUIDE. C C DEPT. OF DATA THEORY, LEIDEN 1985 C C GIFI A. NON LINEAR MULTIVARIATE ANALYSIS. DEPT. OF DATA C C THEORY, LEIDEN 1981 C C GIFI A. NON LINEAR MULTIVARIATE ANALYSIS. DSWO PRESS, C C LEIDEN 1985 C C MEULMAN, J. HOMOGENEITY ANALYSIS OF INCOMPLETE DATA. DSWO C C PRESS, LEIDEN 1982 C C C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC EXTERNAL MAIN2 COMMON/CP1/NOBJ,NOBJA,NVAR,NVARA,MAXCAT,ISUMCA,NW1,NIND,TITLE(20) COMMON/CP2/IWRITE,INPARA C C READER AND PRINTER SPECIFICATION; CHANGE HERE AND IN SUBROUTINE C DECLAR IF OTHER CONVENTIONS ARE USED ON YOUR MACHINE C INPARA=5 IWRITE=6 C CALL COMMEN READ (INPARA,500) TITLE READ (INPARA,501) NOBJ,NOBJA,NVAR,NVARA,MAXCAT,ISUMCA C IF (MAXCAT.GT.NVAR) NW1=MAXCAT IF (MAXCAT.LE.NVAR) NW1=NVAR IF (NOBJ.GT.ISUMCA) NIND=NOBJ+NVAR IF (NOBJ.LE.ISUMCA) NIND=ISUMCA+NVAR C C ** ALLOCATE STORAGE ** C C XY = NOBJ+ISUMCA C DATA = NOBJ*NVAR C DRDC = NOBJ+ISUMCA C DRIDCI = NOBJ+ISUMCA C NCAT = NVAR C ICUMCA = NVAR C CM = NVAR*NVAR C WORK1 = NW1 C WORK2 = NVAR C WORK3 = NVAR+NVAR C IND = NIND C NWORDS = 3*(NOBJ+ISUMCA) + 5*NVAR + NVAR*(NOBJ+NVAR) + NW1 + NIND CALL DECLAR(MAIN2,A,NWORDS) STOP 500 FORMAT(20A4) 501 FORMAT(6I5) END SUBROUTINE DECLAR(SUBRTN,B,NWORDS) C ******************************************************************* C * * C * D E C L A R * C * * C * PURPOSE: A STATIC STORAGE ALLOCATION ROUTINE * C * * C ******************************************************************* DIMENSION A(22482) NDIM=22482 IF (NWORDS.LE.NDIM) GO TO 10 WRITE (6,600) NWORDS RETURN 10 CONTINUE CALL SUBRTN(A,NWORDS) RETURN 600 FORMAT (50H0,ARRAY 'A' IN SUBROUTINE 'DECLAR' IS TO SMALL FOR, X 9H THIS JOB/21H 'A' MUST BE AT LEAST,I8,5H LONG// X 45H EXECUTION OF THE PROGRAM HAS BEEN TERMINATED) END SUBROUTINE MAIN2(A,NWORDS) C ******************************************************************* C * * C * M A I N 2 * C * * C * PURPOSE: COMPUTES THE RELATIVE ADRESSES FOR ALL ARAYS IN THE * C * REST OF THE PROGRAM * C * * C ******************************************************************* DIMENSION A(NWORDS) COMMON/CP1/NOBJ,NOBJA,NVAR,NVARA,MAXCAT,ISUMCA,NW1,NIND,TITLE(20) COMMON/CP2/IWRITE,INPARA C IXY = 1 IDATA = 1+NOBJ+ISUMCA IDRDC = IDATA+NOBJ*NVAR IRICI = IDRDC+NOBJ+ISUMCA INCAT = IRICI+NOBJ+ISUMCA ICUMC = INCAT+NVAR CM = ICUMC+NVAR IWO1 = CM+NVAR*NVAR IWO2 = IWO1+NW1 IWO3 = IWO2+NVAR IIND = IWO3+NVAR+NVAR C NXY = NOBJ+ISUMCA NW3 = NVAR+NVAR C CALL MAIN3(A(IXY),A(IDATA),A(IDRDC),A(IRICI),A(INCAT),A(ICUMC), X A(CM), A(IWO1), A(IWO2), A(IWO3), A(IIND),NXY,NW3) RETURN END SUBROUTINE MAIN3(XY,DATA,DRDC,DRIDCI,NCAT,ICUMCA,CM,WORK1,WORK2, X WORK3,IND,NXY,NW3) C ****************************************************************** C * * C * M A I N 3 * C * * C * PURPOSE: CONTROLS THE PROGRAM * C * * C ****************************************************************** COMMON/CP1/NOBJ,NOBJA,NVAR,NVARA,MAXCAT,ISUMCA,NW1,NIND,TITLE(20) COMMON/CP2/IWRITE,INPARA DIMENSION XY(NXY),DATA(NOBJ,NVAR),DRDC(NXY),DRIDCI(NXY), X NCAT(NVAR),ICUMCA(NVAR),CM(NVAR,NVAR),WORK1(NW1), X WORK2(NVAR),WORK3(NW3),IND(NIND),FORMAT(20),IPRI(10), X IWRI(10),IPLO(11),LUN(5) INTEGER DATA C READ (INPARA,504) (NCAT(J),J=1,NVAR) READ (INPARA,502) NORM,IEIG,NITER,CRIT,MISSOL,MISEST READ (INPARA,503) (LUN(K),K=1,5) READ (INPARA,506) (IPRI(K),K=1,10) READ (INPARA,506) (IPLO(K),K=1,11) READ (INPARA,506) (IWRI(K),K=1,10) READ (INPARA,500) FORMAT C IF (MISSOL.EQ.1) MISEST=0 IF (NORM.EQ.0) NORM=1 IF (NITER.EQ.0) NITER=100 IF (CRIT.EQ.0.0) CRIT=1E-5 IF (CRIT.LT.1E-6) CRIT=1E-6 IF (NOBJA.EQ.0) NOBJA=NOBJ IF (NVARA.EQ.0) NVARA=NVAR C CALL CATINF(NCAT,NVAR,KJ,ICUMCA) CALL QUANT(DRDC(1),DRDC(NOBJ+1),DRIDCI(1),DRIDCI(NOBJ+1),XY(1), X XY(NOBJ+1),DATA,NCAT,FORMAT,TITLE,ICUMCA,CM,WORK1,NW1,WORK2, X NOBJ,NOBJA,NVAR,NVARA,MAXCAT,NORM,IEIG,NITER,CRIT,KJ,IND,NIND, X FIT,WORK3(1),WORK3(NVAR+1),IPRI,IPLO,IWRI,LUN,ISUMCA,MISSOL, X MIS,MISEST) CALL RESUL1(DATA,NOBJ,NVAR,XY(1),XY(NOBJ+1),KJ,WORK1,ICUMCA,NCAT, X DRDC(1),DRDC(NOBJ+1),DRIDCI(1),MAXCAT,NORM,IPRI,IWRI,LUN,MIS) DO 10 J=1,NVAR CALL PREPA1(DATA,NOBJ,NVAR,DRDC(1),DRDC(NOBJ+1),DRIDCI(1), X DRIDCI(NOBJ+1),ICUMCA,NCAT(J),NP,J,XY(1),XY(NOBJ+1),KJ) IF (IPLO(2).EQ.1) WRITE (IWRITE,600) J IF (IPLO(2).EQ.1) CALL PRPLOT(DRDC,DRIDCI, X IND,NP,NOBJ,1,IWRITE,0) NP=NCAT(J) IF (IPLO(3).EQ.1) WRITE (IWRITE,601) J IF (IPLO(3).EQ.1) CALL PRPLOT(DRDC(NOBJ+1),DRIDCI(NOBJ+1), X IND,NP,NP,1,IWRITE,0) 10 CONTINUE CALL RESUL2(DATA,NOBJ,NVAR,XY(1),XY(NOBJ+1),KJ,CM,WORK1,WORK2, X ICUMCA,NCAT,DRDC(NOBJ+1),DRDC(1),DRIDCI(1),DRIDCI(NOBJ+1), X MAXCAT,FIT,NORM,WORK3(1),WORK3(NVAR+1),IPRI,IPLO,IWRI,LUN, X MISEST) NP=NOBJ+NVAR IF (IPLO(8).EQ.1) WRITE (IWRITE,640) IF (IPLO(8).EQ.1) CALL PRPLOT(DRDC,DRIDCI,IND,NP,NOBJ,1,IWRITE,1) NP=NVAR IF (IPLO(7).EQ.1) WRITE (IWRITE,641) IF (IPLO(7).EQ.1) CALL PRPLOT(DRDC(NOBJ+1),DRIDCI(NOBJ+1), X IND,NP,NVAR,1,IWRITE,1) KB=NOBJ+1 KE=NOBJ+NVAR IF (IPLO(11).NE.1) GO TO 50 DO 40 K=KB,KE DRDC(K+NVAR)=DRDC(K) DRIDCI(K+NVAR)=DRIDCI(K) 40 CONTINUE 50 CALL RESUL3(DATA,NOBJ,NVAR,XY(1),XY(NOBJ+1),KJ,CM,WORK1,WORK2, X ICUMCA,NCAT,DRDC(NOBJ+1),DRDC(1),DRIDCI(1),DRIDCI(NOBJ+1), X MAXCAT,FIT,NORM,WORK3(1),WORK3(NVAR+1),IPRI,IPLO,IWRI,LUN, X MISEST) NP=NOBJ+NVAR IF (IPLO(10).EQ.1) WRITE (IWRITE,650) IF (IPLO(10).EQ.1) CALL PRPLOT(DRDC,DRIDCI,IND,NP,NOBJ,1,IWRITE,1) NP=NVAR IF (IPLO(9).EQ.1) WRITE (IWRITE,651) IF (IPLO(9).EQ.1) CALL PRPLOT(DRDC(NOBJ+1),DRIDCI(NOBJ+1), X IND,NP,NVAR,1,IWRITE,1) NP=NVAR+NVAR IF (IPLO(11).EQ.1) WRITE (IWRITE,660) IF (IPLO(11).EQ.1) CALL PRPLOT(DRDC(NOBJ+1),DRIDCI(NOBJ+1), X IND,NP,NVAR,1,IWRITE,1) STOP 500 FORMAT(20A4) 502 FORMAT(3I5,F10.8,2I5) 503 FORMAT(5I5) 504 FORMAT(16I5) 506 FORMAT(11I5) 600 FORMAT(1H1,48H TRANSFORMATION PLOT FOR VARIABLE,I5, X 50H: CATEGORY QUANTIFICATIONS (CHARACTERS) ALONG WITH, X 25H OBJECT SCORES (INTEGERS)) 601 FORMAT(1H1,48H TRANSFORMATION PLOT FOR VARIABLE,I5) 640 FORMAT(1H1,51H PCA SOLUTION BEFORE TRANSFORMATION:, X 51H COMPONENT LOADINGS (CHARACTERS), COMPONENT SCORES , X 10H(INTEGERS)) 641 FORMAT(1H1,51H PCA SOLUTION BEFORE TRANSFORMATION:, X 19H COMPONENT LOADINGS) 650 FORMAT(1H1,50H PCA SOLUTION AFTER TRANSFORMATION:, X 51H COMPONENT LOADINGS (CHARACTERS), COMPONENT SCORES , X 10H(INTEGERS)) 651 FORMAT(1H1,50H PCA SOLUTION AFTER TRANSFORMATION:, X 19H COMPONENT LOADINGS) 660 FORMAT(1H1,42H COMPONENT LOADINGS BEFORE , X 48H(CHARACTERS) AND AFTER TRANSFORMATION (INTEGERS)) END SUBROUTINE QUANT (DR,DC,DRINV,DCINV,XROW,YCOL,DATA,NCAT,FORMAT, X TITLE,ICUMCA,CM,WORK1,NW1,WORK2,NOBJ,NOBJA,NVAR, X NVARA,MAXCAT,NORM,IEIG,NITER,CRIT,KJ,IND,NIND, X FIT,WORK3A,WORK3B,IPRI,IPLO,IWRI,LUN,ISUMCA, X MISSOL,MIS,MISEST) C ****************************************************************** C * * C * Q U A N T * C * * C * PURPOSE : COMPUTATION OF THE ONE DIMENSIONAL HOMOGENEITY * C * ANALYSIS SOLUTION * C * SUBROUTINES CALLED: READIN, DIAGOA, PRINTA, INITIA, ITERA, * C * OUTPUT, TRANYA, TRANXA * C * * C ****************************************************************** DIMENSION DR(NOBJ),DC(KJ),DRINV(NOBJ),DCINV(KJ),XROW(NOBJ), X YCOL(KJ),DATA(NOBJ,NVAR),NCAT(NVAR),FORMAT(20), X ICUMCA(NVAR),WORK1(NW1),CM(NVAR,NVAR),WORK2(NVAR), X TITLE(20),IND(NIND),WORK3A(NVAR),WORK3B(NVAR),IPRI(10), X IPLO(11),IWRI(10),LUN(5) COMMON/CP2/IWRITE,INPARA INTEGER DATA C CALL READIN (DATA,NOBJ,NVAR,LUN(1),FORMAT) CALL DIAGOA(DATA,NOBJ,NOBJ,NVAR,KJ,ICUMCA,DC,DCINV, X DR,DRINV,NCAT,SRC,MIS,IND,MISSOL) C CALL PRINTA(TITLE,DATA,NOBJ,NOBJA,NVAR,NVARA,KJ,DC, X DR,NORM,NITER,CRIT,FORMAT,ICUMCA,NCAT,MAXCAT, X IND(1),IND(NVAR+1),IPRI,IPLO,IWRI,LUN, X ISUMCA,IEIG,MISEST,MISSOL) DO 1 I=1,NOBJ XROW(I)=0.0 1 CONTINUE DO 2 I=1,KJ YCOL(I)=0.0 2 CONTINUE C NVART=NVAR IF (NORM.EQ.1) NVAR=NVARA NOBJT=NOBJ IF (NORM.EQ.2) NOBJ=NOBJA C IF (KJ.GT.ISUMCA) WRITE (IWRITE,900) ISUMCA,KJ IF (KJ.GT.ISUMCA) GO TO 1001 C C ** FOR COMMENT SEE FORMAT STATEMENT IN LINE 900 C MACA=NCAT(1) DO 20 J=2,NVAR IF (NCAT(J).GT.MACA) MACA=NCAT(J) 20 CONTINUE IF (MACA.GT.MAXCAT) WRITE (IWRITE,1900) MAXCAT,MACA IF (MACA.GT.MAXCAT) GO TO 1001 C C ** FOR COMMENT SEE FORMAT STATEMENT IN LINE 1900 C C LMIN=0 IF (IEIG.EQ.1.AND.NORM.EQ.2) LMIN=1 DO 6 I=1,NVAR IF (IND(I).GT.0.AND.LMIN.EQ.1) WRITE (IWRITE,1000) IF (IND(I).GT.0.AND.LMIN.EQ.1) GO TO 1001 6 CONTINUE C C ** FOR COMMENT SEE FORMAT STATEMENT IN LINE 1000 ** C C ** WHEN NOT ALL OBJECTS OR VARIABLES ARE ACTIVE IN THE ANALYSIS, ** C ** NEW COMPUTATION OF THE MARGINAL FREQUENCIES IS NECESSARY ** C IF ((NOBJT.NE.NOBJ).OR.(NVART.NE.NVAR)) X CALL DIAGOA(DATA,NOBJ,NOBJT,NVAR,KJ,ICUMCA, X DC,DCINV,DR,DRINV,NCAT,SRC,MIS,IND,MISSOL) C C ** ARRAY 'IND' CONTAINS THE NUMBER OF MISSING VALUES PER VARIABLE ** C DO 3 I=1,NVAR IF (IND(I).EQ.NOBJ) WRITE (IWRITE,800) IF (IND(I).EQ.NOBJ) GO TO 1001 3 CONTINUE C C ** FOR COMMENT SEE FORMAT STATEMENT IN LINE 800 C C ** WHEN MISSING DATA SHOULD BE TREATED AS MULTIPLE CATEGORIES ** C ** (MISSOL=1) WITH NORMALIZED CATEGORY QUANTIFICATIONS (NORM=2), ** C ** THE FLOW OF THE PROGRAM IS CHANGED BECAUSE OF COMPUTATIONAL ** C ** EFFICIENCY. AFTER CONVERGENCE (LINE 450) THE DESIRED NORMALIZA- ** C ** TION IS RESTORED. C NOR2=0 IF (MISSOL.EQ.1.AND.NORM.EQ.2) NOR2=NORM IF (MISSOL.EQ.1.AND.NORM.EQ.2) NORM=1 C C ** COMPUTATION OF INITIAL VALUES FOR OBJECT SCORES WHEN NORM=1 OR ** C ** FOR CATEGORY QUANTIFICATIONS WHEN NORM=2 ** C IF(NORM.EQ.1) CALL INITIA(NOBJ,NVAR,NOBJ,XROW,DR,SRC,FIT,NORM,MIS, X IEIG,ICUMCA,NCAT) IF(NORM.EQ.2) CALL INITIA(KJ,NVAR,NOBJ,YCOL,DC,SRC,FIT,NORM,MIS, X IEIG,ICUMCA,NCAT) C C ** THE ITERATION PROCES (UPTO LINE 450) ** C OLDFIT=0.0 DO 50 K=1,NITER CALL ITERA(XROW,YCOL,DATA,NOBJ,NOBJT,NVAR,DR,DC, X DRINV,DCINV,SRC,KJ,ICUMCA,FIT,NORM,IEIG, X NCAT,MIS,WORK1,MAXCAT,IND) DIF=FIT-OLDFIT CALL OUTPUT(K,FIT,DIF) IF (ABS(DIF).LT.CRIT) GO TO 450 OLDFIT=FIT 50 CONTINUE C C ** CONVERGENCE HAS BEEN REACHED; IF NECESSARY THE DESIRED NORMALIZA ** C ** TION IS RESTORED ** C 450 IF (NOR2.EQ.2) NORM=2 IF (NOR2.NE.2) GO TO 100 CALL TRANYA(YCOL,XROW,DATA,NOBJ,NOBJT,NVAR,KJ,DCINV, X ICUMCA,NCAT,IEIG,NORM,WORK1,MAXCAT) DO 4 J=1,KJ YCOL(J)=YCOL(J)/(SQRT(FIT)) 4 CONTINUE DO 5 I=1,NOBJ XROW(I)=XROW(I)/(SQRT(FIT)) 5 CONTINUE C C ** IF DESIRED, CATEGORY QUANTIFICATIONS FOR PASSIVE VARIABLES OR ** C ** OBJECT SCORES FOR PASSIVE OBJECTS ARE COMPUTED TOGETHER WITH * C ** FINAL QUANTIFICATIONS OR SCORES FOR ACTIVE VARIABLES OR OBJECTS ** C 100 NVARA=NVAR NVAR=NVART NOBJA=NOBJ NOBJ=NOBJT CALL DIAGOA(DATA,NOBJ,NOBJ,NVAR,KJ,ICUMCA,DC,DCINV, X DR,DRINV,NCAT,SRC,MIS,IND,MISSOL) IF(NORM.EQ.1) CALL TRANYA(YCOL,XROW,DATA,NOBJ,NOBJT,NVAR,KJ,DCINV, X ICUMCA,NCAT,IEIG,NORM,WORK1,MAXCAT) IF(NORM.EQ.2) CALL TRANXA (XROW,YCOL,DATA,NOBJ,NOBJA,NOBJT, X NVAR,KJ,ICUMCA,NCAT,DRINV,MIS,IEIG,NORM) 200 RETURN C 800 FORMAT (1H0//,5X,45H*********************************************, X 52H****************************************************/ X 1H0,4X,47HTHE PROGRAM STOPS BECAUSE ONE OR MORE VARIABLES, X 50H HAVE ONLY MISSING VALUES. CHECK YOUR FORMAT CARD./ X 1H0,4X,45H*********************************************, X 52H****************************************************/) C 900 FORMAT (1H0//,5X,45H*********************************************, X 53H*****************************************************/ X 1H0,4X,45HTHE PROGRAM STOPS BECAUSE THE VALUE OF 'TOTAL, X 53H NUMBER OF CATEGORIES' IN THE DATA SPECIFICATION CARD/ X 1H0,4X,41HIS NOT CORRECT. CHANGE THE VALUE OF THIS , X 14HPARAMETER FROM,I5,5H TO,I5,1H./ X 1H0,4X,45H*********************************************, X 53H*****************************************************/) C 1900 FORMAT (1H0//,5X,45H*********************************************, X 55H*******************************************************/ X 1H0,4X,47HTHE PROGRAM STOPS BECAUSE THE VALUE OF 'MAXIMUM, X 53H NUMBER OF CATEGORIES' IN THE DATA SPECIFICATION CARD/ X 1H0,4X,41HIS NOT CORRECT. CHANGE THE VALUE OF THIS , X 14HPARAMETER FROM,I5,5H TO,I5,1H./ X 1H0,4X,47H***********************************************, X 53H*****************************************************/) C 1000 FORMAT (1H0//,5X,45H*********************************************, X 52H****************************************************/ X 1H0,4X,45HTHE PROGRAM STOPS BECAUSE THE COMBINATION OF , X 52HOPTIONS: NORMALIZATION = 2 AND EIGENVALUE = SMALLEST/ X 1H0,4X,45HIS NOT SUPPORTED WHEN THE DATAMATRIX CONTAINS, X 15H MISSING VALUES/ X 1H0,4X,45H*********************************************, X 52H****************************************************/) 1001 STOP END SUBROUTINE CATINF(NCAT,NVAR,KJ,ICUMCA) C ******************************************************************* C * * C * C A T I N F * C * * C * PURPOSE: CUMULATIVE AND TOTAL NUMBER OF CATEGORIES * C * * C ******************************************************************* DIMENSION ICUMCA(NVAR),NCAT(NVAR) KJ=0 DO 10 J=1,NVAR ICUMCA(J)=KJ KJ=KJ+NCAT(J) 10 CONTINUE RETURN END SUBROUTINE CENORA (A,NN,SRC,FIT,IEIG) C ******************************************************************* C * * C * C E N O R A * C * * C * PURPOSE: OBJECT SCORES IN DEVIATION FROM MEAN * C * AND NORMALIZED TO NUMBER OF OBJECTS * C * IT SUPPOSES NO MISSING VALUES * C * * C ******************************************************************* DIMENSION A(NN) DOUBLE PRECISION CMEAN,SSQ CMEAN=0.0 DO 10 I=1,NN CMEAN=CMEAN+A(I) 10 CONTINUE CMEAN=CMEAN/FLOAT(NN) SSQ=0.0 DO 20 I=1,NN A(I)=A(I)-CMEAN SSQ=SSQ+A(I)*A(I) 20 CONTINUE SSQ=SSQ/FLOAT(NN) STAND=DSQRT (SSQ) DO 30 I=1,NN A(I)=A(I)/STAND 30 CONTINUE IF (IEIG.EQ.0) FIT=STAND IF (IEIG.EQ.1) FIT=1.0-STAND RETURN END SUBROUTINE CENORB (A,NN,NVAR,NOBJ,W,SRC,FIT,IEIG,NORM) C ******************************************************************* C * * C * C E N O R B * C * * C * PURPOSE: QUANTIFICATIONS OR OBJECT SCORES (IF THERE ARE * C * MISSING VALUES) IN WEIGHTED DEVIATION FROM MEAN * C * AND WEIGHTED NORMALIZED * C * * C ******************************************************************* DIMENSION A(NN),W(NN) DOUBLE PRECISION CMEAN,SSQ CMEAN=0.0 DO 10 I=1,NN CMEAN=CMEAN+A(I)*W(I) 10 CONTINUE CMEAN=CMEAN/SRC SSQ=0.0 DO 20 I=1,NN A(I)=A(I)-CMEAN SSQ=SSQ+A(I)*A(I)*W(I) 20 CONTINUE IF (NORM.EQ.1) SSQ=SSQ/ FLOAT (NN*NVAR) IF (NORM.EQ.2) SSQ=SSQ/ FLOAT (NOBJ*NVAR) STAND=DSQRT (SSQ) DO 30 I=1,NN A(I)=A(I)/STAND 30 CONTINUE IF (IEIG.EQ.0) FIT=STAND IF (IEIG.EQ.1) FIT=1.0-STAND RETURN END SUBROUTINE CENORC (A,NN,NVAR,NOBJ,W,SRC,FIT,IEIG,ICUMCA,NCAT) C ******************************************************************* C * * C * C E N O R C * C * * C * PURPOSE: QUANTIFICATIONS IN WEIGHTED DEVIATION FROM MEAN * C * AND WEIGHTED NORMALIZED IN MINIMIZING SMALLEST * C * EIGENVALUE * C * * C ******************************************************************* DIMENSION A(NN),W(NN),ICUMCA(NVAR),NCAT(NVAR) DOUBLE PRECISION CMEAN,SSQT,T,SSQ SSQT=0.0 DO 10 J=1,NVAR T=0.0 CMEAN=0.0 KE=ICUMCA(J)+1 KB=ICUMCA(J)+NCAT(J) DO 20 K=KE,KB T=T+W(K) CMEAN=CMEAN+A(K)*W(K) 20 CONTINUE CMEAN=CMEAN/T SSQ=0.0 DO 30 K=KE,KB A(K)=A(K)-CMEAN SSQ=SSQ+A(K)*A(K)*W(K) 30 CONTINUE SSQT=SSQT+SSQ 10 CONTINUE SSQT=SSQT/FLOAT (NVAR*NOBJ) STAND=DSQRT (SSQT) DO 40 I=1,NN A(I)=A(I)/STAND 40 CONTINUE IF (IEIG.EQ.0) FIT=STAND IF (IEIG.EQ.1) FIT=1.0-STAND RETURN END SUBROUTINE COMMEN C ****************************************************************** C * * C * C O M M E N * C * * C * PURPOSE : PRINTING OF GENERAL COMMENTS ABOUT THE * C * PRESENT VERSION OF THE PROGRAM * C * * C * * C ****************************************************************** COMMON/CP2/IWRITE,INPARA WRITE (IWRITE,600) WRITE (IWRITE,601) 600 FORMAT (///////// X15X,55H* * * * * * * * * * * * * * * * * * * * * * * * * * * */ X15X,55H* ..................DECEMBER 1984................... */ X15X,55H* P R I M A L S..........................VERSION 3.0 */ X15X,55H* .................................................. */ X15X,55H* SOME NOTABLE CHANGES THAT HAVE BEEN PERFORMED */ X15X,55H* COMPARED TO THE SECOND VERSION: */ X15X,55H* */ X15X,55H* - A DIFFERENT INITIALIZATION FOR THE OBJECT SCORES */ X15X,55H* OR FOR THE CATEGORY QUANTIFICATIONS (RANDOM); */ X15X,55H* - A DIFFERENT NORMALIZATION OF THE COMPONENT */ X15X,55H* SCORES (THE SUM OF SQUARES USED TO BE ONE, NOW */ X15X,55H* IT EQUALS THE NUMBER OF OBJECTS); */ X15X,55H* - THE ORIGIN HAS BEEN ADDED TO THE PLOTS FOR THE */ X15X,55H* PRINCIPAL COMPONENTS SOLUTIONS. * ) 601 FORMAT(/ X15X,55H* REFERENCES: */ X15X,55H* VAN DE GEER, J.P. & MEULMAN, J. PRIMALS USER'S */ X15X,55H* GUIDE. DEPT. OF DATA THEORY, LEIDEN 1985. */ X15X,55H* GIFI A. NONLINEAR MULTIVARIATE ANALYSIS. DEPT. */ X15X,55H* OF DATA THEORY, LEIDEN 1981. */ X15X,55H* GIFI A. NONLINEAR MULTIVARIATE ANALYSIS. DSWO */ X15X,55H* PRESS, LEIDEN 1985. */ X15X,55H* MEULMAN, J. HOMOGENEITY ANALYSIS OF INCOMPLETE */ X15X,55H* DATA. DSWO PRESS, LEIDEN 1982. */ X15X,55H* */ X15X,55H*******************************************************) RETURN END SUBROUTINE CORMA (A,B,N,NVAR,RHO,J,K,CMEAN,SSQ,NOJ,NOK,M) C ******************************************************************* C * * C * C O R M A * C * * C * PURPOSE: CALCULATE CELLS OF THE CORRELATION-MATRIX * C * * C ******************************************************************* DIMENSION A(N),B(N),CMEAN(NVAR),SSQ(NVAR) REAL NOJI,NOKI IF (J.GT.1) GO TO 50 CMEAN(J)=0.0 CMEAN(K)=0.0 NOJI=1.0/FLOAT(NOJ) IF (NOJ.EQ.0) NOJI=0.0 IF (NOJ.EQ.0) GO TO 1 NOKI=1.0/FLOAT(NOK) 1 IF (NOK.EQ.0) NOKI=0.0 IF (NOK.EQ.0) GO TO 2 2 DO 10 I=1,N CMEAN(J)=CMEAN(J)+A(I) CMEAN(K)=CMEAN(K)+B(I) 10 CONTINUE CMEAN(J)=CMEAN(J)*NOJI CMEAN(K)=CMEAN(K)*NOKI SSQ(J)=0.0 SSQ(K)=0.0 SSCP=0.0 DO 20 I=1,N IF (A(I).EQ.0.0.AND.M.EQ.1) A(I)=CMEAN(J) A(I)=A(I)-CMEAN(J) IF (B(I).EQ.0.0.AND.M.EQ.1) B(I)=CMEAN(K) B(I)=B(I)-CMEAN(K) SSQ(J)=SSQ(J)+A(I)*A(I) SSQ(K)=SSQ(K)+B(I)*B(I) SSCP=SSCP+A(I)*B(I) 20 CONTINUE IF (SSCP.EQ.0.0) RHO=0.0 IF (SSCP.EQ.0.0) GO TO 30 RHO=SSCP/SQRT (SSQ(J)*SSQ(K)) C 30 GO TO 700 30 IF (J.EQ.1) GO TO 700 50 SSCP=0.0 DO 60 I=1,N IF (A(I).EQ.0.0.AND.M.EQ.1) A(I)=CMEAN(J) A(I)=A(I)-CMEAN(J) IF (B(I).EQ.0.0.AND.M.EQ.1) B(I)=CMEAN(K) B(I)=B(I)-CMEAN(K) SSCP=SSCP+A(I)*B(I) 60 CONTINUE IF (SSCP.EQ.0.0) RHO=0.0 IF (SSCP.EQ.0.0) GO TO 700 RHO=SSCP/SQRT (SSQ(J)*SSQ(K)) 700 RETURN END SUBROUTINE DIAGOA(DATA,NOBJ,NOBJT,NVAR,KJ,ICUMCA,DC,DCINV,DR, X DRINV,NCAT,SRC,MIS,NMIS,MISSOL) C ****************************************************************** C * * C * D I A G O A * C * * C * PURPOSE: COMPUTATION OF MARGINAL FREQUENCIES OF THE ROWS * C * AND COLUMNS, TOTAL NUMBER OF OBJECTS AND NUMBER OF * C * MISSING VALUES PER VARIABLE * C * * C ****************************************************************** DIMENSION DATA(NOBJT,NVAR),DC(KJ),DCINV(KJ),ICUMCA(NVAR),DR(NOBJ), X NCAT(NVAR),DRINV(NOBJ),NMIS(NVAR) INTEGER DATA SRC=0.0 DO 10 K=1,KJ DC(K)=0.0 DCINV(K)=0.0 10 CONTINUE DO 20 J=1,NVAR NMIS(J)=0 DO 30 I=1,NOBJ IF(DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) NMIS(J)=NMIS(J)+1 IF(DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) GO TO 30 IND=DATA(I,J)+ICUMCA(J) DC(IND)=DC(IND)+1.0 30 CONTINUE 20 CONTINUE DO 40 K=1,KJ IF (DC(K).NE.0) DCINV(K)=1.0/DC(K) 40 CONTINUE DO 50 I=1,NOBJ DR(I)=0.0 DRINV(I)=0.0 DO 60 J=1,NVAR IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) GO TO 60 DR(I)=DR(I)+1.0 60 CONTINUE IF (DR(I).EQ.0) GO TO 50 SRC=SRC+DR(I) DRINV(I)=1.0/DR(I) 50 CONTINUE IF (MISSOL.EQ.1) SRC=FLOAT(NOBJ)*FLOAT(NVAR) IF (SRC.LT. (FLOAT(NOBJ)*FLOAT(NVAR)) ) MIS=1 IF (SRC.EQ. (FLOAT(NOBJ)*FLOAT(NVAR)) ) MIS=0 RETURN END SUBROUTINE IMTQL2(NM,N,D,E,Z,IERR) C ****************************************************************** C * * C * I M T Q L 2 * C * * C * PURPOSE: DETERMINES THE EIGENVALUES AND EIGENVECTORS OF A * C * SYMMETRIC TRIDIAGONAL MATRIX USING THE IMPLICIT * C * QL METHOD. * C * * C * SUBROUTINES CALLED: NONE * C * * C * ADAPTED FROM EISPACK * C * * C ****************************************************************** DIMENSION D(N),E(N),Z(NM,N) REAL MACHEP C C ********** MACHEP IS A MACHINE DEPENDENT PARAMETER SPCIFYING C THE RELATIVE PRECISION OF FLOATING POINT ARITHMETIC. C C ********** C MACHEP = 2.0**(-20) C IERR = 0 IF (N .EQ. 1) GO TO 1001 C DO 100 I = 2, N 100 E(I-1) = E(I) C E(N) = 0.0 C DO 240 L = 1, N J = 0 C ********** LOOK FOR SMALL SUB-DIAGONAL ELEMENT ********** 105 DO 110 M = L, N IF (M .EQ. N) GO TO 120 IF (ABS(E(M)) .LE. MACHEP* (ABS(D(M)) + ABS(D(M+1)))) X GO TO 120 110 CONTINUE C 120 P = D(L) IF (M .EQ. L) GO TO 240 IF (J .EQ.30) GO TO 1000 J = J + 1 C ********** FORM SHIFT ********** G = (D(L+1) - P) / (2.0 * E(L)) R = SQRT(G*G+1.0) G = D(M) - P + E(L) / (G + SIGN(R,G)) S = 1.0 C = 1.0 P = 0.0 NML = M - L C ********** FOR I=M-1 STEP -1 UNTIL L DO -- ********** DO 200 II = 1, NML I = M - II F = S * E(I) B = C * E(I) IF (ABS(F) .LT. ABS(G)) GO TO 150 C = G / F R = SQRT(C*C+1.0) E(I+1) = F * R S = 1.0 / R C = C * S GO TO 160 150 S = F / G R = SQRT(S*S+1.0) E(I+1) = G* R C = 1.0 / R S = S * C 160 G = D(I+1) - P R = (D(I) - G) * S + 2.0 * C * B P = S * R D(I+1) = G + P G = C * R - B C ********** FORM VECTOR ********** DO 180 K = 1, N F = Z(K,I+1) Z(K,I+1) = S * Z(K,I) + C * F Z(K,I) = C * Z(K,I) - S * F 180 CONTINUE C 200 CONTINUE C D(L) = D(L) - P E(L) = G E(M) = 0.0 GO TO 105 240 CONTINUE C ********** ORDER EIGENVALUES AND EIGENVECTORS ********** DO 300 II = 2, N I = II - 1 K = I P = D(I) C DO 260 J = II, N IF (D(J) .LE. P ) GO TO 260 K = J P = D(J) 260 CONTINUE C IF (K .EQ. I) GO TO 300 D(K) = D(I) D(I) = P C DO 280 J = 1, N P = Z(J,I) Z(J,I) = Z(J,K) Z(J,K) = P 280 CONTINUE C 300 CONTINUE C GO TO 1001 C ********** SET ERROR -- NO CONVERGENCE TO AN C EIGENVALUE AFTER 30 ITERATIONS ********** 1000 IERR = L 1001 RETURN END SUBROUTINE INITIA(NN,NVAR,NOBJ,A,W,SRC,FIT,NORM,MIS,IEIG,ICUMCA, X NCAT) C ****************************************************************** C * * C * I N I T I A * C * * C * PURPOSE:INITIAL VALUES FOR OBJECT SCORES OR * C * CATEGORY QUANTIFICATIONS * C * * C ****************************************************************** DIMENSION A(NN),W(NN),ICUMCA(NVAR),NCAT(NVAR) K=15000 DO 10 I=1,NN IND=1.+ FLOAT(NN)*RANDOM(K) A(I)=IND 10 CONTINUE C IF (NORM.EQ.1.AND.MIS.EQ.0) CALL CENORA(A,NN,SRC,FIT,IEIG) IF (NORM.EQ.1.AND.MIS.EQ.1) X CALL CENORB(A,NN,NVAR,NOBJ,W,SRC,FIT,IEIG,NORM) IF (NORM.EQ.2.AND.MIS.EQ.1) X CALL CENORB(A,NN,NVAR,NOBJ,W,SRC,FIT,IEIG,NORM) IF (NORM.EQ.2.AND.MIS.EQ.0) X CALL CENORC(A,NN,NVAR,NOBJ,W,SRC,FIT,IEIG,ICUMCA,NCAT) RETURN END C FUNCTION RANDOM(N) M=2796203 N=125*N I=N/M N=N-I*M RANDOM=FLOAT(N)/FLOAT(M) RETURN END SUBROUTINE ITERA(XROW,YCOL,DATA,NOBJ,NOBJT,NVAR,DR,DC,DRINV,DCINV, X SRC,KJ,ICUMCA,FIT,NORM,IEIG,NCAT,MIS,WORK,MAXCAT,IND) C ****************************************************************** C * * C * I T E R A * C * * C * PURPOSE: CALL OF SELECTED SUBROUTINES IN ITERATION-PROCES * C * DEPENDING ON NORMALIZATION AND OCCURENCE OF * C * MISSING DATA * C * SUBROUTINES CALLED:TRANXA,TRANYA,CENORA,CENORB * C * * C ****************************************************************** DIMENSION DATA(NOBJT,NVAR),XROW(NOBJ),YCOL(KJ),DR(NOBJ), X DC(KJ),DRINV(NOBJ),DCINV(KJ),ICUMCA(NVAR),NCAT(NVAR), X WORK(MAXCAT),IND(NVAR) IF(NORM.EQ.1) CALL TRANYA(YCOL,XROW,DATA,NOBJ,NOBJT,NVAR,KJ,DCINV, X ICUMCA,NCAT,IEIG,NORM,WORK,MAXCAT) CALL TRANXA(XROW,YCOL,DATA,NOBJ,NOBJ,NOBJT,NVAR,KJ,ICUMCA, X NCAT,DRINV,MIS,IEIG,NORM) IF(NORM.EQ.2) CALL TRANYA(YCOL,XROW,DATA,NOBJ,NOBJT,NVAR,KJ,DCINV, X ICUMCA,NCAT,IEIG,NORM,WORK,MAXCAT) IF (NORM.EQ.1.AND.MIS.EQ.0) X CALL CENORA (XROW,NOBJ,SRC,FIT,IEIG) IF (NORM.EQ.1.AND.MIS.EQ.1) X CALL CENORB (XROW,NOBJ,NVAR,NOBJ,DR,SRC,FIT,IEIG,NORM) IF (NORM.EQ.2.AND.MIS.EQ.1) X CALL CENORB (YCOL, KJ, NVAR,NOBJ,DC,SRC,FIT,IEIG,NORM) IF (NORM.EQ.2.AND.MIS.EQ.0) X CALL CENORC (YCOL, KJ, NVAR,NOBJ,DC,SRC,FIT,IEIG,ICUMCA,NCAT) RETURN END SUBROUTINE OUTPUT(K,FIT,DIF) C ****************************************************************** C * * C * O U T P U T * C * * C * PURPOSE: PRINTING OF ITERATION HISTORY * C * * C ****************************************************************** COMMON/CP2/IWRITE,INPARA IF (K.GT.1) GO TO 100 WRITE (IWRITE,605) WRITE (IWRITE,610) WRITE (IWRITE,615) WRITE (IWRITE,620) WRITE (IWRITE,625) K,FIT GO TO 700 100 WRITE (IWRITE,625) K,FIT,DIF 605 FORMAT (1H1//,5X,17HITERATION HISTORY) 610 FORMAT (1H ,4X,17H*****************) 615 FORMAT (1H0//,49X,15HDIFFERENCE WITH) 620 FORMAT (1H ,21X,46HITERATION FIT PRECEDING ITERATION) 625 FORMAT (1H0,26X,I4,1X,F11.6,3X,F11.6) 700 RETURN END SUBROUTINE PCA(DATA,NOBJ,NVAR,A,B,C,D,WORK1,WORK2,CM,ICUMCA,NCAT, X XROW,YCOL,KJ,CMEAN,SSQ,KRAW,FIT,NORM,IPRI,IPLO, X IWRI,LUN,IPCAB,MISEST) C ******************************************************************* C * * C * P C A * C * * C * PURPOSE: PRINCIPAL COMPONENT SOLUTION FOR INDIVIDUAL SCORES * C * AND VARIABLES * C * * C ******************************************************************* DIMENSION DATA(NOBJ,NVAR),WORK1(NVAR),WORK2(NVAR),A(NOBJ),B(NOBJ), X C(KJ),D(KJ),CM(NVAR,NVAR),ICUMCA(NVAR),NCAT(NVAR), X XROW(NOBJ),YCOL(KJ),CMEAN(NVAR),SSQ(NVAR),IPRI(10), X IPLO(11),IWRI(10),LUN(5) COMMON/CP2/IWRITE,INPARA INTEGER DATA,ALL REAL MAX1,MAX2 DATA S/1H*/ C IF (IPCAB.EQ.0) GO TO 100 C C ** IPCAB INDICATES PRINCIPAL COMPONENTS SOLUTION BEFORE TRANSFORMATION C C ** ARRAY WORK1 CONTAINS EIGENVALUES, ARRAY CM CONTAINS EIGENVECTORS, ** C ** ARRAY WORK2 IS USED FOR ROWS OF RAW DATA C ICOSC=IPRI(8)+IPLO(8)+IWRI(8) C C ** ICOSC INDICATES THE COMPUTATION OF COMPONENT SCORES B.T. C IF (ICOSC.EQ.0) GO TO 500 DO 10 I=1,NOBJ A(I)=0.0 B(I)=0.0 DO 20 J=1,NVAR WORK2(J)=DATA(I,J) IF(WORK2(J).EQ.0.0.OR.WORK2(J).GT.NCAT(J)) WORK2(J)=CMEAN(J) K=1 WORK2(J)=WORK2(J)-CMEAN(J) IF (SSQ(J).EQ.0.0) GO TO 16 WORK2(J)=WORK2(J)/SQRT (SSQ(J)) 16 A(I)=A(I)+WORK2(J)*CM(J,K) K=2 B(I)=B(I)+WORK2(J)*CM(J,K) 20 CONTINUE IF (WORK1(1).EQ.0.0) GO TO 15 A(I)=A(I)*SQRT(FLOAT(NOBJ))/SQRT (WORK1(1)) 15 IF (WORK1(2).EQ.0.0) GO TO 10 B(I)=B(I)*SQRT(FLOAT(NOBJ))/SQRT (WORK1(2)) 10 CONTINUE C C ** COMPUTE COMPONENT LOADINGS BEFORE TRANSFORMATION C 500 DO 40 J=1,NVAR K=1 C(J)=CM(J,K)*SQRT (WORK1(K)) K=2 D(J)=CM(J,K)*SQRT (WORK1(K)) 40 CONTINUE GO TO 700 C C ** COMPUTE OBJECT-SCORES AND COMPONENT-LOADINGS FOR TRANSFORMED DATA** C C ARRAY WORK2 IS USED TO REPLACE CATEGORY-LABELS BY QUANTIFICATIONS. ** C 100 ICOSC=IPRI(10)+IPLO(10)+IWRI(10) IF (ICOSC.EQ.0) GO TO 510 IF (MISEST.NE.2) GO TO 501 C C** IF MISEST.EQ.2 THE AVERAGES OF INDIVIDUALS WITH MISSING DATA ARE C COMPUTED (ARRAY C) AND SUBSTITUTED (ARRAY WORK2) C DO 101 J=1,NVAR SUM=0.0 KTEL=0 DO 102 I=1,NOBJ IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) SUM=SUM+XROW(I) IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) KTEL=KTEL+1 102 CONTINUE IF (KTEL.EQ.0) GO TO 101 AVE=SUM/FLOAT (KTEL) C(J)=AVE 101 CONTINUE C C ** COMPUTE COMPONENT SCORES AFTER TRANSFORMATION (ARRAYS A AND B) C 501 DO 110 I=1,NOBJ A(I)=0.0 B(I)=0.0 DO 120 J=1,NVAR MMM=0 IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) MMM=1 IF (MMM.NE.1) GO TO 149 IF (MISEST.EQ.0) WORK2(J)=XROW(I) IF (MISEST.EQ.1) WORK2(J)=0.0 IF (MISEST.EQ.2) WORK2(J)=C(J) IF (MISEST.EQ.0.AND.NORM.EQ.2) WORK2(J)=WORK2(J)/FIT IF (MMM.EQ.1) GO TO 150 149 IND=DATA(I,J)+ICUMCA(J) WORK2(J)=YCOL(IND) 150 K=1 WORK2(J)=WORK2(J)-CMEAN(J) IF (SSQ(J).EQ.0.0) GO TO 116 WORK2(J)=WORK2(J)/SQRT (SSQ(J)) 116 A(I)=A(I)+WORK2(J)*CM(J,K) K=2 B(I)=B(I)+WORK2(J)*CM(J,K) 120 CONTINUE IF (WORK1(1).EQ.0.0) GO TO 115 A(I)=A(I)*SQRT(FLOAT(NOBJ))/SQRT (WORK1(1)) 115 IF (WORK1(2).EQ.0.0) GO TO 110 B(I)=B(I)*SQRT(FLOAT(NOBJ))/SQRT (WORK1(2)) 110 CONTINUE C C ** COMPUTE COMPONENT LOADINGS AFTER TRANSFORMATION (ARRAYS C AND D) C 510 DO 140 J=1,NVAR K=1 C(J)=CM(J,K)*SQRT (WORK1(K)) K=2 D(J)=CM(J,K)*SQRT (WORK1(K)) 140 CONTINUE C 700 MAX1=C(1) MAX2=D(1) DO 160 J=2,NVAR IF (ABS( MAX1 ).LT.ABS( C(J) )) MAX1=C(J) IF (ABS( MAX2 ).LT.ABS( D(J) )) MAX2=D(J) 160 CONTINUE IF (MAX1.GT.0.0) GO TO 850 IF (ICOSC.EQ.0) GO TO 171 DO 170 I=1,NOBJ A(I)=A(I)*(-1) 170 CONTINUE 171 DO 180 J=1,NVAR C(J)=C(J)*(-1) 180 CONTINUE C 850 IF (MAX2.GT.0.0) GO TO 900 IF (ICOSC.EQ.0) GO TO 191 DO 190 I=1,NOBJ B(I)=B(I)*(-1) 190 CONTINUE 191 DO 200 J=1,NVAR D(J)=D(J)*(-1) 200 CONTINUE C 900 IF (IPCAB.EQ.0) GO TO 1900 IF (IPRI(8).EQ.0) GO TO 1000 IF (IPRI(5).EQ.0) WRITE (6,610) IF (IPRI(5).NE.0) WRITE (6,600) WRITE (6,601) (K,K=1,2) KK=18 WRITE (6,602) (S,J=1,KK) DO 46 I=1,NOBJ WRITE (6,605) I,A(I),B(I) 46 CONTINUE 1000 IF (IPRI(7).EQ.0) GO TO 1100 IF (IPRI(5).EQ.0.AND.IPRI(8).EQ.0) WRITE (6,613) IF (IPRI(5).NE.0.OR.IPRI(8).NE.0) WRITE (6,603) WRITE (6,601) (K,K=1,2) KK=18 WRITE (6,602) (S,J=1,KK) DO 50 J=1,NVAR WRITE (6,605) J,C(J),D(J) 50 CONTINUE 1100 IF (IWRI(8).EQ.1) CALL WRITE2(LUN(5),A,B,NOBJ) IF (IWRI(7).EQ.1) CALL WRITE2(LUN(5),C,D,NVAR) GO TO 3000 C 1900 IF (IPRI(10).EQ.0) GO TO 2000 IF (IPRI(6).EQ.0) WRITE (6,620) IF (IPRI(6).NE.0) WRITE (6,600) WRITE (6,601) (K,K=1,2) KK=18 WRITE (6,602) (S,J=1,KK) DO 146 I=1,NOBJ WRITE (6,605) I,A(I),B(I) 146 CONTINUE 2000 IF (IPRI(9).EQ.0) GO TO 2100 IF (IPRI(6).EQ.0.AND.IPRI(10).EQ.0) WRITE (6,623) IF (IPRI(6).NE.0.OR.IPRI(10).NE.0) WRITE (6,603) WRITE (6,601) (K,K=1,2) KK=18 WRITE (6,602) (S,J=1,KK) DO 151 J=1,NVAR WRITE (6,605) J,C(J),D(J) 151 CONTINUE 2100 IF (IWRI(10).EQ.1) CALL WRITE2(LUN(5),A,B,NOBJ) IF (IWRI(9).EQ.1) CALL WRITE2(LUN(5),C,D,NVAR) 600 FORMAT (1H0//,5X,16HCOMPONENT SCORES) 610 FORMAT (1H0//,5X,38HCOMPONENT SCORES BEFORE TRANSFORMATION) 620 FORMAT (1H0//,5X,37HCOMPONENT SCORES AFTER TRANSFORMATION) 601 FORMAT (1H0,6X,2I8) 602 FORMAT (1H ,7X,18A1) 603 FORMAT (1H0//,5X,18HCOMPONENT LOADINGS) 613 FORMAT (1H0//,5X,40HCOMPONENT LOADINGS BEFORE TRANSFORMATION) 623 FORMAT (1H0//,5X,39HCOMPONENT LOADINGS AFTER TRANSFORMATION) 605 FORMAT (1H ,1X,I5,1X,1H*,2F8.3) 3000 RETURN END SUBROUTINE PREPA1 (DATA,NOBJ,NVAR,DR,DC,DRI,DCI,ICUMCA,NCAT,NP,J, X X,Y,KJ) C ****************************************************************** C * * C * P R E P A 1 * C * * C * PURPOSE: PREPARE INPUT TO PLOT QUANTIFICATIONS VERSUS ORI- * C * GINAL CATEGORY LABELS, TOGETHER WITH INDIVIDUAL * C * SCORES * C * * C ****************************************************************** DIMENSION DATA(NOBJ,NVAR),DR(NOBJ),DC(NCAT),ICUMCA(NVAR), X DRI(NOBJ),DCI(NCAT),X(NOBJ),Y(KJ) INTEGER DATA DO 20 I=1,NOBJ DR(I)=DATA(I,J) IF (DR(I).GT.NCAT) DR(I)=0.0 20 CONTINUE DO 30 K=1,NCAT DC(K)=FLOAT(K) 30 CONTINUE NP=NOBJ+NCAT DO 40 I=1,NOBJ DRI(I)=X(I) 40 CONTINUE KB=ICUMCA(J)+1 KE=ICUMCA(J)+NCAT DO 50 K=KB,KE IND=K-ICUMCA(J) DCI(IND)=Y(K) 50 CONTINUE RETURN END SUBROUTINE PRIINT(X,M,N,ALL,MISOB,DR) C ****************************************************************** C * * C * P R I I N T * C * * C * PURPOSE: PRINTS A M*N ARRAY WITH ROW- AND COLUMN IDENTIFI- * C * CATION, EACH LINE CONTAINING AT MOST 25 VALUES (I4) * C * * C * SUBROUTINES CALLED: NONE * C * * C ****************************************************************** DIMENSION X(M,N),DR(M) COMMON/CP2/IWRITE,INPARA INTEGER X,ALL DATA A/1H*/ C C C CHECK HOW MANY TIMES THE NUMBER OF COLUMNS EXCEEDS 25 C C C MM=M IF (ALL.EQ.0) MM=10 IF (M.LT.10) MM=M IF (ALL.NE.0) MM=M NTRU=N/25 NREM=N-NTRU*25 IF (NREM.GT.0) GO TO 1 NREM=25 NTRU=NTRU-1 1 NREP=NTRU+1 C C C PRINT NREP TIMES A BLOCK OF ELEMENTS C C C NK=25 KB=1 DO 2 I=1,NREP IF (I.EQ.NREP) NK=NREM KE=KB+NK-1 WRITE (IWRITE,621) IF (MISOB.EQ.0) WRITE (IWRITE,622) (K,K=KB,KE) IF (MISOB.EQ.1) WRITE (IWRITE,620) (K,K=KB,KE) KK=NK*4+6 WRITE (IWRITE,623) (A,J=1,KK) DO 3 L=1,MM IF (MISOB.EQ.1) MV=N-DR(L) IF (MISOB.EQ.0) WRITE (IWRITE,624) (X(L,K),K=KB,KE) 3 IF (MISOB.EQ.1) WRITE (IWRITE,625) L,MV,(X(L,K),K=KB,KE) KB=KB+NK 2 CONTINUE 621 FORMAT (1H0/,15X,9HVARIABLES) 622 FORMAT(1H0,19X,25I4) 620 FORMAT(1H0,18X,1HM,25I4) 623 FORMAT (1H ,14X,106A1) 624 FORMAT(1H ,14X,1H*,4X,25I4) 625 FORMAT(1H ,8X,I5,1X,1H*,26I4) RETURN END SUBROUTINE PRINTA (TITLE,DATA,NOBJ,NOBJA,NVAR,NVARA,KJ,DC,DR,NORM, X NITER,CRIT,FORMAT,ICUMCA,NCAT,MAXCAT,NMIS,IDC, X IPRI,IPLO,IWRI,LUN,ISUMCA,IEIG,MISEST,MISSOL) C ******************************************************************* C * * C * P R I N T A * C * * C * PURPOSE: PRINTING OF DECK-SETUP AND DATA-MATRIX * C * SUBROUTINES CALLED: PRIINT * C * * C ******************************************************************* DIMENSION DATA(NOBJ,NVAR),TITLE(20),FORMAT(20),ICUMCA(NVAR), X DC(KJ),IDC(KJ),DR(NOBJ),NCAT(NVAR),NMIS(NVAR),IPRI(10), X IPLO(11),IWRI(10),LUN(5) COMMON/CP2/IWRITE,INPARA INTEGER DATA DATA A/1H*/,AA/1H-/ INT=4.+(FLOAT(NVAR))/16 C WRITE (IWRITE,600) WRITE (IWRITE,601) WRITE (IWRITE,602) WRITE (IWRITE,6602) WRITE (IWRITE,603) WRITE (IWRITE,604) WRITE (IWRITE,6604) WRITE (IWRITE,605) TITLE WRITE (IWRITE,606) WRITE (IWRITE,607) WRITE (IWRITE,610) NOBJ NOPAS=NOBJ-NOBJA IF (NOPAS.EQ.0.OR.NORM.EQ.2) WRITE (IWRITE,611) NOBJA IF (NOPAS.NE.0.AND.NORM.EQ.1) WRITE (IWRITE,612) NOBJ,A WRITE (IWRITE,615) NVAR NVPAS=NVAR-NVARA IF (NVPAS.EQ.0.OR.NORM.EQ.1) WRITE (IWRITE,614) NVARA IF (NVPAS.NE.0.AND.NORM.EQ.2) WRITE (IWRITE,613) NVAR,A WRITE (IWRITE,616) MAXCAT WRITE (IWRITE,617) ISUMCA IF (NOPAS.NE.0.AND.NORM.EQ.1) WRITE (IWRITE,712) A IF (NVPAS.NE.0.AND.NORM.EQ.2) WRITE (IWRITE,713) A WRITE (IWRITE,618) WRITE (IWRITE,619) WRITE (IWRITE,630) NORM WRITE (IWRITE,631) WRITE (IWRITE,632) WRITE (IWRITE,650) IEIG WRITE (IWRITE,651) WRITE (IWRITE,652) WRITE (IWRITE,633) NITER WRITE (IWRITE,634) CRIT WRITE (IWRITE,839) MISSOL WRITE (IWRITE,840) WRITE (IWRITE,841) WRITE (IWRITE,835) MISEST WRITE (IWRITE,836) WRITE (IWRITE,837) WRITE (IWRITE,838) WRITE (IWRITE,655) WRITE (IWRITE,656) WRITE (IWRITE,635) LUN(1),FORMAT IWR=IWRI(2)+IWRI(3) IF (LUN(2).NE.0.AND.IWR.NE.0) WRITE (IWRITE,636) LUN(2) IF (LUN(2).EQ.0.OR. IWR.EQ.0) WRITE (IWRITE,736) AA IF (LUN(3).NE.0.AND.IWRI(4).NE.0) WRITE (IWRITE,637) LUN(3) IF (LUN(3).EQ.0.OR. IWRI(4).EQ.0) WRITE (IWRITE,737) AA IWR=IWRI(5)+IWRI(6) IF (LUN(4).NE.0.AND.IWR.NE.0) WRITE (IWRITE,638) LUN(4) IF (LUN(4).EQ.0.OR. IWR.EQ.0) WRITE (IWRITE,738) AA IWR=IWRI(7)+IWRI(8)+IWRI(9)+IWRI(10) IF (LUN(5).NE.0.AND.IWR.NE.0) WRITE (IWRITE,639) LUN(5) IF (LUN(5).EQ.0.OR. IWR.EQ.0) WRITE (IWRITE,739) AA WRITE (IWRITE,660) WRITE (IWRITE,661) IPRI(1),IPLO(1),IWRI(1) WRITE (IWRITE,662) IPRI(2),IPLO(2),IWRI(2) WRITE (IWRITE,663) IPRI(3),IPLO(3),IWRI(3) WRITE (IWRITE,664) IPRI(4),IPLO(4),IWRI(4) WRITE (IWRITE,665) IPRI(5),IPLO(5),IWRI(5) WRITE (IWRITE,666) IPRI(6),IPLO(6),IWRI(6) WRITE (IWRITE,667) IPRI(7),IPLO(7),IWRI(7) WRITE (IWRITE,668) IPRI(8),IPLO(8),IWRI(8) WRITE (IWRITE,669) IPRI(9),IPLO(9),IWRI(9) WRITE (IWRITE,670) IPRI(10),IPLO(10),IWRI(10) WRITE (IWRITE,671) IPLO(11) WRITE (IWRITE,1100) INT WRITE (IWRITE,1101) NOBJ,NOBJA,NVAR,NVARA,MAXCAT,ISUMCA WRITE (IWRITE,1102) (NCAT(J),J=1,NVAR) WRITE (IWRITE,1103) NORM,IEIG,NITER,CRIT,MISSOL,MISEST IF (KJ.GT.ISUMCA) RETURN IF (IPRI(1).EQ.0) N=10 IF (NOBJ.LT.10) N=NOBJ IF (IPRI(1).NE.0) N=NOBJ WRITE (IWRITE,620) N CALL PRIINT(DATA,NOBJ,NVAR,IPRI(1),1,DR) WRITE (IWRITE,621) CALL PRIINT(NCAT,1,NVAR,1,0,DR) DO 20 J=1,KJ IDC(J)=DC(J) 20 CONTINUE WRITE (IWRITE,640) WRITE (IWRITE,641) IF (MAXCAT.GT.18) MM=18 IF (MAXCAT.LE.18) MM=MAXCAT WRITE (IWRITE,642) (J,J=1,MM) KK=MM*6+8 WRITE (IWRITE,643) (A,J=1,KK) DO 50 J=1,NVAR K=ICUMCA(J)+1 L=NCAT(J)+ICUMCA(J) WRITE (IWRITE,644) J,NMIS(J),(IDC(I),I=K,L) 50 CONTINUE 600 FORMAT (1H1//,5X,17H*****************,76X,11HVERSION 3.0, X 23H RELEASED DECEMBER 1984) 601 FORMAT (1H0,4X,17H* P R I M A L S *, X 20X,36HONE-DIMENSIONAL HOMOGENEITY ANALYSIS, X 20X,32HJACQUELINE MEULMAN & ALBERT GIFI) 602 FORMAT (1H0,4X,17H*****************,76X,19HDEPT. OF DATATHEORY) 6602 FORMAT (1H0,97X,20HUNIVERSITY OF LEIDEN) 603 FORMAT (1H0,97X,17HMIDDELSTEGRACHT 4) 604 FORMAT (1H0,97X,14H2312 TW LEIDEN) 6604 FORMAT (1H0,97X,15HTHE NETHERLANDS) 605 FORMAT (1H0,4X,8HTITLE : ,20A4) 606 FORMAT (1H0//,5X,19HDATA SPECIFICATIONS) 607 FORMAT (1H ,4X,19H*******************) 610 FORMAT (1H0,4X,32HNUMBER OF OBJECTS OR INDIVIDUALS,12X,I5) 611 FORMAT (1H0/,5X,40HNUMBER OF OBJECTS ACTIVE IN THE ANALYSIS, X 4X,I5) 612 FORMAT (1H0/,5X,40HNUMBER OF OBJECTS ACTIVE IN THE ANALYSIS, X 4X,I5,1X,1A1) 712 FORMAT (1H0/,5X,1A1,1X,39HBECAUSE THE CATEGORIES HAD TO BE IN THE, X 54H CENTROID OF THE OBJECTS, ALL OBJECTS HAVE BEEN FORCED, X 23H ACTIVE IN THE ANALYSIS) 615 FORMAT (1H0/,5X,19HNUMBER OF VARIABLES,25X,I5) 614 FORMAT (1H0/,5X,42HNUMBER OF VARIABLES ACTIVE IN THE ANALYSIS, X 2X,I5) 613 FORMAT (1H0/,5X,42HNUMBER OF VARIABLES ACTIVE IN THE ANALYSIS, X 2X,I5,1X,1A1) 713 FORMAT (1H0/,5X,1A1,1X,36HBECAUSE THE OBJECTS HAD TO BE IN THE, X 59H CENTROID OF THE CATEGORIES, ALL VARIABLES HAVE BEEN FORCED, X 23H ACTIVE IN THE ANALYSIS) 616 FORMAT (1H0/,5X,28HMAXIMUM NUMBER OF CATEGORIES,16X,I5) 617 FORMAT (1H0/,5X,26HTOTAL NUMBER OF CATEGORIES,18X,I5) 618 FORMAT (1H0//,5X,23HANALYSIS SPECIFICATIONS) 619 FORMAT (1H ,4X,23H***********************) 620 FORMAT (1H1//,5X,10HLISTING OF,I5,23H ROWS OF THE DATAMATRIX) 621 FORMAT (1H0//,5X,20HNUMBER OF CATEGORIES) 630 FORMAT (1H0,4X,13HNORMALIZATION,31X,I5) 631 FORMAT (1H0,6X,45H1 : CATEGORIES ARE IN THE CENTROID OF OBJECTS) 632 FORMAT (1H0,6X,45H2 : OBJECTS ARE IN THE CENTROID OF CATEGORIES) 633 FORMAT (1H0/,5X,28HMAXIMUM NUMBER OF ITERATIONS,16X,I5) 634 FORMAT (1H0/,5X,21HCONVERGENCE CRITERION,18X,F10.8) 655 FORMAT (1H0//,5X,27HINPUT/OUTPUT SPECIFICATIONS) 656 FORMAT (1H ,4X,27H***************************) 635 FORMAT ((1H0,5X,35H1 THE DATA HAVE BEEN READ FROM UNIT,I3, X 13H WITH FORMAT ,19A4)/(49X,16A4)) 636 FORMAT (1H0,5X,47H2 THE OBJ.SCORES/CAT.QUANTIFICATIONS HAVE BEEN, X 16H WRITTEN TO UNIT,I3,24H WITH FORMAT (I10,F12.7)) 736 FORMAT (1H0,5X,47H2 THE OBJ.SCORES/CAT.QUANTIFICATIONS HAVE BEEN, X 16H WRITTEN TO UNIT,2X,1A1) 637 FORMAT (1H0,5X,28H3 THE TRANSFORMED DATAMATRIX,11X,8HHAS BEEN, X 16H WRITTEN TO UNIT,I3,21H WITH FORMAT (10F8.4)) 737 FORMAT (1H0,5X,28H3 THE TRANSFORMED DATAMATRIX,11X,8HHAS BEEN, X 16H WRITTEN TO UNIT,2X,1A1) 638 FORMAT (1H0,5X,35H4 THE COR.MATRIX B/A TRANSFORMATION,4X, X 8HHAS BEEN, X 16H WRITTEN TO UNIT,I3,21H WITH FORMAT (10F8.4)) 738 FORMAT (1H0,5X,35H4 THE COR.MATRIX B/A TRANSFORMATION,4X, X 8HHAS BEEN, X 16H WRITTEN TO UNIT,2X,1A1) 639 FORMAT (1H0,5X,47H5 THE PCA SOLUTION B/A TRANSFORMATION HAS BEEN, X 16H WRITTEN TO UNIT,I3,25H WITH FORMAT (I10,2F12.7)) 739 FORMAT (1H0,5X,47H5 THE PCA SOLUTION B/A TRANSFORMATION HAS BEEN, X 16H WRITTEN TO UNIT,2X,1A1) 660 FORMAT (1H0//,5X,28HOUTPUT OPTIONS 0:NO, 1:YES,18X,5HPRINT,2X, X 4HPLOT,1X,5HWRITE) 661 FORMAT (1H0/,5X,13H 1 DATAMATRIX,31X,I5,1X,I5,1X,I5) 662 FORMAT (1H0,4X,16H 2 OBJECT SCORES,28X,I5,1X,I5,1X,I5) 663 FORMAT (1H0,4X,27H 3 CATEGORY QUANTIFICATIONS,17X,I5,1X,I5,1X,I5) 664 FORMAT (1H0,4X,25H 4 TRANSFORMED DATAMATRIX,19X,I5,1X,I5,1X,I5) 665 FORMAT (1H0,4X,43H 5 CORRELATION MATRIX BEFORE TRANSFORMATION, X 1X,I5,1X,I5,1X,I5) 666 FORMAT (1H0,4X,43H 6 CORRELATION MATRIX AFTER TRANSFORMATION, X 1X,I5,1X,I5,1X,I5) 667 FORMAT (1H0,4X,43H 7 COMPONENT LOADINGS BEFORE TRANSFORMATION, X 1X,I5,1X,I5,1X,I5) 668 FORMAT (1H0,4X,43H 8 COMPONENT SCORES BEFORE TRANSFORMATION, X 1X,I5,1X,I5,1X,I5) 669 FORMAT (1H0,4X,43H 9 COMPONENT LOADINGS AFTER TRANSFORMATION, X 1X,I5,1X,I5,1X,I5) 670 FORMAT (1H0,4X,43H10 COMPONENT SCORES AFTER TRANSFORMATION, X 1X,I5,1X,I5,1X,I5) 671 FORMAT (1H0,4X,43H11 COMPONENT LOADINGS B & A TRANSFORMATION, X 7X,I5) 640 FORMAT (1H0//,5X,20HMARGINAL FREQUENCIES) 641 FORMAT (1H0/,15X,10HCATEGORIES) 642 FORMAT (1H0,20X,1HM,18I6) 643 FORMAT (1H ,4X,9HVARIABLES,1X,116A1) 644 FORMAT (1H ,8X,I5,1X,1H*,19I6/(15X,1H*,6X,18I6)) 650 FORMAT (1H0/,5X,10HEIGENVALUE,34X,I5) 651 FORMAT (1H0,6X,33H0 : MAXIMIZING LARGEST EIGENVALUE) 652 FORMAT (1H0,6X,34H1 : MINIMIZING SMALLEST EIGENVALUE) 835 FORMAT (1H0/,5X,29HSUBSTITUTION FOR MISSING DATA,15X,I5) 836 FORMAT (1H0,6X,17H0 : OBJECT SCORES) 837 FORMAT (1H0,6X,7H1 : 0.0) 838 FORMAT (1H0,6X,41H2 : CENTROID OF OBJECTS WITH MISSING DATA) 839 FORMAT (1H0/,5X,25HTREATMENT OF MISSING DATA,19X,I5) 840 FORMAT (1H0,6X,45H0 : MISSING DATA PASSIVE CATEGORIES (DELETED)) 841 FORMAT (1H0,6X,36H1 : MISSING DATA MULTIPLE CATEGORIES) 1100 FORMAT (1H0/////,5X,28HECHO OF PARAMETER CARDS (2 -,I3,1H)) 1101 FORMAT (1H0,4X,6I5) 1102 FORMAT (1H0,4X,16I5) 1103 FORMAT (1H0,4X,3I5,F10.8,2I5) RETURN END SUBROUTINE PRIREA(X,M,N,ALL) C ****************************************************************** C * * C * P R I R E A * C * * C * PURPOSE: PRINTS A M*N ARRAY WITH ROW- AND COLUMN IDENTIFI- * C * CATION, EACH LINE CONTAINING AT MOST 15 VALUES (F8.3)* C * * C * SUBROUTINES CALLED: NONE * C * * C ****************************************************************** DIMENSION X(M,N) COMMON/CP2/IWRITE,INPARA DATA A/1H*/ INTEGER ALL C C C CHECK HOW MANY TIMES THE NUMBER OF COLUMNS EXCEEDS 15 C C C NN=N IF (ALL.EQ.0) NN=10 IF (N.LT.10) NN=N IF (ALL.NE.0) NN=N NTRU=NN/15 NREM=NN-NTRU*15 IF (NREM.GT.0) GO TO 1 NREM=15 NTRU=NTRU-1 1 NREP=NTRU+1 C C C PRINT NREP TIMES A BLOCK OF ELEMENTS C C C NK=15 KB=1 DO 2 I=1,NREP IF (I.EQ.NREP) NK=NREM KE=KB+NK-1 C WRITE(IWRITE,621) WRITE(IWRITE,60) (K,K=KB,KE) KK=NK*8+2 WRITE(IWRITE,61)(A,J=1,KK) DO 3 L=1,M 3 WRITE(IWRITE,62) L,(X(L,K),K=KB,KE) KB=KB+NK 2 CONTINUE RETURN C 621 FORMAT(1H0/,17X,9HVARIABLES) 60 FORMAT(1H0,6X,15I8) 61 FORMAT(1H ,4X,3X,122A1) 62 FORMAT(1H ,1X,I5,1X,1H*,15F8.3) END SUBROUTINE PRPLOT(X,Y,IND,NITEM,NRC,IDENT,IWRITE,NOR) C ****************************************************************** C * * C * P R P L O T * C * * C * PURPOSE: PRINTPLOT OF Y VERSUS X * C * IF IDENT.GT.0, THE FIRST NRC POINTS ARE IDENTIFIED BY* C * INTEGER NUMBERS 1-9,0,1-9... , THE REST OF THE POINTS* C * BY CHARACTERS A-I,J,A-I... * C * IF IDENT.EQ.0, ALL POINTS ARE PRINTED AS STARS * C * IF IDENT.LT.0, THE FIRST NRC POINTS ARE IDENTIFIED BY* C * STARS, THE OTHERS BY A 'D' * C * * C * SUBROUTINES CALLED: NONE * C * * C * AUTHOR: WILLEM HEISER RELEASED: MARCH 1978 * C * ADAPTED : NOVEMBER 1982 * C * ADAPTED : DECEMBER 1984 * C * * C ****************************************************************** C C TO ADJUST THE PLOT CHANGE THE CONSTANTS UNDER WHICH C ## IS PLACED ACCORDING TO THE COMMENTS C DIMENSION LINE(71),XPR(11),X(NITEM),Y(NITEM),LABEL(20),IND(NITEM) C ##=NN+1 (NN=INT(NC/7)) C ##=NUMBER OF COLUMNS PER LINE=NC INTEGER BLANK,STAR DATA BLANK,STAR,MORE/1H ,1H*,1HM/,NL/56/,NC/71/,NN/10/ C ##=NC ##=NN C ##=NUMBER OF LINES PER PAGE DATA LABEL/1H0,1H1,1H2,1H3,1H4,1H5,1H6,1H7,1H8,1H9, X 1HJ,1HA,1HB,1HC,1HD,1HE,1HF,1HG,1HH,1HI/ C 1 FORMAT(16X,3H.--,10(7H+------),5H+---.) C ##=NN 2 FORMAT(8X,F7.3,1X,1HI,2X,71A1,3X,1HI) C ##=NC C ##=NC+5 4 FORMAT(15X,11F7.3) C ##=NN+1 C DO 10 I=1,NITEM 10 IND(I)=I C C THE INDICES ARE ORDERED SUCH THAT THEY WOULD ORDER Y IN DESCENDING C ORDER C M=NITEM 20 M=M/2 IF (M.LT.1) GO TO 40 K=NITEM-M J=1 41 I=J 49 L=I+M IF (Y(IND(I)).GE.Y(IND(L))) GO TO 60 II=IND(I) IND(I)=IND(L) IND(L)=II I=I-M IF (I.GE.1) GO TO 49 60 J=J+1 IF (J-K) 41,41,20 40 CONTINUE C C THE PLOT IS ADJUSTED TO THE RANGE OF THE 'LONGEST' AXIS C XMAX=X(1) XMIN=X(1) DO 70 I=2,NITEM IF (X(I).GT.XMAX) XMAX=X(I) IF (X(I).LT.XMIN) XMIN=X(I) 70 CONTINUE IF (XMIN.GT. 0.0.AND.NOR.EQ.1) XMIN=0.0 IF (XMAX.LT. 0.0.AND.NOR.EQ.1) XMAX=0.0 SPANX=XMAX-XMIN YMAX=Y(IND(1)) IF (YMAX.LT. 0.0.AND.NOR.EQ.1) YMAX=0.0 YMIN=Y(IND(NITEM)) IF (YMIN.GT. 0.0.AND.NOR.EQ.1) YMIN=0.0 SPANY=YMAX-YMIN IF (SPANY.LE.SPANX) GO TO 75 XMIN=XMIN-(SPANY-SPANX)/2.0 SPANX=SPANY 75 STEPX=SPANX/(7.0*NN) STEPY=(STEPX*NC)/(NL+0.) HSTEP=STEPY/2.0 TOP=(YMIN+YMAX+SPANX)/2.0 C C THE POINTS ARE PLOTTED LINE AFTER LINE C WRITE(IWRITE,1) L=1 I=1 YPR=TOP+STEPY IF (NOR.EQ.1) IP=(0.0-XMIN)/STEPX+1.0 80 YPR=YPR-STEPY 90 DO 95 J=1,NC 95 LINE(J)=BLANK IF ((YPR-HSTEP.LE.0.0.AND.YPR+HSTEP.GT.0.0).AND.NOR.EQ.1) X LINE(IP)=STAR IF (I.GT.NITEM) GO TO 110 IF ((YPR-Y(IND(I))).GT.HSTEP) GO TO 110 100 JP=(X(IND(I))-XMIN)/STEPX+1.0 IF (LINE(JP).EQ.BLANK) GO TO 101 LINE(JP)=MORE GO TO 103 101 IF (IDENT.EQ.0) GO TO 102 IF (IDENT.GT.0) GO TO 104 IF (IND(I).LE.NRC) LINE(JP)=STAR IF (IND(I).GT.NRC) LINE(JP)=LABEL(15) GO TO 103 104 IF (IND(I).LE.NRC) LINE(JP)=LABEL(MOD(IND(I),10)+1) IF (IND(I).GT.NRC) X LINE(JP)=LABEL(MOD((IND(I)-NRC),10)+11) GO TO 103 102 LINE(JP)=STAR 103 I=I+1 IF (I.GT.NITEM) GO TO 105 IF ((YPR-Y(IND(I))).LT.HSTEP) GO TO 100 105 WRITE(IWRITE,2) YPR,LINE GO TO 115 110 WRITE(IWRITE,2) YPR,LINE 115 L=L+1 IF (L.LE.NL) GO TO 80 WRITE(IWRITE,1) C C VALUES OF X-AXIS ARE PRINTED C XPR(1)=XMIN DO 130 J=1,10 130 XPR(J+1)=XPR(J)+STEPX*7.0 WRITE(IWRITE,4) XPR RETURN END SUBROUTINE READIN(DATA,NOBJ,NVAR,LUNI,FORMAT) C ******************************************************************* C * * C * R E A D I N * C * * C * PURPOSE:READING OF DATA * C * * C ******************************************************************* DIMENSION DATA(NOBJ,NVAR),FORMAT(20) INTEGER DATA DO 10 I=1,NOBJ READ(LUNI,FORMAT) (DATA(I,J),J=1,NVAR) 10 CONTINUE RETURN END SUBROUTINE RESUL1 (DATA,NOBJ,NVAR,XROW,YCOL,KJ,WORK1,ICUMCA,NCAT, X DR,DC,WORK,MAXCAT,NORM,IPRI,IWRI,LUN,MIS) C ****************************************************************** C * * C * R E S U L 1 * C * * C * PURPOSE: PRINTING OF OBJECTS-SCORES, CATEGORY-QUANTIFICA- * C * TIONS, AND DISCRIMINATION-MEASURES. * C * * C ****************************************************************** DIMENSION DATA(NOBJ,NVAR),XROW(NOBJ),YCOL(KJ),ICUMCA(NVAR), X WORK1(NVAR),DR(NOBJ),DC(KJ),NCAT(NVAR),WORK(NOBJ),IPRI(10), X IWRI(10),LUN(5) COMMON/CP2/IWRITE,INPARA INTEGER DATA DATA S/1H*/ C ** C ** COMPUTE DISCRIMINATION-MEASURES FOR VARIABLES OR OBJECTS C ** IF (NORM.EQ.2) GO TO 500 DO 10 J=1,NVAR WORK1(J)=0.0 NMISJ=0 DM=0.0 DO 11 I=1,NOBJ IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) NMISJ=NMISJ+1 11 CONTINUE IF (NMISJ.NE.0.AND.MIS.EQ.0) GO TO 15 K=ICUMCA(J)+1 L=ICUMCA(J)+NCAT(J) DO 20 I=K,L WORK1(J)=WORK1(J)+YCOL(I)*YCOL(I)*DC(I)/NOBJ 20 CONTINUE IF (NMISJ.EQ.0.OR.MIS.EQ.1) GO TO 10 15 DO 30 I=1,NOBJ MMJ=0 IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) MMJ=1 IF (MMJ.NE.1) GO TO 35 A=XROW(I) IF (MMJ.EQ.1) GO TO 36 35 INDJ=DATA(I,J)+ICUMCA(J) A=YCOL(INDJ) 36 DM=DM+A*A 30 CONTINUE WORK1(J)=DM/FLOAT (NOBJ) 10 CONTINUE C C ** WRITE DISCRIMINATION-MEASURES ** C WRITE (IWRITE,614) WRITE (IWRITE,622) CALL PRIREA(WORK1,1,NVAR,1) GO TO 550 C 500 DO 110 I=1,NOBJ IF (MIS.EQ.1) WORK(I)=XROW(I)*XROW(I)*DR(I)/FLOAT (NVAR) IF (MIS.EQ.0) WORK(I)=XROW(I)*XROW(I) 110 CONTINUE C C ** WRITE DISCRIMINATION-MEASURES ** C WRITE (IWRITE,614) WRITE (IWRITE,623) CALL PRIREA(WORK,1,NOBJ,1) C C ** WRITE OBJECT SCORES ** C 550 IF (IPRI(2).EQ.1) WRITE (IWRITE,605) IF (IPRI(2).EQ.1) CALL PRIREA(XROW,1,NOBJ,1) IF (IWRI(2).EQ.1) CALL WRITE1(LUN(2),XROW,NOBJ) C C ** WRITE CATEGORY QUANTIFICATIONS ** C IF ((IPRI(3).NE.1.AND.IPRI(3).NE.11).AND.IWRI(3).NE.11) GO TO 560 IF (IPRI(3).EQ.1) WRITE (IWRITE,640) IF (IPRI(3).EQ.11) WRITE (IWRITE,740) IF (IPRI(3).EQ.1.OR.IPRI(3).EQ.11) WRITE (IWRITE,641) IF (MAXCAT.GT.14) MM=14 IF (MAXCAT.LE.14) MM=MAXCAT IF (IPRI(3).EQ.1.OR.IPRI(3).EQ.11) WRITE (IWRITE,642) (J,J=1,MM) KK=MM*8+2 IF (IPRI(3).EQ.1.OR.IPRI(3).EQ.11) WRITE (IWRITE,643) (S,J=1,KK) DO 5 J=1,NVAR K=ICUMCA(J)+1 L=NCAT(J)+ICUMCA(J) IF (IPRI(3).EQ.1) WRITE (IWRITE,644) J,(YCOL(I),I=K,L) IF (IPRI(3).NE.11.AND.IWRI(3).NE.11) GO TO 5 DO 4 I=K,L IF (NORM.EQ.1) YCOL(I)=YCOL(I)/SQRT (WORK1(J)) 4 CONTINUE IF (IPRI(3).EQ.11) WRITE (IWRITE,644) J,(YCOL(I),I=K,L) 5 CONTINUE 560 IF (IWRI(3).EQ.1.OR.IWRI(3).EQ.11) CALL WRITE1(LUN(2),YCOL,KJ) 605 FORMAT(1H0////,5X,13HOBJECT SCORES) 614 FORMAT(1H1//,5X,23HDISCRIMINATION MEASURES) 622 FORMAT (1H0,13X,9HVARIABLES) 623 FORMAT (1H0,13X,7HOBJECTS) 640 FORMAT (1H0////,5X,24HCATEGORY QUANTIFICATIONS) 740 FORMAT (1H0////,5X,37HSTANDARDIZED CATEGORY QUANTIFICATIONS) 641 FORMAT (1H0/,17X,10HCATEGORIES) 642 FORMAT (1H0,15X,15I8) 643 FORMAT (1H ,4X,9HVARIABLES,3X,120A1) 644 FORMAT (1H ,10X,I5,1X,1H*,14F8.3/(17X,1H*,14F8.3)) RETURN END SUBROUTINE RESUL2 (DATA,NOBJ,NVAR,XROW,YCOL,KJ,CM,WORK1,WORK2, X ICUMCA,NCAT,DC,A,B,D,MAXCAT,FIT,NORM,WORK3A, X WORK3B,IPRI,IPLO,IWRI,LUN,MISEST) C ****************************************************************** C * * C * R E S U L 2 * C * * C * PURPOSE: COMPUTE AND PRINT CORRELATION MATRIX, PRINCIPAL * C * COMPONENTS SOLUTION BEFORE TRANSFORMATION * C * * C ****************************************************************** DIMENSION DATA(NOBJ,NVAR),XROW(NOBJ),YCOL(KJ),ICUMCA(NVAR), X CM(NVAR,NVAR),WORK1(NVAR),WORK2(NVAR),DC(KJ),A(NOBJ), X B(NOBJ),D(KJ),WORK3A(NVAR),WORK3B(NVAR),NCAT(NVAR), X IPRI(10),IPLO(11),IWRI(10),LUN(5) COMMON/CP2/IWRITE,INPARA C C ** A AND B ARE WORK-ARRAYS, THE ARRAYS PREVIOUSLY FILLED WITH DR AND C ** DRINV ARE USED FOR THEM C INTEGER DATA C C ** ARRAY WORK1 IS FILLED WITH THE ROWS OF THE DATA-MATRIX, WHILE C ** THE LABELS ARE REPLACED BY THE CATEGORY-QUANTIFICATIONS. C IF (IPRI(4).NE.1.AND.IWRI(4).NE.1) GO TO 510 C C ** PRINT TRANSFORMED DATA MATRIX OR WRITE TO LOGICAL UNIT NUMBER C IF (MISEST.NE.2) GO TO 20 C C ** THE MEAN VALUES OF THE NON-MISSING RELEVANT C ** OBJECT SCORES ARE SUBSTITUTED FOR MISSING DATA C DO 100 J=1,NVAR SUM=0.0 KTEL=0 DO 110 I=1,NOBJ IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) SUM=SUM+XROW(I) IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) KTEL=KTEL+1 110 CONTINUE IF (KTEL.EQ.0) GO TO 100 AVE=SUM/ FLOAT (KTEL) WORK2(J)=AVE 100 CONTINUE C C ** WRITE TRANSFORMED DATAMATRIX TO LOGICAL UNIT NUMBER C 20 IF (IPRI(4).EQ.1) WRITE (IWRITE,620) DO 25 I=1,NOBJ DO 26 J=1,NVAR MMM=0 IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) MMM=1 IF (MMM.NE.1) GO TO 24 IF (MISEST.EQ.0) WORK1(J)=XROW(I) IF (MISEST.EQ.1) WORK1(J)=0.0 IF (MISEST.EQ.2) WORK1(J)=WORK2(J) IF (MISEST.EQ.0.AND.NORM.EQ.2) WORK1(J)=WORK1(J)/FIT IF (MMM.EQ.1) GO TO 26 24 IND=DATA(I,J)+ICUMCA(J) WORK1(J)=YCOL(IND) 26 CONTINUE IF (IPRI(4).EQ.1) WRITE (IWRITE,610) (WORK1(J),J=1,NVAR) IF (IWRI(4).EQ.1) LUNO=LUN(3) IF (IWRI(4).EQ.1) WRITE (LUNO,615) (WORK1(J),J=1,NVAR) 25 CONTINUE C ** C ** COMPUTE CORRELATIONS RAW DATA C ** 510 ICCB=IPRI(5)+IPRI(7)+IPRI(8)+IPLO(7)+IPLO(8)+IWRI(5)+IWRI(7) X +IWRI(8) IF (ICCB.EQ.0) GO TO 700 JE=NVAR-1 CM(NVAR,NVAR)=1.0 DO 27 J=1,JE KB=J+1 CM(J,J)=1.0 DO 28 K=KB,NVAR NOJ=0 NOK=0 DO 29 I=1,NOBJ A(I)=DATA(I,J) IF (A(I).GT.NCAT(J)) A(I)=0.0 IF (A(I).NE.0.0) NOJ=NOJ+1 B(I)=DATA(I,K) IF (B(I).GT.NCAT(K)) B(I)=0.0 IF (B(I).NE.0.0) NOK=NOK+1 29 CONTINUE CALL CORMA(A,B,NOBJ,NVAR,RHO,J,K,WORK3A,WORK3B,NOJ,NOK,1) CM(J,K)=RHO CM(K,J)=RHO 28 CONTINUE 27 CONTINUE C ** C ** WRITE CORRELATIONS RAW DATA C ** IF (IPRI(5).EQ.1) WRITE (IWRITE,625) IF (IPRI(5).EQ.1) CALL PRIREA(CM,NVAR,NVAR,1) IF (IWRI(5).NE.1) GO TO 11 DO 10 J=1,NVAR LUNO=LUN(4) WRITE (LUNO,615) (CM(J,K),K=1,NVAR) 10 CONTINUE C ** C ** COMPUTE EIGENVALUES AND EIGENVECTORS OF CORRELATION-MATRIX C ** 11 CALL TRED2(NVAR,NVAR,WORK1,WORK2,CM) CALL IMTQL2(NVAR,NVAR,WORK1,WORK2,CM,IERR) IF (IERR.GT.0) GO TO 1000 IF (IPRI(5).EQ.1) WRITE (IWRITE,630) IF (IPRI(5).NE.1) WRITE (IWRITE,631) WRITE (IWRITE,610) (WORK1(J),J=1,NVAR) C IPCAB=(ICCB-IPRI(5))-IWRI(5) C C ** ICCB INDICATED THE NECESSITY TO COMPUTE THE CORRELATION MATRIX C ** BEFORE TRANSFORMATION. IPCAB.NE.0 INDICATES THE COMPUTATION OF C ** THE PRINCIPAL COMPONENTS SOLUTION C IF (IPCAB.NE.0) CALL PCA(DATA,NOBJ,NVAR,A,B,DC,D,WORK1,WORK2,CM, X ICUMCA,NCAT,XROW,YCOL,KJ,WORK3A,WORK3B,KRAW,FIT,NORM, X IPRI,IPLO,IWRI,LUN,IPCAB,MISEST) GO TO 700 C 1000 WRITE (IWRITE,1001) 1001 FORMAT (1H0//,5X,14HNO CONVERGENCE) 610 FORMAT(1H0,8X,15F8.3) 615 FORMAT(10F8.4) 620 FORMAT(1H1//,5X,23HTRANSFORMED DATA MATRIX) 625 FORMAT(1H1//,5X,40HCORRELATION MATRIX BEFORE TRANSFORMATION) 630 FORMAT (1H0//,5X,11HEIGENVALUES) 631 FORMAT (1H0//,5X,37HEIGENVALUES CORRELATION MATRIX BEFORE, X 15H TRANSFORMATION) 700 RETURN END SUBROUTINE RESUL3 (DATA,NOBJ,NVAR,XROW,YCOL,KJ,CM,WORK1,WORK2, X ICUMCA,NCAT,DC,A,B,D,MAXCAT,FIT,NORM, X WORK3A,WORK3B,IPRI,IPLO,IWRI,LUN,MISEST) C ****************************************************************** C * * C * R E S U L 3 * C * * C * PURPOSE: COMPUTE AND PRINT CORRELATION MATRIX AND PRINCI- * C * PAL COMPONENTS SOLUTION AFTER TRANSFORMATION * C * * C ****************************************************************** DIMENSION DATA(NOBJ,NVAR),XROW(NOBJ),YCOL(KJ),ICUMCA(NVAR), X CM(NVAR,NVAR),WORK1(NVAR),WORK2(NVAR),DC(KJ),A(NOBJ), X B(NOBJ),D(KJ),WORK3A(NVAR),WORK3B(NVAR),NCAT(NVAR), X IPRI(10),IPLO(11),IWRI(10),LUN(5) COMMON/CP2/IWRITE,INPARA C C ** A AND B ARE WORK-ARRAYS, THE ARRAYS PREVIOUSLY FILLED WITH DR AND C ** DRINV ARE USED FOR THEM C INTEGER DATA ICCA=IPRI(6)+IPRI(9)+IPRI(10)+IPLO(9)+IPLO(10)+IWRI(6)+IWRI(9) X +IWRI(10) IF (ICCA.EQ.0) GO TO 700 IF (MISEST.NE.2) GO TO 500 DO 100 J=1,NVAR SUM=0.0 KTEL=0 DO 110 I=1,NOBJ IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) SUM=SUM+XROW(I) IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) KTEL=KTEL+1 110 CONTINUE IF (KTEL.EQ.0) GO TO 100 AVE=SUM/FLOAT (KTEL) WORK2(J)=AVE 100 CONTINUE C C ** C ** COMPUTE CORRELATIONS AFTER TRANSFORMATIONS C ** 500 JE=NVAR-1 CM(NVAR,NVAR)=1.0 DO 30 J=1,JE KB=J+1 CM(J,J)=1.0 DO 40 K=KB,NVAR DO 50 I=1,NOBJ MMJ=0 IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) MMJ=1 IF (MMJ.NE.1) GO TO 54 IF (MISEST.EQ.0) A(I)=XROW(I) IF (MISEST.EQ.1) A(I)=0.0 IF (MISEST.EQ.2) A(I)=WORK2(J) IF (MISEST.EQ.0.AND.NORM.EQ.2) A(I)=A(I)/FIT IF (MMJ.EQ.1) GO TO 55 54 INDJ=DATA(I,J)+ICUMCA(J) A(I)=YCOL(INDJ) 55 MMK=0 IF (DATA(I,K).EQ.0.OR.DATA(I,K).GT.NCAT(K)) MMK=1 IF (MMK.NE.1) GO TO 49 IF (MISEST.EQ.0) B(I)=XROW(I) IF (MISEST.EQ.1) B(I)=0.0 IF (MISEST.EQ.2) B(I)=WORK2(K) IF (MISEST.EQ.0.AND.NORM.EQ.2) B(I)=B(I)/FIT IF (MMK.EQ.1) GO TO 50 49 INDK=DATA(I,K)+ICUMCA(K) B(I)=YCOL(INDK) 50 CONTINUE CALL CORMA(A,B,NOBJ,NVAR,RHO,J,K,WORK3A,WORK3B,NOBJ,NOBJ,0) CM(J,K)=RHO CM(K,J)=RHO 40 CONTINUE 30 CONTINUE IF (IPRI(6).EQ.1) WRITE (IWRITE,625) IF (IPRI(6).EQ.1) CALL PRIREA(CM,NVAR,NVAR,1) IF (IWRI(6).NE.1) GO TO 121 LUNO=LUN(4) DO 120 J=1,NVAR WRITE (LUNO,615) (CM(J,K),K=1,NVAR) 120 CONTINUE C 121 CALL TRED2(NVAR,NVAR,WORK1,WORK2,CM) CALL IMTQL2(NVAR,NVAR,WORK1,WORK2,CM,IERR) IF (IERR.GT.0) GO TO 1000 IF (IPRI(6).EQ.1) WRITE (IWRITE,630) IF (IPRI(6).NE.1) WRITE (IWRITE,631) WRITE (IWRITE,610) (WORK1(J),J=1,NVAR) C IPCAB=0 IPCAA=(ICCA-IPRI(6))-IWRI(6) IF (IPCAA.NE.0) CALL PCA(DATA,NOBJ,NVAR,A,B,DC,D,WORK1,WORK2,CM, X ICUMCA,NCAT,XROW,YCOL,KJ,WORK3A,WORK3B,KRAW,FIT, X NORM,IPRI,IPLO,IWRI,LUN,IPCAB,MISEST) GO TO 700 C C 1000 WRITE (IWRITE,1001) 1001 FORMAT (1H0//,5X,14HNO CONVERGENCE) 610 FORMAT(1H0,8X,15F8.3) 615 FORMAT(10F8.4) 625 FORMAT(1H1//,5X,39HCORRELATION MATRIX AFTER TRANSFORMATION) 630 FORMAT (1H0//,5X,11HEIGENVALUES) 631 FORMAT (1H0//,5X,45HEIGENVALUES CORRELATION MATRIX AFTER TRANSFOR, X 6HMATION) 700 RETURN END SUBROUTINE TRANXA(XNEW,YOLD,DATA,NOBJ,NOB,NOBJT,NVAR,KJ,ICUMCA, X NCAT,DRINV,MIS,IEIG,NORM) C ******************************************************************* C * * C * T R A N X A * C * * C * PURPOSE: NEW VALUES FOR OBJECT-SCORES * C * OPTIONS: LARGEST OR SMALLEST EIGENVALUE (IEIG) * C * : NORMALIZE X OR Y (NORM) * C * : MISSING DATA PASSIVE = MISSING DATA (MISS=1) * C * : MISSING DATA MULTIPLE= NO MISSING DATA (MISS=0) * C * * C * DURING THE ITERATION PROCES NOB HAS THE SAME VALUE AS NOBJ, SO * C * * C * THE VALUE OF MISS INDICATES WHETHER MISSING DATA SHOULD BE * C * * C * DELETED OR SHOULD BE TREATED AS MULTIPLE CATEGORIES. WHEN THERE * C * * C * HAPPEN TO BE PASSIVE OBJECTS WITH MISSING DATA, OPTION 'MISSING * C * * C * DATA DELETED' IS REQUIRED TO ENSURE THAT THE PASSIVE OBJECTS * C * * C * ARE IN THE CENTROID OF THE NON MISSING CATEGORIES. THAT IS WHY * C * * C * DURING THE LAST COMPUTATION FOR OBJECT SCORES, WHEN NOB INDI- * C * * C * CATES THE NUMBER OF ACTIVE OBJECTS AND NOBJ INDICATES THE TOTAL * C * * C * NUMBER OF OBJECTS, THE VALUE OF MISS IS SET TO 1 (PASSIVE) FOR * C * * C * THE PASSIVE OBJECTS. * C ******************************************************************* DIMENSION DATA(NOBJT,NVAR),YOLD(KJ),XNEW(NOBJ),ICUMCA(NVAR), X NCAT(NVAR),DRINV(NOBJ) INTEGER DATA DO 10 I=1,NOBJ IF (I.LE.NOB) MISS=MIS IF (I.GT.NOB) MISS=1 IF (DRINV(I).EQ.0.0.AND.MISS.EQ.1) OLD=0.0 IF (DRINV(I).NE.0.0.OR.MISS.EQ.0) OLD=XNEW(I) 15 XNEW(I)=0.0 DO 20 J=1,NVAR IF ((DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)).AND.MISS.EQ.0) X XNEW(I)=XNEW(I)+OLD IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) GO TO 20 IND=DATA(I,J)+ICUMCA(J) XNEW(I)=XNEW(I)+YOLD(IND) 20 CONTINUE IF ((NORM.EQ.2.OR.IEIG.EQ.0).AND.MISS.EQ.0) X XNEW(I)=XNEW(I)/NVAR IF ((NORM.EQ.2.OR.IEIG.EQ.0).AND.MISS.EQ.1) X XNEW(I)=XNEW(I)*DRINV(I) IF ((NORM.EQ.1.AND.IEIG.EQ.1).AND.MISS.EQ.0) X XNEW(I)=OLD-(XNEW(I)/NVAR) IF ((NORM.EQ.1.AND.IEIG.EQ.1).AND.MISS.EQ.1) X XNEW(I)=OLD-(XNEW(I)*DRINV(I)) 10 CONTINUE RETURN END SUBROUTINE TRANYA(YNEW,XOLD,DATA,NOBJ,NOBJT,NVAR,KJ,DCINV,ICUMCA, X NCAT,IEIG,NORM,OLD,MAXCAT) C ******************************************************************* C * * C * T R A N Y A * C * * C * PURPOSE: NEW VALUES FOR CATEGORY-QUANTIFICATIONS * C * * C ******************************************************************* DIMENSION DATA(NOBJT,NVAR),DCINV(KJ),YNEW(KJ),XOLD(NOBJ), X ICUMCA(NVAR),NCAT(NVAR),OLD(MAXCAT) INTEGER DATA DO 10 J=1,NVAR KB=ICUMCA(J)+1 KE=ICUMCA(J)+NCAT(J) L=0 DO 20 K=KB,KE IF (NORM.EQ.1.OR.IEIG.EQ.0) GO TO 15 L=L+1 IF (DCINV(K).EQ.0.0) OLD(L)=0.0 IF (DCINV(K).NE.0.0) OLD(L)=YNEW(K) 15 YNEW(K)=0.0 20 CONTINUE DO 30 I=1,NOBJ IF (DATA(I,J).EQ.0.OR.DATA(I,J).GT.NCAT(J)) GO TO 30 IND=DATA(I,J)+ICUMCA(J) YNEW(IND)=YNEW(IND)+XOLD(I) 30 CONTINUE L=0 DO 40 K=KB,KE L=L+1 IF(NORM.EQ.1.OR.IEIG.EQ.0) YNEW(K)=YNEW(K)*DCINV(K) C IF(NORM.EQ.2.AND.IEIG.EQ.1) CALL CENORC(YNEW,KJ,NVAR,DC,SRC, C X FIT,IEIG,ICUMCA,NCAT) IF(NORM.EQ.2.AND.IEIG.EQ.1) X YNEW(K)=OLD(L)-(YNEW(K)*DCINV(K)) 40 CONTINUE 10 CONTINUE RETURN END SUBROUTINE TRED2(NM,N,D,E,Z) C ****************************************************************** C * * C * T R E D 2 * C * * C * PURPOSE: REDUCES A REAL SYMMETRIC MATRIX TO A SYMMETRIC TRI- * C * DIAGONAL MATRIX USING ORTHOGONAL SIMILARITY TRANS- * C * FORMATIONS; THE ACCUMULATED ORTHOGONAL TRANSFORMA- * C * TIONS ARE RETAINED IN Z. * C * * C * SUBROUTINES CALLED: NONE * C * * C * ADAPTED FROM EISPACK * C * * C ****************************************************************** DIMENSION Z(NM,N),D(N),E(N) C IF (N .EQ. 1) GO TO 320 C ********** FOR I=N STEP -1 UNTIL 2 DO -- ********** DO 300 II =2, N I = N + 2 - II L = I - 1 H = 0.0 SCALE = 0.0 IF (L .LT.2) GO TO 130 C ********** SCALE ROW (ALGOL TOL THEN NOT NEEDED) ********** DO 120 K =1, L 120 SCALE = SCALE + ABS(Z(I,K)) C IF (SCALE .NE. 0.0) GO TO 140 130 E(I) = Z(I,L) GO TO 290 C 140 DO 150 K = 1, L Z(I,K) = Z(I,K) / SCALE H = H + Z(I,K) * Z(I,K) 150 CONTINUE C F = Z(I,L) G = -SIGN(SQRT(H),F) E(I) = SCALE * G H = H - F * G Z(I,L) = F - G F = 0.0 C DO 240 J = 1, L Z(J,I) = Z(I,J) / (SCALE * H) G = 0.0 C ********** FORM ELEMENT OF A*U ********** DO 180 K = 1, J 180 G = G + Z(J,K) * Z(I,K) C JP1 = J + 1 IF (L .LT. JP1) GO TO 220 C DO 200 K = JP1, L 200 G = G + Z(K,J) * Z(I,K) C ********** FORM ELEMENT OF P ********** 220 E(J) = G / H F = F + E(J) * Z(I,J) 240 CONTINUE C HH = F / (H + H) C ********** FORM REDUCED A ********** DO 260 J =1, L F = Z(I,J) G = E(J) - HH * F E(J) = G C DO 260 K = 1, J Z(J,K) = Z(J,K) - F * E(K) -G * Z(I,K) 260 CONTINUE C DO 280 K = 1, L 280 Z(I,K) = SCALE * Z(I,K) C 290 D(I) = H 300 CONTINUE C 320 D(1) = 0.0 E(1) = 0.0 C ********** ACCUMULATION OF TRANSFORMATION MATRICES ********** DO 500 I = 1, N L = I - 1 IF (D(I) .EQ. 0.0) GO TO 380 C C DO 360 J = 1, L G = 0.0 DO 340 K =1, L 340 G = G + Z(I,K) * Z(K,J) C DO 360 K = 1, L Z(K,J) = Z(K,J) - G * Z(K,I) 360 CONTINUE C 380 D(I) = Z(I,I) Z(I,I) = 1.0 IF (L .LT. 1) GO TO 500 C DO 400 J =1, L Z(I,J) = 0.0 Z(J,I) = 0.0 400 CONTINUE C 500 CONTINUE C RETURN END SUBROUTINE WRITE1(LUNO,A,N) C ******************************************************************* C * * C * W R I T E 2 * C * * C * PURPOSE: WRITE OUTPUT TO LOGICAL UNIT NUMBER * C * * C ******************************************************************* DIMENSION A(N) DO 10 I=1,N WRITE(LUNO,600) I,A(I) 10 CONTINUE 600 FORMAT (I10,1F12.7) RETURN END SUBROUTINE WRITE2(LUNO,A,B,N) C ******************************************************************* C * * C * W R I T E 2 * C * * C * PURPOSE: WRITE OUTPUT TO LOGICAL UNIT NUMBER * C * * C ******************************************************************* DIMENSION A(N),B(N) DO 10 I=1,N WRITE(LUNO,600) I,A(I),B(I) 10 CONTINUE 600 FORMAT (I10,2F12.7) RETURN END