Trailing-Edge
-
PDP-10 Archives
-
decuslib10-09
-
43,50466/stp4.stp
There is 1 other file named stp4.stp in the archive. Click here to see a list.
C *** STAT PACK ***
C SUBROUTINE FOR KENDALL TAU CORRELATION.
C CALLING SEQUENCE: CALL STKTAU(NV,NC,MV,MC,DATA,AV,NAMES)
C WHERE NV - IS THE NUMBER OF COLUMNS ACTUALLY FILLED(VARIABLES)
C NC - IS THE NUMBER OF ROWS ACTUALLY FILLED(OBSERVATIONS,CASES)
C MV - IS MAXIMUM NUMBER OF COLUMNS, AS SPECIFIED IN MAIN.
C MC - IS MAXIMUM NUMBER OF ROWS, AS SPECIFIED IN MAIN.
C DATA - STORAGE FOR DATA DIMENSIONED FOR MAXIMUM MATRIX.
C AV - IS A VECTOR DIMENSIONED AT LEAST NV.
C NAMES - IS A VECTOR CONTAINING VARIABLE NAMES
C
C PROGRAM OUTPUTS ENTIRE KENDALL TAU CORRELATION MATRIX FOR ALL DATA
C
SUBROUTINE STKTAU(NV,NC,MV,MC,DATA,AV,NAMES)
DIMENSION DATA(MC,MV),AV(1),NAMES(1)
COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
COMMON /PRNT/ LINPP,ICOPS,RUNPRG
COMMON/EXTRA/HEDR(70),NSZ
IF(IOUT.NE.21) WRITE(IOUT,5566)(HEDR(I),I=1,NSZ)
5566 FORMAT('1',70A1)
IF(IOUT.EQ.21) CALL PRNTHD
ISQ=7
IF(IOUT.EQ.21) ISQ=15
SWTCH=0
WRITE(IOUT,10)
10 FORMAT('0',10X,'***** KENDALL TAU CORRELATION MATRIX *****')
LINES=4
14 DO 1 I=1,NV
DO 2 J=1,I
IF(I.NE.J) GO TO 3
AV(J)=1.0
GO TO 2
3 AM=0.
TX=0.
TY=0.
DO 4 K=1,NC-1
SX=DATA(K,I)
SY=DATA(K,J)
DO 4 L=K+1,NC
SM=(SX-DATA(L,I))*(SY-DATA(L,J))
IF(SM.EQ.0) GO TO 5
IF(SM.LT.0) AM=AM-1.
IF(SM.GT.0) AM=AM+1.
GO TO 4
5 IF((SX-DATA(L,I)).EQ.0)TX=TX+1.
IF((SY-DATA(L,J)).EQ.0)TY=TY+1.
4 CONTINUE
TC=NC
XNGO=.5*(TC*(TC-1.))
DEN=(XNGO-TX)*(XNGO-TY)
IF(DEN.GT.0) GOTO 6
SWTCH=1.
AV(J)=100**4
GO TO 2
6 AV(J)=AM/SQRT(DEN)
2 CONTINUE
IF(IOUT.NE.21) GO TO 17
KK=(I+ISQ-1)/ISQ+1
LINES=LINES+KK
IF(LINES.LE.(LINPP-KK)) GO TO 17
WRITE(IOUT,15)
DO 16 J=1,I-1,ISQ
NEND=J-ISQ-1
IF(NEND.GT.(I-1)) NEND=I-1
16 WRITE(IOUT,8)(NAMES(K),K=J,NEND)
CALL PRNTHD
LINES=2+KK
17 DO 1 J=1,I,ISQ
NEND=J+ISQ-1
IF(NEND.GT.I) NEND=I
IF(J.EQ.1) WRITE(IOUT,7) NAMES(I),(AV(K),K=J,NEND)
IF(J.NE.1) WRITE(IOUT,12)(AV(K),K=J,NEND)
7 FORMAT('0',A5,2X,15F8.4)
12 FORMAT(8X,15F8.4)
1 CONTINUE
WRITE(IOUT,15)
15 FORMAT(1X)
DO 13 J=1,NV,ISQ
NEND=J+ISQ-1
IF(NEND.GT.NV) NEND=NV
13 WRITE(IOUT,8)(NAMES(K),K=J,NEND)
8 FORMAT(8X,15(3X,A5))
IF(SWTCH.NE.0) WRITE(IOUT,9)
9 FORMAT('0NOTE: ALL STARS COULD NOT BE CALCULATED BECAUSE '
1,'ONE OF THE S.D.S WAS EQUAL TO ZERO')
RETURN
END
C *** STAT PACK ***
C SUBROUTINE FOR STEPWISE REGRESSION
C CALLING SEQUENCE: CALL STSTRG(NV,NC,MV,MC,DATA,COR,VMN,STD,CO,NAMES)
C WHERE NV - IS THE NUMBER OF COLUMNS ACTUALLY FILLED(VARIABLES)
C NC - IS THE NUMBER OF ROWS ACTUALLY FILLED(OBSERVATIONS,CASES)
C MV - IS MAXIMUM NUMBER OF COLUMNS, AS SPECIFIED IN MAIN.
C MC - IS MAXIMUM NUMBER OF ROWS, AS SPECIFIED IN MAIN.
C DATA - STORAGE FOR DATA DIMENSIONED FOR MAXIMUM MATRIX.
C COR - IS THE CORRELATION MATRIX.
C VMN - IS A VECTOR CONTAINING VARIABLE MEANS.
C STD - IS A VECTOR CONTAINING VARIABLE STANDARD DEVIATIONS.
C CO - IS A VECTOR DIMENSIONED AT LEAST MV
C NAMES - IS A VECTOR CONTAINING VARIABLE NAMES
C
C USER HAS MANY OPTIONS AVAILABLE IN STEPWISE REGRESSION,
C AMONG WHICH ARE OUTPUT FOR ANALYSIS OF VARIANCE AND THE DURBIN-WATSON
C TEST FOR AUTO CORRELATION. LIKE THE MULTIPLE REGRESSION, OPTIONS
C ARE AVAILABLE FOR STORING THE RESIDUALS.
C PROGRAM WAS ORIGINALLY OBTAINED FROM WAYNE STATE UNIVERSITY,
C MODIFICATIONS MADE FOR W.M.U., AND FINAL MODIFICATION FOR STP.
C ORIGINAL SOURCE DOCUMENT FOR STEPWISE REGRESSION WAS:
C RALSTON AND WILF "MATHEMATICAL METHODS FOR DIGITAL COMPUTERS".
C
SUBROUTINE STSTRG(NV,NC,MV,MC,DATA,COR,VMN,STD,CO,NAMES)
DIMENSION DATA(MC,MV),COR(MV,MV),VMN(1),IV(20),DCOR(20,20)
DIMENSION SP(20),SIGCO(20),STUDT(20),BETA(20),STD(1)
DIMENSION OPT(6),AVG(20),INDEX(20),SIGMA(20),COEN(20)
DIMENSION NAMES(1),IVA(20),CO(1)
COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
COMMON /PRNT/ LINPP,ICOPS,RUNPRG
COMMON/EXTRA/HEDR(70),NSZ
C
C ************************************************************
C
C OPTIONS ARE ENTERED AND INTERUPTED
C
21 TOL=.0001
COE=0
EFIN=0
EFOUT=0
NOFRC=0
OBSER=NC
IF(ICC.NE.2) WRITE(IDLG,16)
16 FORMAT('0LIST THE OPTIONS YOU WISH SEPARATED BY COMMAS'/)
READ (ICC,17)(SP(I),I=1,14)
17 FORMAT(14(A5,1X))
IF(SP(1).EQ.'!') RETURN
DO 22 I=1,6
22 OPT(I)=0
DO 18 I=1,14
EXIT=0
40 IF(SP(I).EQ.' ') GO TO 41
IF(SP(I).NE.'HELP') GO TO 19
WRITE(IDLG,20)
20 FORMAT('0POSSIBLE OPTIONS ARE:'/
1' "TOLER" - INDICATE A TOLERANCE OTHER THAN .0001'/
2' "F-VAL" - INDICATE F-VALUES FOR ENTERING OR OMITTING',
3' A VARIABLE'/
4' "FORCE" - INDICATE VARIABLES TO BE FORCED INTO THE REGRESSION'
5/' "DURWT" - OUTPUT DURBIN-WATSON TEST FOR AUTOCORRELATION'/
6' "ANOVA" - ANALYSIS OF VARIANCE'/
7' "RESID" - STORE RESIDUALS'/
8'0IF NO OPTIONS ARE DESIRED TYPE A RETURN')
GO TO 21
19 IF(SP(I).NE.'TOLER') GO TO 25
IF(OPT(1).NE.0) WRITE(IDLG,23)
23 FORMAT('0OPTION "TOLER" LISTED MORE THAN ONCE')
OPT(1)=1.
GO TO 18
25 IF(SP(I).NE.'F-VAL') GO TO 27
IF(OPT(2).NE.0) WRITE(IDLG,26)
26 FORMAT('0OPTION "F-VAL" LISTED MORE THAN ONCE')
OPT(2)=1
GO TO 18
27 IF(SP(I).NE.'FORCE') GO TO 30
IF(OPT(3).NE.0) WRITE(IDLG,28)
28 FORMAT('0OPTION "FORCE" LISTED MORE THAN ONCE')
OPT(3)=1
GO TO 18
30 IF(SP(I).NE.'DURWT') GO TO 32
IF(OPT(4).NE.0) WRITE(IDLG,31)
31 FORMAT('0OPTION "DURWT" LISTED MORE THAN ONCE')
OPT(4)=1
GO TO 18
32 IF(SP(I).NE.'ANOVA') GO TO 34
IF(OPT(5).NE.0) WRITE(IDLG,33)
33 FORMAT('0OPTION "ANOVA" LISTED MORE THAN ONCE')
OPT(5)=1
GO TO 18
34 IF(SP(I).NE.'RESID') GO TO 36
IF(OPT(6).NE.0) WRITE(IDLG,35)
35 FORMAT('0OPTION "RESID" LISTED MORE THAN ONCE')
OPT(6)=1
GO TO 18
36 WRITE(IDLG,37) SP(I)
37 FORMAT('0OPTION "',A5,'" DOES NOT EXIST')
GO TO 21
18 CONTINUE
C
41 IF(OPT(1).NE.1) GO TO 44
IF(ICC.NE.2) WRITE(IDLG,42)
42 FORMAT(' ENTER TOLERANCE ',$)
READ(ICC,43)TOL
C
C ************************************************************
C
43 FORMAT(F)
44 IF(OPT(2).NE.1) GO TO 60
IF(ICC.NE.2) WRITE(IDLG,45)
45 FORMAT(' ENTER F-VALUE FOR ENTERING A VARIABLE? ',$)
READ (ICC,43)EFIN
IF(ICC.NE.2) WRITE(IDLG,46)
46 FORMAT(' ENTER F-VALUE FOR OMITTING A VARIABLE? ',$)
READ(ICC,43) EFOUT
C
60 IF(OPT(6).NE.1) GOTO 3
57 IF(NV.EQ.MV) GO TO 56
62 IF(ICC.NE.2) WRITE(IDLG,61)
61 FORMAT('0WHICH VARIABLE DO YOU WANT THE RESIDUALS STORED UNDER? ',
1$)
GO TO 58
56 WRITE(IDLG,63)
63 FORMAT('0RESIDUALS CANNOT BE FIT IN WITHOUT OVERLAYING A VARIABLE'
1/' PLEASE INDICATE WHICH ONE',$)
58 CALL ALPHA(INDXKK,1,J,IRET,IHELP,IERR,NAMES,MV)
IF(IERR.EQ.1) GO TO 57
IF(IRET.EQ.1) RETURN
IF(IHELP.NE.1) GO TO 320
WRITE(IDLG,59)
59 FORMAT(' ENTER VARIABLE NAME OR NUMBER WHERE RESIDUALS CAN BE'/
1' STORED. NAME ASSIGNED TO RESIDUALS WILL BE "RESID"')
GO TO 57
320 IF(INDXKK.GT.0) GO TO 3
WRITE(IDLG,38)
GO TO 57
3 IF(ICC.NE.2) WRITE(IDLG,1)
1 FORMAT('0LIST THE INDEPENDENT VARIABLES?'/)
CALL ALPHA(IVA,19,NIND,IRET,IHELP,IERR,NAMES,NV)
IF(IRET.EQ.1) RETURN
IF(IERR.EQ.1) GO TO 3
IF(IHELP.EQ.1) GO TO 3
4 K=1
DO 5 I=1,NIND
IV(I)=IVA(I)
IF(IVA(I).GT.0) GO TO 5
IV(I)=K
K=K+1
5 CONTINUE
DO 6 I=1,NIND-1
IF(IVA(I).LT.0) GO TO 6
DO 7 J=I+1,NIND
IF(IVA(J).LT.0) GO TO 7
IF(IVA(I).NE.IVA(J)) GO TO 7
WRITE(IDLG,8)
8 FORMAT(' THE SAME INDEP. VARIABLE HAS BEEN LISTED TWICE')
GO TO 3
7 CONTINUE
6 CONTINUE
IF((NIND+1).LE.NV) GO TO 11
WRITE(IDLG,9)
9 FORMAT(' YOU HAVE LISTED MORE VARIABLES THAN POSSIBLE'/
1' WITH THE DATA SET AVAILABLE')
GO TO 3
11 NIND=NIND+1
13 IF(ICC.NE.2) WRITE(IDLG,10)
10 FORMAT(' WHICH IS THE DEPENDENT VARIABLE? ',$)
CALL ALPHA(J,1,I,IRET,IHELP,IERR,NAMES,NV)
IF(IRET.EQ.1) RETURN
IF(IERR.EQ.1) GOTO 13
IF(IHELP.EQ.1) GOTO 13
IVA(NIND)=J
IV(NIND)=IVA(NIND)
IF(IVA(NIND).LT.0) IV(NIND)=1
IF(IVA(NIND).LT.0) GO TO 50
DO 12 J=1,NIND-1
IF(IVA(J).LT.0) GO TO 12
IF(IVA(J).NE.IVA(NIND)) GO TO 12
WRITE(IDLG,2)
2 FORMAT(' THE DEPENDENT VARIABLE ALSO EXISTS AS AN INDEP. VAR.')
GO TO 3
12 CONTINUE
50 IF(OPT(3).NE.1) GO TO 318
55 IF(ICC.NE.2) WRITE(IDLG,51)
51 FORMAT(' ENTER THE VARIABLES TO BE FORCED INTO',
1' THE REGRESSION'/)
CALL ALPHA(INDEX,20,NOFRC,IRET,IHELP,IERR,NAMES,NV)
IF(IRET.EQ.1) RETURN
IF(IERR.EQ.1) GO TO 55
IF(IHELP.NE.1) GO TO 47
WRITE(IDLG,48)
48 FORMAT(' ENTER THE VARIABLE NAME OR NUMBER TO BE FORCED INTO'/
1' THE REGRESSION. ONLY THOSE VARIABLES WHCIH WERE LISTED AS'/
2' INDEPENDENT VARIABLES MAY BE USED. *, ?, ALL ARE NOT'/
3' ACCEPTABLE')
GO TO 55
47 DO 49 K=1,NOFRC
IF(INDEX(I).GT.0) GO TO 39
WRITE(IDLG,38)
38 FORMAT(' ?, *, ALL MAY NOT BE USED HERE')
GO TO 55
39 DO 53 J=1,NIND-1
IF(INDEX(K).EQ.IVA(J)) GO TO 52
53 CONTINUE
WRITE(IDLG,54) NAMES(INDEX(K))
54 FORMAT(' VARIABLE "',A5,'" WAS NOT LISTED AS A DEP. VAR.')
GO TO 55
52 INDEX(K)=J
49 CONTINUE
GO TO 318
C
C *************************************************************
C
C HERE THE DATA IS SET UP REPLACEING ALL *,?, AND ALL WITH VALID
C VARIABLES
C
300 J=NIND-1
314 IF(IVA(J).GT.0) GO TO 315
IV(J)=IV(J)+1
IF(IV(J).LE.NV) GO TO 316
315 J=J-1
IF(J.GE.1) GO TO 314
IF(IVA(NIND).GT.0) RETURN
IV(NIND)=IV(NIND)+1
IF(IV(NIND).GT.NV) RETURN
K=1
DO 312 I=1,NIND-1
IF(IVA(I).GT.0) GO TO 312
IV(I)=K
K=K+1
312 CONTINUE
J=NIND
GO TO 318
316 K=IV(J)
IF(J.EQ.NIND-1) GO TO 318
DO 317 I=J+1,NIND-1
IF(IVA(I).GT.0) GO TO 317
K=K+1
IF(K.GT.NV) GO TO 315
IV(I)=K
317 CONTINUE
318 DO 313 I=1,NIND-1
DO 313 K=I+1,NIND
IF(IV(I).EQ.IV(K)) GO TO 300
313 CONTINUE
C
C ********************************************************
C
14 DO 15 I=1,NIND
K=IV(I)
AVG(I)=VMN(K)
SIGMA(I)=STD(K)*SQRT(NC-1.)
DO 15 J=1,I
L=IV(J)
DCOR(I,J)=COR(K,L)
15 DCOR(J,I)=DCOR(I,J)
65 LDVAR=IV(NIND)
ISTRE=O
KDV1=0
IDW=1
SDEV2=0
IDIF=0
KTR=0
SDIF2=0.0
IF(IOUT.NE.21) WRITE(IOUT,5566)(HEDR(I),I=1,NSZ)
5566 FORMAT('1',70A1)
IF(IOUT.EQ.21) CALL PRNTHD
WRITE(IOUT,67)
67 FORMAT('0',10X,'***** STEPWISE REGRESSION *****')
WRITE(IOUT,66)NIND,NAMES(IV(NIND))
66 FORMAT('0',I2,' VARIABLES; VARIABLE: ',A5,' IS DEPENDENT')
LINES=6
DEFR=NC-1.
IF(OPT(3).NE.1) GO TO 76
LINES=LINES+(NOFRC+4)/5
WRITE(IOUT,70)NOFRC,(NAMES(IV(INDEX(I))),I=1,NOFRC)
70 FORMAT('0',I2,' VARIABLES FORCED INTO REGRESSION: ',
15(A5,1X)/(38X,5(A5,1X)))
C
C SOLVE SYSTEM OF EQUATIONS FOR FORCED VARIABLES
C
DO 71 L=1,NOFRC
K=INDEX(L)
DIAG=DCOR(K,K)
DO 72 I=1,NIND
IF(I.EQ.K) GO TO 72
FACTR=DCOR(I,K)/DIAG
DO 73 J=1,NIND
IF(J.EQ.K) GO TO 73
DCOR(I,J)=DCOR(I,J)-(FACTR*DCOR(K,J))
73 CONTINUE
72 CONTINUE
DO 74 I=1,NIND
74 DCOR(I,K)=-DCOR(I,K)/DIAG
DO 75 J=1,NIND
75 DCOR(K,J)=DCOR(K,J)/DIAG
71 DCOR(K,K)=1.0/DIAG
FACTR=NOFRC
DEFR=DEFR-FACTR
FLEVL=(DCOR(K,NIND)**2/DIAG)*DEFR/DCOR(NIND,NIND)
NOENT=K
NOSTP=NOFRC-1
GO TO 77
76 NOSTP=-1
C
C SELECT VARIABLES BY ANAL OF VAR AND SOLVE EQUATIONS
C
77 NOSTP=NOSTP+1
IF(DCOR(NIND,NIND).LT.0.0) WRITE(IOUT,78)
78 FORMAT('0LAST DIAGONAL ELEMENT IS NEGATIVE')
IF(DCOR(NIND,NIND).EQ.0.0) WRITE(IOUT,79)
79 FORMAT('0LAST DIAGONAL ELEMENT IS ZERO')
IF(DCOR(NIND,NIND).LE.0.0) GO TO 300
SIGY=SIGMA(NIND)*SQRT(DCOR(NIND,NIND)/DEFR)
DEFR=DEFR-1.0
IF(DEFR.GT.0.0) GO TO 80
WRITE(IDLG,81)DEFR
81 FORMAT('0DEGREES OF FREEDOM IS',F8.2)
GO TO 300
80 VMIN=-10000.0
VMAX=0.0
NOIN=0
DO 85 I=1,NIND-1
DIAG=DCOR(I,I)
IF(DIAG.GT.0.0) GO TO 83
IF(DIAG+TOL.GE.0) GO TO 85
WRITE(IDLG,82)NAMES(IV(I))
82 FORMAT('0NEGATIVE DIAGONAL ELEMENT BELOW TOLERANCE VARIABLE: ',A5)
GO TO 300
83 IF((DIAG-TOL).LT.0.0) GO TO 85
VAR=DCOR(I,NIND)*DCOR(NIND,I)/DIAG
SIG=SIGMA(I)
IF(VAR.EQ.0.0) GO TO85
IF(VAR.GT.0.0) GO TO 89
IF(SIG.LT.0.0)GO TO 86
IF(SIG.GT.0.0) GO TO 87
90 WRITE(IDLG,91)SIG,NAMES(IV(I))
91 FORMAT('0STANDARD DEVIATION IS',G,' FOR VAR ',A5)
GO TO 300
86 SIG=-SIG
GO TO 88
87 IF(VAR.LE.VMIN) GO TO 88
VMIN=VAR
NOMIN=I
88 NOIN=NOIN+1
INDEX(NOIN)=I
COEN(NOIN)=DCOR(I,NIND)*SIGMA(NIND)/SIG
SIGCO(NOIN)=(SIGY/SIG)*SQRT(DIAG)
GO TO 85
89 IF(SIG.LE.0.0) GO TO 90
IF(VAR.LE.VMAX) GO TO 85
VMAX=VAR
NOMAX=I
85 CONTINUE
IF(NOIN.LT.0) GO TO 110
IF(NOIN.GT.0) GO TO 94
KTR=KTR+1
IF(IOUT.NE.21) GO TO 400
LINES=LINES+2
IF(LINES.LE.LINPP) GO TO 400
CALL PRNTHD
LINES=4
400 WRITE(IOUT,92) SIGY
92 FORMAT('0STANDARD ERROR OF Y =',G)
IF(KTR.LE.1) GO TO 96
113 IF(IOUT.NE.21) GO TO 401
LINES=LINES+2
IF(LINES.LE.LINPP) GO TO 401
CALL PRNTHD
LINES=4
401 WRITE(IOUT,93)
93 FORMAT('0ERROR, NO VARIABLE CAN ENTER THE EQUATION')
RETURN
94 CNST=AVG(NIND)
DO 95 J=1,NOIN
L=INDEX(J)
95 CNST=CNST-(COEN(J)*AVG(L))
FLEV=FLEVL
NOEN=NOENT
K1=K
IF(VMIN.EQ.-10000.0) GO TO 96
FLEVL=-VMIN*(DEFR+1.0)/DCOR(NIND,NIND)
IF(EFOUT.LE.FLEVL) GO TO 96
K=NOMIN
DEFR=DEFR+2.0
NOENT=0
GO TO 97
96 FLEVL=VMAX*DEFR/(DCOR(NIND,NIND)-VMAX)
IF(EFIN.GE.FLEVL) GO TO 112
K=NOMAX
NOENT=K
97 COEFDT=1.0-DCOR(NIND,NIND)
COD=COEFDT-COE
COE=COEFDT
RRR=SQRT(COEFDT)
IF(NOIN.LT.0) GO TO 110
IF(NOIN.EQ.0) GO TO 104
IF(NOEN.LE.0) GO TO 99
IF(IOUT.NE.21) GO TO 402
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 402
CALL PRNTHD
LINES=8
402 WRITE(IOUT,98)NOSTP,NAMES(IV(K1))
98 FORMAT(///'0STEP NO.',I3/' ENTERING VARIABLE: ',A5)
ISTRE=K1
GO TO 102
99 IF(K1.NE.ISTRE) GO TO 100
IF(K1.NE.0) GO TO 119
100 IF(IOUT.NE.21) GO TO 403
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 403
CALL PRNTHD
LINES=8
403 WRITE(IOUT,101)NOSTP,NAMES(IV(K1))
101 FORMAT(///'0STEP NO.',I3/' VARIABLE: ',A5,' REMOVED')
102 NDGN=1
NDGD=NC-NOIN+1
PROB=FISHER(NDGN,NDGD,FLEV)
IF(IOUT.NE.21) GO TO 404
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 404
CALL PRNTHD
LINES=8
404 WRITE(IOUT,103)FLEV,PROB,SIGY,COEFDT,RRR,COD,CNST
103 FORMAT(3X,'F-LEVEL',2X,G,' WITH PROB. ',F7.4/
13X,'STANDARD ERROR OF',' ESTIMATE',G/
23X,'COEFFICIENT OF DETERMINATION =',G/
23X,'COEFFICIENT OF MULTIPLE REGRESSION =',G/
33X,'INCREASE IN COEFFICIENT OF DETERMINATION =',G/
43X,'CONSTANT',G13.5)
IF(IOUT.NE.21) GO TO 407
LINES=LINES+2
IF(LINES.LE.(LINPP-5)) GO TO 407
CALL PRNTHD
LINES=4
407 WRITE(IOUT,501)
501 FORMAT('0VARIABLE',3X,'COEFFICIENT',5X,'STD ERROR OF COEFF.')
DO 502 J=1,NOIN
IF(IOUT.NE.21) GO TO 408
LINES=LINES+1
IF(LINES.LE.LINPP) GO TO 408
CALL PRNTHD
WRITE(IOUT,501)
LINES=5
408 WRITE(IOUT,503) NAMES(IV(INDEX(J))),COEN(J),SIGCO(J)
503 FORMAT(2X,A5,5X,G15.7,1X,G15.7)
502 CONTINUE
104 DIAG=DCOR(K,K)
DO 105 I=1,NIND
IF(I.EQ.K) GO TO 105
FACTR=DCOR(I,K)/DIAG
DO 106 J=1,NIND
IF(J.EQ.K) GO TO 106
DCOR(I,J)=DCOR(I,J)-(FACTR*DCOR(K,J))
106 CONTINUE
105 CONTINUE
DO 107 I=1,NIND
107 DCOR(I,K)=-DCOR(I,K)/DIAG
DO 108 J=1,NIND
108 DCOR(K,J)=DCOR(K,J)/DIAG
DCOR(K,K)=1.0/DIAG
IF(NOSTP.LT.250) GO TO 77
WRITE(IDLG,109) NOSTP
109 FORMAT('0ERROR, NUMBER OF STEPS =',I4,' CHECK F-VALUES')
GO TO 300
110 WRITE(IDLG,111)NOIN
111 FORMAT('0ERROR, NOIN IS',I7)
GO TO 300
C
C FINAL OUTPUT
C
112 IF(NOIN.EQ.0) GO TO 113
DO 114 I=1,NOIN
J=INDEX(I)
STUDT(I)=COEN(I)/SIGCO(I)
114 BETA(I)=DCOR(J,NIND)
IF(NOEN.GE.0) GO TO 116
WRITE(IDLG,115)NOEN
115 FORMAT('0ERROR, NOEN=',I6)
GO TO 300
116 IF(NOEN.EQ.0) GO TO 118
IF(IOUT.NE.21) GO TO 405
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 405
CALL PRNTHD
LINES=8
405 WRITE(IOUT,98)NOSTP,NAMES(IV(K1))
ISTRE=K1
GO TO 123
118 IF(K1.NE.ISTRE) GO TO 121
IF(L1.EQ.0) GO TO 121
119 WRITE(IDLG,120)
120 FORMAT(' VARIABLE REMOVED SAME AS VARIABLE ENTERED')
GO TO 300
121 IF(IOUT.NE.21) GO TO 406
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 406
CALL PRNTHD
LINES=8
406 WRITE(IOUT,101) NOSTP,NAMES(IV(K1))
123 COEFDT=1.0-DCOR(NIND,NIND)
COD=COEFDT-COE
RRR=SQRT(COEFDT)
NDGN=1
NDGD=NC-NOIN+1
PROB=FISHER(NDGN,NDGD,FLEV)
IF(IOUT.NE.21) GO TO 409
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 409
CALL PRNTHD
LINES=8
409 WRITE(IOUT,103)FLEV,PROB,SIGY,COEFDT,RRR,COD,CNST
IF(IOUT.NE.21) GO TO 410
LINES=LINES+3
IF(LINES.LE.(LINPP-5)) GO TO 410
CALL PRNTHD
LINES=5
410 WRITE(IOUT,504)
504 FORMAT('0',24X,'STD.ERR.',5X,'STANDARDIZED'/
1' VARIABLE',4X,'COEFF.',6X,'OF COEFF.',4X,'COEFFICIENT',
25X,'T-VALUE',5X,'PROB.')
DO 1100 I=1,NOIN
NDGN=1
NDGD=NC-NOIN-1
PROB=STUDT(I)**2
PROB=FISHER(NDGN,NDGD,PROB)
IF(IOUT.NE.21) GO TO 1100
LINES=LINES+1
IF(LINES.LE.LINPP) GO TO 1100
CALL PRNTHD
WRITE(IOUT,504)
LINES=6
1100 WRITE(IOUT,1101)NAMES(IV(INDEX(I))),COEN(I),SIGCO(I),BETA(I),
1STUDT(I),PROB
1101 FORMAT(2X,A5,4X,G12.4,1X,G12.4,1X,G,1X,G11.4,1X,F6.3)
IF((OPT(4).NE.1).AND.(OPT(5).NE.1).AND.(OPT(6).NE.1)) GO TO 300
IF(OPT(6).NE.1) GO TO 135
DO 200 I=1,NV
200 CO(I)=0
IF(INDXKK.LE.NV) GO TO 135
INDXKK=NV+1
NV=INDXKK
135 NAMES(INDXKK)='RESID'
SSS=0
SST=0
DO 125 J=1,NC
YPRED=CNST
SSS=DATA(J,IV(NIND))+SSS
SST=SST+DATA(J,IV(NIND))**2
DO 126 I=1,NOIN
K=IV(INDEX(I))
126 YPRED=YPRED+COEN(I)*DATA(J,K)
DEV=DATA(J,IV(NIND))-YPRED
SDEV2=SDEV2+DEV**2
IF(IDIF.EQ.0) GO TO 127
DIF=DEV-DEVP
SDIF2=SDIF2+DIF*DIF
GO TO 128
127 IDIF=1
128 DEVP=DEV
IF(OPT(6).NE.1)GO TO 125
DATA(J,INDXKK)=DEV
DO 201 I=1,NV
201 CO(I)=DATA(J,I)*DEV+CO(I)
125 CONTINUE
IF(OPT(6).NE.1) GO TO 202
VMN(INDXKK)=0
STD(INDXKK)=SQRT(CO(INDXKK)/(NC-1.))
DO 203 I=1,NV
IF(I.EQ.INDXKK) GO TO 205
IF((STD(I)*STD(INDXKK)).EQ.0) GO TO 204
COR(I,INDXKK)=CO(I)/((NC-1.)*STD(I)*STD(INDXKK))
GO TO 203
204 COR(I,INDXKK)=0
GO TO 203
205 COR(I,INDXKK)=1.
203 COR(INDXKK,I)=COR(I,INDXKK)
202 IF(OPT(5).NE.1) GO TO 132
SST=SST-SSS*SSS/NC
SSR=SST-SDEV2
NERR=NC-(NOIN+1)
NTOTL=NERR+NOIN
ERMS=SDEV2/NERR
RMS=SSR/NOIN
F=RMS/ERMS
PROB=FISHER(NOIN,NERR,F)
IF(IOUT.NE.21) GO TO 411
LINES=LINES+6
IF(LINES.LE.LINPP) GO TO 411
CALL PRNTHD
LINES=8
411 WRITE(IOUT,130)
130 FORMAT('0ANALYSIS OF VARIANCE')
WRITE(IOUT,131)SSR,NOIN,RMS,F,PROB,SDEV2,NERR,ERMS,SST,NTOTL
131 FORMAT(4X,'SOURCE',6X,'SUM OF SQ.',6X,'DF',5X,
1'MEAN SQ.',6X,'F',7X,'PROB.'/' REGRESSION',F17.6,I6,F15.5,F9.3
2,2X,F7.4/6X,'ERROR',F17.6,I6,F15.5/6X,'TOTAL',F17.6,I6)
132 IF(OPT(4).NE.1) GO TO 300
DWS=SDIF2/SDEV2
DW1=4.-DWS
IF(IOUT.NE.21) GO TO 412
LINES=LINES+7
IF(LINES.LE.LINPP) GO TO 412
CALL PRNTHD
LINES=9
412 WRITE(IOUT,133)SDIF2,DWS,DW1
133 FORMAT(//14X,'DURBIN-WATSON TEST FOR AUTOCORRRELATION'//
110X,'SUM OF SQUARED DIFFERENCES =',G/10X,'DWS=',G/
210X,'4 - DWS =',G)
GO TO 300
END
C *** STAT PACK ***
C SUBROUTINE FOR PARTIAL CORRELATIONS
C CALLING SEQUENCE: CALL PCORR(NV,NC,MV,MC,COR,SP,NAMES)
C WHERE NV - IS THE NUMBER OF COLUMNS ACTUALLY FILLED(VARIABLES)
C NC - IS THE NUMBER OF ROWS ACTUALLY FILLED(OBSERVATIONS,CASES)
C MV - IS MAXIMUM NUMBER OF COLUMNS, AS SPECIFIED IN MAIN.
C MC - IS MAXIMUM NUMBER OF ROWS, AS SPECIFIED IN MAIN.
C COR - IS THE CORRELATION MATRIX.
C SP - IS A VECTOR DIMENSIONED AT LEAST NV.
C NAMES - IS A VECTOR CONTAINING VARIABLE NAMES
C
C GIVES PARTIAL CORRELATIONS MATRIX FOR SPECIFIED VARIABLES.
C
SUBROUTINE PCORR(NV,NC,MV,MC,COR,SP,NAMES)
COMMON /DEV/ICC,IDATA,IOUT,IDLG,IDSK
COMMON /PRNT/ LINPP,ICOPS,RUNPRG
COMMON/EXTRA/HEDR(70),NSZ
DIMENSION COR(MV,MV),SP(1),DCOR(20,20),XVER(20,20)
DIMENSION IV(20),NAMES(1),IVA(20)
ISQ=7
IF(IOUT.EQ.21) ISQ=15
1 IF(ICC.NE.2) WRITE(IDLG,2)
2 FORMAT('0WHICH VARIABLES? ',$)
CALL ALPHA(IVA,20,NN,IRET,IHELP,IERR,NAMES,NV)
IF(IRET.EQ.1) RETURN
IF(IERR.EQ.1) GO TO 1
IF(IHELP.EQ.1) GO TO 1
K=1
DO 3 I=1,NN
IV(I)=IVA(I)
IF(IVA(I).GT.0) GO TO 3
IV(I)=K
K=K+1
3 CONTINUE
DO 4 I=1,NN-1
IF(IVA(I).LT.0) GOTO 4
DO 5 J=I+1,NN
IF(IVA(I).LT.0) GO TO 5
IF(IVA(I).NE.IVA(J)) GO TO 5
WRITE(IDLG,6)
6 FORMAT(' THE SAME VARIABLE HAS BEEN LISTED TWICE')
GO TO 1
5 CONTINUE
4 CONTINUE
IF(NN.LE.NV) GO TO 58
WRITE(IDLG,7)
7 FORMAT(' YOU HAVE LISTED MORE VARIABLES THAN IS POSSIBLE'/
1' WITH THE DATA SET AVAILABLE')
GO TO 1
50 J=NN
54 IF(IVA(J).GT.0) GOTO 55
IV(J)=IV(J)+1
IF(IV(J).LE.NV) GO TO 56
55 J=J-1
IF(J.GE.1) GO TO 54
RETURN
56 K=IV(J)
IF(J.EQ.NN) GOTO 58
DO 57 I=J+1,NN
IF(IVA(I).GT.0) GO TO 57
K=K+1
IF(K.GT.NV) GO TO 55
IV(I)=K
57 CONTINUE
58 DO 53 I=1,NN-1
DO 53 K=I+1,NN
IF(IV(I).EQ.IV(K)) GO TO 50
53 CONTINUE
DO 8 I=1,NN
DO 8 J=1,NN
8 DCOR(I,J)=COR(IV(I),IV(J))
C
C *****************************************************************
C INVERSE BY LINEAR ROW TRANSFORMATION
C
DO 9 I=1,NN
DO 10 J=1,NN
10 XVER(I,J)=0
9 XVER(I,I)=1.0
DO 11 I=1,NN
IF(((DCOR(I,I)+100.)-100.).NE.0.0) GO TO 15
IF(I.EQ.NN) GO TO 30
DO 12 J=I+1,NN
IF(((DCOR(J,I)+100.)-100.).NE.0.0) GO TO 13
12 CONTINUE
30 WRITE(IDLG,31)
31 FORMAT('0THE INVERSE DOES NOT EXIST')
GO TO 50
13 DO 14 K=1,NN
DCOR(I,K)=DCOR(I,K)+DCOR(J,K)
14 XVER(I,K)=XVER(I,K)+XVER(J,K)
15 G=DCOR(I,I)
DO 16 J=1,NN
DCOR(I,J)=DCOR(I,J)/G
16 XVER(I,J)=XVER(I,J)/G
DO 17 L=1,NN
IF(L.EQ.I) GO TO 17
G=DCOR(L,I)
DO 18 J=1,NN
DCOR(L,J)=DCOR(L,J)-G*DCOR(I,J)
18 XVER(L,J)=XVER(L,J)-G*XVER(I,J)
17 CONTINUE
11 CONTINUE
C
C ***************************************************************
C
IF(IOUT.NE.21) WRITE(IOUT,5566)(HEDR(I),I=1,NSZ)
5566 FORMAT('1',70A1)
IF(IOUT.EQ.21) CALL PRNTHD
WRITE(IOUT,19)
19 FORMAT('0',10X,'***** PARTIAL CORRELATION MATRIX *****')
LINES=4
DO 20 I=2,NN
DO 21 J=1,I-1
21 SP(J)=-XVER(I,J)/SQRT(XVER(I,I)*XVER(J,J))
IF(IOUT.NE.21) GO TO 27
L=(I+ISQ-2)/ISQ+1
LINES=LINES+L
IF(LINES.LE.(LINPP-L-1)) GO TO 27
WRITE(IOUT,25)
DO 28 J=1,I-2,ISQ
NEND=J+ISQ-1
IF(NEND.GT.(I-2)) NEND=I-2
28 WRITE(IOUT,23)(NAMES(IV(J)),K=J,NEND)
CALL PRNTHD
LINES=2+L
27 DO 20 K=1,I-1,ISQ
NEND=K+ISQ-1
IF(NEND.GT.I-1) NEND=I-1
IF(K.EQ.1)WRITE(IOUT,22) NAMES(IV(I)),(SP(J),J=K,NEND)
22 FORMAT('0',A5,1X,15F8.4)
IF(K.NE.1) WRITE(IOUT,26) (SP(J),J=K,NEND)
26 FORMAT(7X,15F8.4)
20 CONTINUE
WRITE(IOUT,25)
25 FORMAT(1X)
DO 24 I=1,NN,ISQ
NEND=I+ISQ-1
IF(NEND.GT.NN-1) NEND=NN-1
24 WRITE(IOUT,23)(NAMES(IV(J)),J=I,NEND)
23 FORMAT(7X,15(3X,A5))
GO TO 50
END
C SUBROUTINE FOR WILCOXON MATCHED PAIRS SIGNED-RANK TEST.
C CALLING SEQUENCE: CALL WILCX(NV,NC,MV,MC,DATA,IV,H,NAMES)
C WHERE NV - IS THE NUMBER OF COLUMNS ACTUALLY FILLED(VARIABLES)
C NC - IS THE NUMBER OF ROWS ACTUALLY FILLED(OBSERVATIONS,CASES)
C MV - IS MAXIMUM NUMBER OF COLUMNS, AS SPECIFIED IN MAIN.
C MC - IS MAXIMUM NUMBER OF ROWS, AS SPECIFIED IN MAIN.
C DATA - STORAGE FOR DATA DIMENSIONED FOR MAXIMUM MATRIX.
C IV - IS A VECTOR DIMENSIONED AT LEAST NC.
C H - IS A VECTOR DIMENSIONED FOR AT LEAST NC.
C NAMES - IS A VECTOR CONTAINING VARIABLE NAMES
C
C PROGRAM AUTOMATICALLY OUTPUTS ENTIRE MATRIX, FORMAING IT ONE LINE
C AT A TIME. MAKES USE OF A PARTITIONING SORT.
C
SUBROUTINE WILCX(NV,NC,MV,MC,DATA,IV,H,NAMES)
DIMENSION DATA(MC,MV),IV(1),H(1),NAMES(1)
C FOLLOWING DIMENSION IS FOR USE IN THE SORT
DIMENSION IU(16),IL(16)
COMMON/DEV/ICC,IDATA,IOUT,IDLG,IDSK
COMMON /PRNT/ LINPP,ICOPS,RUNPRG
COMMON/EXTRA/HEDR(70),NSZ
INTEGER BASE
IF(IOUT.NE.21) WRITE(IOUT,5566) (HEDR(I),I=1,NSZ)
5566 FORMAT('1',70A1)
IF(IOUT.EQ.21) CALL PRNTHD
WRITE(IOUT,1)
1 FORMAT('0',10X,'***** WILCOXON MATCHED-PAIRS SIGNED-RANK TEST',
1' *****'/'0VAR. VS VAR.',7X,'MEAN',10X,'S.D.',4X,'N',6X
2,'Z',10X,'T')
LINES=6
DO 2 I=1,NV-1
DO 3 JJ=I+1,NV
L=1
DO 4 K=1,NC
IF(DATA(K,I)-DATA(K,JJ))5,4,6
5 IV(L)=1
H(L)=DATA(K,JJ)-DATA(K,I)
GO TO 7
6 IV(L)=0
H(L)=DATA(K,I)-DATA(K,JJ)
7 L=L+1
4 CONTINUE
C
C SORT BY PARTITION - ACM
C
NN=L-1
IF(NN.LT.2) GO TO 27
M=1
II=1
J=NN
71 IF(II.GE.J) GO TO 78
72 K=II
IJ=(J+II)/2
T=H(IJ)
IF(H(II).LE.T) GO TO 73
H(IJ)=H(II)
H(II)=T
T=H(IJ)
ISAV=IV(IJ)
IV(IJ)=IV(II)
IV(II)=ISAV
73 LL=J
IF(H(J).GE.T) GO TO 75
H(IJ)=H(J)
H(J)=T
T=H(IJ)
ISAV=IV(IJ)
IV(IJ)=IV(J)
IV(J)=ISAV
IF(H(II).LE.T) GO TO 75
H(IJ)=H(II)
H(II)=T
T=H(IJ)
ISAV=IV(IJ)
IV(IJ)=IV(II)
IV(II)=ISAV
GO TO 75
74 H(LL)=H(K)
H(K)=TT
ISAV=IV(LL)
IV(LL)=IV(K)
IV(K)=ISAV
75 LL=LL-1
IF(H(LL).GT.T) GO TO 75
TT=H(LL)
76 K=K+1
IF(H(K).LT.T) GO TO 76
IF(K.LE.LL) GO TO 74
IF((LL-II).LE.(J-K)) GO TO 77
IL(M)=II
IU(M)=LL
II=K
M=M+1
GO TO 79
77 IL(M)=K
IU(M)=J
J=LL
M=M+1
GO TO 79
78 M=M-1
IF(M.EQ.0) GO TO 20
II=IL(M)
J=IU(M)
79 IF((J-II).GE.11) GO TO 72
IF(II.EQ.1) GO TO 71
II=II-1
80 II=II+1
IF(II.EQ.J) GO TO 78
T=H(II+1)
ISAV=IV(II+1)
IF(H(II).LE.T) GO TO 80
K=II
81 H(K+1)=H(K)
IV(K+1)=IV(K)
K=K-1
IF(T.LT.H(K)) GO TO 81
H(K+1)=T
IV(K+1)=ISAV
GO TO 80
C
C END SORT BEGIN RANKING
C
20 SAV=H(1)
SUM=1
IVV=1
DO 21 K=2,NN
IF(H(K).EQ.SAV) GO TO 23
H(K-1)=SUM/IVV
DO 22 L=K-IVV,K-1
22 H(L)=H(K-1)
SAV=H(K)
SUM=K
IVV=1
GO TO 21
23 IVV=IVV+1
SUM=SUM+K
21 CONTINUE
H(NN)=SUM/IVV
DO 24 L=NN-IVV+1,NN
24 H(L)=H(NN)
C
C RANKING COMPLETE
C
C FIND NUMBER OF ORIGINAL NEGATIVES
C
27 SD=0
Z=0
SAV=0
IVV=0
SUM=0
IF(NN.LT.2) GO TO 28
DO 25 L=1,NN
IF(IV(L).EQ.0) GO TO 25
SUM=SUM+H(L)
IVV=IVV+1
25 CONTINUE
SAV=(NN*(NN+1))/2
IF((IVV*2).LT.NN) GO TO 26
C
C IF IVE COUNTED THE WROND SIGN ADJUST THE SUM BY KNOWING
C THE TOTAL OF THE RANKS IS N(N+1)/2
C
SUM=SAV-SUM
IVV=NN-IVV
26 SAV=SAV/2
SD=SQRT(((2.*NN+1.)*SAV)/6.)
Z=(SUM-SAV)/SD
28 IF(IOUT.NE.21) GO TO 31
LINES=LINES+1
IF(LINES.LE.LINPP) GO TO 31
CALL PRNTHD
WRITE(IOUT,1)
LINES=7
31 WRITE(IOUT,30) NAMES(I),NAMES(JJ),SAV,SD,NN,Z,SUM
30 FORMAT(1X,A5,2X,A5,F11.2,F14.2,I5,F8.2,F11.2)
3 CONTINUE
2 CONTINUE
RETURN
END