1 1 GENLIN: TESTING THE GENERAL LINEAR MODEL DONALD E. RAMIREZ DEPARTMENT OF MATHEMATICS MATH-ASTRO BUILDING UNIVERSITY OF VIRGINIA CHARLOTTESVILLE, VIRGINIA 22903 TECHNICAL REPORT (MAY 12, 1988) 1 INTRODUCTION THIS TECHNICAL REPORT DESCRIBES THE PROGRAM GENLIN WHICH IS PRESENTED IN THE PAPER: GENLIN: TESTING THE GENERAL LINEAR MODEL, MULTIVARIATE BEHAVIORAL RESEARCH (TO APPEAR) GENLIN IS A FORTRAN 77 PROGRAM WHICH SOLVES THE MULTIVARIATE LINEAR REGRESSION MODEL OF THE FORM: Y = X*B + E WHERE E(Y) = X*B AND COV(Y) = IxS. THE USER CAN DELETE COLUMNS OF X TO FORM A NULL HYPOTHESES TO TEST AGAINST THE FULL MODEL. THE OUTPUT CAN INCLUDE INV(X'*X), S^, B^, AND COV(B^). THE TEST STATISTIC FOR TESTING THE NULL HYPOTHESIS IS WILKS' LAMBDA CRITERION. THE DESCRIPTIVE LEVEL IS ALSO GIVEN. THE USER IS GIVEN A CHANCE TO TEST RESTRICTIVE HYPOTHESES OF THE FORM: H: C*B = 0. INDEX APPENDIX A CONTAINS THE RAW DATA SET WITH N=37 ,Q=7, AND P=3. APPENDIX B CONTAINS THE TERMINAL OUTPUT FROM GENLIN APPENDIX C CONTAINS THE FORTRAN CODE FOR GENLIN APPENDIX D CONTAINS THE CALL STATEMENTS FROM GENLIN $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ APPENDIX A: THE ORIGINAL DATA 1 1 2 6 12 16 49 48 8 1 5 14 14 30 27 47 76 13 1 0 10 21 16 16 11 40 13 1 0 2 5 17 8 9 52 9 1 2 7 11 26 17 69 63 15 1 2 15 21 34 25 35 82 14 1 0 1 20 23 18 6 71 21 1 0 0 10 19 14 8 68 8 1 0 0 7 16 13 49 74 11 1 3 2 21 26 25 8 70 15 1 8 16 15 35 24 47 70 15 1 5 4 7 15 14 6 61 11 1 1 12 13 27 21 14 54 12 1 2 1 12 20 17 30 55 13 1 3 12 20 26 22 4 54 10 1 0 2 5 14 8 24 40 14 1 7 12 21 35 27 19 66 13 1 0 6 6 14 16 45 54 10 1 12 8 19 27 26 22 64 14 1 3 9 15 18 10 16 47 16 1 0 7 9 14 18 32 48 16 1 4 6 20 26 26 37 52 14 1 4 9 14 23 23 47 74 19 1 0 2 4 11 8 5 57 12 1 0 1 16 15 17 6 57 10 1 3 8 18 28 21 60 80 11 1 1 18 19 34 23 58 78 13 1 2 11 9 23 11 6 70 16 1 0 10 7 12 8 16 47 14 1 8 10 28 32 32 45 94 19 1 2 12 5 25 14 9 63 11 1 7 11 18 29 21 69 76 16 1 2 5 10 23 24 35 59 11 1 0 1 14 19 12 19 55 8 1 1 0 10 18 18 58 74 14 1 6 4 23 31 26 58 71 17 1 0 6 6 15 14 79 54 14 N=37, Columns of X = 6, Columns of Y = 3 REFERENCE TIMM, N. (1975) MULTIVARIATE ANALYSIS WITH APPLICATIONS IN EDUCATION AND PSYCHOLOGY. MONTREREY, CALIFORNIA: BROOKS/COLE PUBLISHING COMPANY. $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ APPENDIX B: TERMINAL OUTPUT FROM GENLIN OK, SEG GENLIN ENTER THE NAME OF THE DATA FILE (WITH X BEFORE Y) TIMM_P314.DATA ENTER THE NUMBER OF CASES : <= 200 37 ENTER THE NUMBER OF COLUMNS OF X : <= 10 6 ENTER THE NUMBER OF COLUMNS OF Y : <= 10 3 DO YOU WANT TO SEE THE ALTERNATIVE DATA? (1=YES) 1 EFFECT: 1 2 3 4 5 6 ROW 1 1.00 1.00 2.00 6.00 12.00 16.00 2 1.00 5.00 14.00 14.00 30.00 27.00 3 1.00 0.00 10.00 21.00 16.00 16.00 4 1.00 0.00 2.00 5.00 17.00 8.00 5 1.00 2.00 7.00 11.00 26.00 17.00 6 1.00 2.00 15.00 21.00 34.00 25.00 7 1.00 0.00 1.00 20.00 23.00 18.00 8 1.00 0.00 0.00 10.00 19.00 14.00 9 1.00 0.00 0.00 7.00 16.00 13.00 10 1.00 3.00 2.00 21.00 26.00 25.00 11 1.00 8.00 16.00 15.00 35.00 24.00 12 1.00 5.00 4.00 7.00 15.00 14.00 13 1.00 1.00 12.00 13.00 27.00 21.00 14 1.00 2.00 1.00 12.00 20.00 17.00 15 1.00 3.00 12.00 20.00 26.00 22.00 16 1.00 0.00 2.00 5.00 14.00 8.00 17 1.00 7.00 12.00 21.00 35.00 27.00 18 1.00 0.00 6.00 6.00 14.00 16.00 19 1.00 12.00 8.00 19.00 27.00 26.00 20 1.00 3.00 9.00 15.00 18.00 10.00 21 1.00 0.00 7.00 9.00 14.00 18.00 22 1.00 4.00 6.00 20.00 26.00 26.00 23 1.00 4.00 9.00 14.00 23.00 23.00 24 1.00 0.00 2.00 4.00 11.00 8.00 25 1.00 0.00 1.00 16.00 15.00 17.00 26 1.00 3.00 8.00 18.00 28.00 21.00 27 1.00 1.00 18.00 19.00 34.00 23.00 28 1.00 2.00 11.00 9.00 23.00 11.00 29 1.00 0.00 10.00 7.00 12.00 8.00 30 1.00 8.00 10.00 28.00 32.00 32.00 31 1.00 2.00 12.00 5.00 25.00 14.00 32 1.00 7.00 11.00 18.00 29.00 21.00 33 1.00 2.00 5.00 10.00 23.00 24.00 34 1.00 0.00 1.00 14.00 19.00 12.00 35 1.00 1.00 0.00 10.00 18.00 18.00 36 1.00 6.00 4.00 23.00 31.00 26.00 37 1.00 0.00 6.00 6.00 15.00 14.00 THERE ARE 6 COLUMNS IN THE ALTERNATIVE HYPOTHESIS. HOW MANY DO YOU WISH TO DELETE FOR THE NULL (AT LEAST 1)? 3 WHICH 3 COLUMNS ARE TO BE DELETED? 2 3 4 WHICH INV(X"X) MATRICES DO YOU WANT TO SEE? (0=NEITHER, 1=ALTERNATIVE, 2=NULL, 3=BOTH) 3 WHICH COV(BETA) MATRICES DO YOU WISH TO SEE ? (0=NEITHER, 1=ALTERNATIVE, 2=NULL, 3=BOTH) 3 NULL HYPOTHESIS ==== ========== BETA ESTIMATES: 10.155 35.864 9.746 -0.033 0.963 0.079 1.189 0.285 0.094 SIGMA ESTIMATES (UNBIASED--N = 37, P = 3 ): 459.405 51.641 7.127 51.641 89.824 4.791 7.127 4.791 8.837 INV(X"X): 0.306 -0.008 -0.005 -0.008 0.001 -0.001 -0.005 -0.001 0.002 COVARIANCE BETA: B11 B12 B13 B21 B22 B23 B31 B32 B33 140.602 15.805 2.181 -3.811 -0.428 -0.059 -2.335 -0.262 -0.036 15.805 27.491 1.466 -0.428 -0.745 -0.040 -0.262 -0.456 -0.024 2.181 1.466 2.705 -0.059 -0.040 -0.073 -0.036 -0.024 -0.045 -3.811 -0.428 -0.059 0.661 0.074 0.010 -0.598 -0.067 -0.009 -0.428 -0.745 -0.040 0.074 0.129 0.007 -0.067 -0.117 -0.006 -0.059 -0.040 -0.073 0.010 0.007 0.013 -0.009 -0.006 -0.012 -2.335 -0.262 -0.036 -0.598 -0.067 -0.009 0.855 0.096 0.013 -0.262 -0.456 -0.024 -0.067 -0.117 -0.006 0.096 0.167 0.009 -0.036 -0.024 -0.045 -0.009 -0.006 -0.012 0.013 0.009 0.016 ALTERNATIVE HYPOTHESIS =========== ========== BETA ESTIMATES: 4.151 33.005 11.173 -0.609 -0.081 0.211 -0.050 -0.721 0.065 -1.732 -0.298 0.214 0.495 1.470 -0.037 2.248 0.324 -0.052 SIGMA ESTIMATES (UNBIASED--N = 37, P = 6 ): 449.338 49.371 14.748 49.371 89.185 6.893 14.748 6.893 8.654 INV(X"X): 0.424 0.024 0.005 0.005 -0.016 -0.011 0.024 0.006 0.000 0.000 -0.001 -0.001 0.005 0.000 0.002 0.000 -0.001 0.000 0.005 0.000 0.000 0.002 -0.001 -0.001 -0.016 -0.001 -0.001 -0.001 0.002 -0.001 -0.011 -0.001 0.000 -0.001 -0.001 0.003 COVARIANCE BETA: B11 B12 B13 B21 B22 B23 B31 B32 B33 B41 B42 B43 B51 B52 B53 B61 B62 B63 190.394 20.920 6.249 10.811 1.188 0.355 2.417 0.266 0.079 2.091 0.230 0.069 -7.036 -0.773 -0.231 -5.071 -0.557 -0.166 20.920 37.790 2.921 1.188 2.146 0.166 0.266 0.480 0.037 0.230 0.415 0.032 -0.773 -1.396 -0.108 -0.557 -1.007 -0.078 6.249 2.921 3.667 0.355 0.166 0.208 0.079 0.037 0.047 0.069 0.032 0.040 -0.231 -0.108 -0.136 -0.166 -0.078 -0.098 10.811 1.188 0.355 2.793 0.307 0.092 -0.051 -0.006 -0.002 0.030 0.003 0.001 -0.333 -0.037 -0.011 -0.572 -0.063 -0.019 1.188 2.146 0.166 0.307 0.554 0.043 -0.006 -0.010 -0.001 0.003 0.006 0.000 -0.037 -0.066 -0.005 -0.063 -0.113 -0.009 0.355 0.166 0.208 0.092 0.043 0.054 -0.002 -0.001 -0.001 0.001 0.000 0.001 -0.011 -0.005 -0.006 -0.019 -0.009 -0.011 2.417 0.266 0.079 -0.051 -0.006 -0.002 0.886 0.097 0.029 0.112 0.012 0.004 -0.543 -0.060 -0.018 0.121 0.013 0.004 0.266 0.480 0.037 -0.006 -0.010 -0.001 0.097 0.176 0.014 0.012 0.022 0.002 -0.060 -0.108 -0.008 0.013 0.024 0.002 0.079 0.037 0.047 -0.002 -0.001 -0.001 0.029 0.014 0.017 0.004 0.002 0.002 -0.018 -0.008 -0.010 0.004 0.002 0.002 2.092 0.230 0.069 0.030 0.003 0.001 0.112 0.012 0.004 0.829 0.091 0.027 -0.271 -0.030 -0.009 -0.439 -0.048 -0.014 0.230 0.415 0.032 0.003 0.006 0.000 0.012 0.022 0.002 0.091 0.165 0.013 -0.030 -0.054 -0.004 -0.048 -0.087 -0.007 0.069 0.032 0.040 0.001 0.000 0.001 0.004 0.002 0.002 0.027 0.013 0.016 -0.009 -0.004 -0.005 -0.014 -0.007 -0.008 -7.036 -0.773 -0.231 -0.333 -0.037 -0.011 -0.543 -0.060 -0.018 -0.271 -0.030 -0.009 1.075 0.118 0.035 -0.477 -0.052 -0.016 -0.773 -1.396 -0.108 -0.037 -0.066 -0.005 -0.060 -0.108 -0.008 -0.030 -0.054 -0.004 0.118 0.213 0.016 -0.052 -0.095 -0.007 -0.231 -0.108 -0.136 -0.011 -0.005 -0.006 -0.018 -0.008 -0.010 -0.009 -0.004 -0.005 0.035 0.016 0.021 -0.016 -0.007 -0.009 -5.071 -0.557 -0.166 -0.572 -0.063 -0.019 0.121 0.013 0.004 -0.439 -0.048 -0.014 -0.477 -0.052 -0.016 1.212 0.133 0.040 -0.557 -1.007 -0.078 -0.063 -0.113 -0.009 0.013 0.024 0.002 -0.048 -0.087 -0.007 -0.052 -0.095 -0.007 0.133 0.241 0.019 -0.166 -0.078 -0.098 -0.019 -0.009 -0.011 0.004 0.002 0.002 -0.014 -0.007 -0.008 -0.016 -0.007 -0.009 0.040 0.019 0.023 DETERMINANT FOR FULL MODEL: 173507. DETERMINANT FOR NULL MODEL: 255701. WILKS LAMBDA( 3,31, 3) : 0.6786 DESCRIPTIVE LEVEL: 0.2232 ===> CANNOT REJECT NULL ========================================================= DO YOU WANT TO TEST A RESTRICTION OF THE FORM H: CB=0 ? (1=YES) 1 IS THE TEST ON THE ALTERNATIVE (1) OR THE NULL (0)? 1 HOW MANY RESTRICTIONS? 3 INPUT C, 3 X 6 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 DETERMINANT FOR FULL MODEL: 173835. DETERMINANT FOR NULL MODEL: 256153. WILKS LAMBDA( 3,31, 3) : 0.6786 DESCRIPTIVE LEVEL: 0.2234 ===> CANNOT REJECT NULL DO YOU WANT TO TEST A RESTRICTION OF THE FORM H: CB=0 ? (1=YES) 0 DO YOU WISH TO MAKE ANOTHER RUN (1=YES) ? 0 **** STOP $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ APPENDIX C: THE FORTRAN PROGRAM GENLIN PROGRAM GENLIN 0000001 C THIS PROGRAM WILL READ IN DATA FROM A DATA FILE. 0000002 C THE USER MAY THEN CHOOSE COLUMNS TO BE DELETED FROM 0000003 C THE DESIGN MATRIX THUS FORMING THE NULL HYPOTHESIS. WILKS 0000004 C LAMBDA TEST WILL BE CALCULATED TO TEST THIS HYPOTHESIS 0000005 C AGAINST THE ALTERNATIVE (FULL DESIGN). 0000006 C AFTER MAKING THIS TEST, THE USER IS ALLOWED TO TEST 0000007 C RESTRICTION HYPOTHESES OF THE FORM H:CB=0. THE USER MAY 0000008 C CHOOSE EITHER THE NULL OR ALTERNATIVE DESIGN MATRICES FROM 0000009 C THE PREVIOUS TEST TO BE USED FOR THESE TESTS. 0000010 C IF THE USER WISHES TO TEST MORE HYPOTHESES BY DELETING 0000011 C MORE OR DIFFERENT COLUMNS HE MAY DO SO. HE IS GIVEN 0000012 C THE OPTION OF SUBSTITUTING THE NULL HYPOTHESIS FROM THE 0000013 C PREVIOUS TEST FOR THE ALTERNATIVE IN THE NEW TEST. 0000014 C DEVELOPED BY DONALD E. RAMIREZ 0000015 C MATH DEPT 0000016 C MATH-ASTRO BUILDING 0000017 C UNIVERSITY OF VIRGINIA 0000018 C CHARLOTTESVILLE, VIRGINIA 22903 0000019 C VERSION 05/12/88 0000020 C PROGRAMMING ASSISTANTS: BRENDA WILSON, LARRY SCHUMANN, DAVID AGARD0000021 INTEGER IA,IB,N,PA,PY,I,J,K,NUMDEL,PO,OCOL,L 0000022 INTEGER P,M,Q,ANSWER,G,MODEL,QT,UN 0000023 C IF IA OR IB ARE CHANGED HERE THEN CHANGE IRB ALSO 0000024 C IN THE SUBROUTINES RTEST, BETSIG, AND WILKS 0000025 PARAMETER (IA=200, IB=10) 0000026 INTEGER DEL(IB) 0000027 REAL XA(IA,IB),XO(IA,IB),BETAA(IB,IB),BETAO(IB,IB),XXA(IB,IB) 0000028 REAL SIGMAA(IB,IB),SIGMAO(IB,IB),LAMBDA,Y(IA,IB),DESCRL 0000029 REAL XXO(IB,IB),C(IB,IB),WORK(IA,IB) 0000030 LOGICAL COPY,XPXA,XPXO,RUN1,CBA,CBO,NEW,OLD 0000031 CHARACTER FN*60,TITLE*80 0000032 EXTERNAL BETSIG,BPRINT,WILKS,RTEST,GTOPFN 0000033 C READ ALTERNATIVE X MATRIX IN FROM DATA FILE AND PLACE IT IN XA. 0000034 UN=5 0000035 NEW=.FALSE. 0000036 OLD=.TRUE. 0000037 TITLE=' ENTER THE NAME OF THE DATA FILE (WITH X BEFORE Y)' 0000038 CALL GTOPFN(UN,FN,NEW,OLD,TITLE) 0000039 PRINT*,'ENTER THE NUMBER OF CASES : <= ',IA 0000040 READ*,N 0000041 PRINT*,'ENTER THE NUMBER OF COLUMNS OF X : <= ',IB 0000042 READ*,PA 0000043 PRINT*,'ENTER THE NUMBER OF COLUMNS OF Y : <= ',IB 0000044 READ*,PY 0000045 DO 10 I=1,N 0000046 READ(5,*) (XA(I,J),J=1,PA),(Y(I,K),K=1,PY) 0000047 10 CONTINUE 0000048 RUN1 = .TRUE. 0000049 20 CONTINUE 0000050 PRINT*,'DO YOU WANT TO SEE THE ALTERNATIVE DATA? (1=YES)' 0000051 READ*,ANSWER 0000052 IF(ANSWER .EQ. 1) THEN 0000053 PRINT*,' ' 0000054 PRINT*,' ' 0000055 PRINT 30, (I,I=1,PA) 0000056 30 FORMAT (' EFFECT: ',10(I1,8X)) 0000057 PRINT*,' ' 0000058 PRINT*,' ROW' 0000059 PRINT*,' ' 0000060 DO 50 I=1,N 0000061 PRINT 40, I, (XA(I,J),J=1,PA) 0000062 40 FORMAT(I3,3X,10(F8.2,1X)) 0000063 50 CONTINUE 0000064 PRINT*,' ' 0000065 PRINT*,' ' 0000066 ENDIF 0000067 C USER CHOOSES COLUMNS TO BE DELETED IN NULL HYPOTHESIS. 0000068 PRINT*,'THERE ARE ',PA,' COLUMNS IN THE ALTERNATIVE HYPOTHESIS.' 0000069 55 CONTINUE 0000070 PRINT*,'HOW MANY DO YOU WISH TO DELETE FOR THE NULL (AT LEAST 1)?'0000071 READ*,NUMDEL 0000072 IF (NUMDEL .LT. 1) GOTO 55 0000073 PRINT*,'WHICH ',NUMDEL,' COLUMNS ARE TO BE DELETED?' 0000074 READ*,(DEL(I),I=1,NUMDEL) 0000075 C PRINT INV(X'X) ? 0000076 PRINT*,'WHICH INV(X"X) MATRICES DO YOU WANT TO SEE?' 0000077 PRINT*,'(0=NEITHER, 1=ALTERNATIVE, 2=NULL, 3=BOTH)' 0000078 READ*,ANSWER 0000079 IF ((ANSWER.EQ.1).OR.(ANSWER.EQ.3)) THEN 0000080 XPXA = .TRUE. 0000081 ELSE 0000082 XPXA = .FALSE. 0000083 ENDIF 0000084 IF ((ANSWER .EQ.2).OR.(ANSWER.EQ.3)) THEN 0000085 XPXO = .TRUE. 0000086 ELSE 0000087 XPXO = .FALSE. 0000088 ENDIF 0000089 PRINT*,'WHICH COV(BETA) MATRICES DO YOU WISH TO SEE ?' 0000090 PRINT*,'(0=NEITHER, 1=ALTERNATIVE, 2=NULL, 3=BOTH)' 0000091 READ*,ANSWER 0000092 IF ((ANSWER.EQ.1).OR.(ANSWER.EQ.3)) THEN 0000093 CBA = .TRUE. 0000094 ELSE 0000095 CBA = .FALSE. 0000096 ENDIF 0000097 IF ((ANSWER.EQ.2).OR.(ANSWER.EQ.3)) THEN 0000098 CBO = .TRUE. 0000099 ELSE 0000100 CBO=.FALSE. 0000101 ENDIF 0000102 C FORM THE NULL X MATRIX, XO 0000103 PO=PA-NUMDEL 0000104 OCOL=0 0000105 DO 120 I=1,PA 0000106 COPY=.TRUE. 0000107 L=NUMDEL 0000108 DO 100 J=1,L 0000109 IF (I.EQ.DEL(J)) THEN 0000110 COPY = .FALSE. 0000111 L=J 0000112 ENDIF 0000113 100 CONTINUE 0000114 IF (COPY) THEN 0000115 OCOL=OCOL+1 0000116 DO 110 J=1,N 0000117 XO(J,OCOL)=XA(J,I) 0000118 110 CONTINUE 0000119 ENDIF 0000120 120 CONTINUE 0000121 C CALCULATE THE BETA MATRICES: BETAA(PA,PY), BETAO(PO.PY) 0000122 C AND THE SIGMA MATRICES: SIGMAA(PY,PY), SIGMAO(PY,PY) 0000123 IF (RUN1) THEN 0000124 CALL BETSIG(Y,XA,N,PY,PA,BETAA,SIGMAA,IA,IB,XXA,WORK) 0000125 ENDIF 0000126 CALL BETSIG(Y,XO,N,PY,PO,BETAO,SIGMAO,IA,IB,XXO,WORK) 0000127 C PRINT OUT BETA AND SIGMA MATRICES 0000128 PRINT*,' ' 0000129 PRINT*,' ' 0000130 PRINT*,'NULL HYPOTHESIS' 0000131 PRINT*,'==== ==========' 0000132 CALL BPRINT(BETAO,SIGMAO,PO,PY,N,IB,XXO,XPXO,CBO) 0000133 PRINT*,' ' 0000134 PRINT*,' ' 0000135 PRINT*,'ALTERNATIVE HYPOTHESIS' 0000136 PRINT*,'=========== ==========' 0000137 CALL BPRINT(BETAA,SIGMAA,PA,PY,N,IB,XXA,XPXA,CBA) 0000138 C CALCULATE WILKS LAMBDA AND THE DESCRIPTIVE LEVEL. WILKS HAS 0000139 C THE PARAMETERS P,M,Q. 0000140 P=PY 0000141 M=N-PA 0000142 Q=NUMDEL 0000143 CALL WILKS(SIGMAA,SIGMAO,P,M,Q,PY,IB,LAMBDA,DESCRL) 0000144 PRINT*,' ' 0000145 PRINT*,' ' 0000146 PRINT 130,P,M,Q,LAMBDA 0000147 130 FORMAT(' WILKS LAMBDA(',I2,',',I2,',',I2,') : ',F8.4) 0000148 IF (DESCRL .GT. 0.05) THEN 0000149 PRINT 140,DESCRL 0000150 140 FORMAT(/' DESCRIPTIVE LEVEL: ',F7.4,' ===> CANNOT REJECT NULL') 0000151 ELSE 0000152 PRINT 150,DESCRL 0000153 150 FORMAT(/' DESCRIPTIVE LEVEL: ',F7.4,' ===> REJECT NULL') 0000154 ENDIF 0000155 PRINT*,' ' 0000156 PRINT*,'========================================================='0000157 PRINT*,' ' 0000158 PRINT*,' ' 0000159 157 CONTINUE 0000160 C OFFER USER THE CHANCE TO TEST HYPOTHESES OF THE FORM H:CB=0 0000161 PRINT*,'DO YOU WANT TO TEST A RESTRICTION OF THE FORM' 0000162 PRINT*,'H: CB=0 ? (1=YES)' 0000163 READ*,ANSWER 0000164 IF (ANSWER .NE. 1) GOTO 158 0000165 PRINT*,'IS THE TEST ON THE ALTERNATIVE (1) OR THE NULL (0)?' 0000166 READ*,MODEL 0000167 IF (MODEL .NE. 0) MODEL=1 0000168 QT=MODEL*PA + (1-MODEL)*PO 0000169 PRINT*,'HOW MANY RESTRICTIONS?' 0000170 READ*,G 0000171 PRINT*,'INPUT C, ',G,' X ',QT 0000172 READ*,((C(I,J),J=1,QT),I=1,G) 0000173 IF (MODEL .EQ. 1) THEN 0000174 CALL RTEST(C,BETAA,Y,XXA,XA,PA,PY,N,G,IA,IB,LAMBDA,DESCRL,WORK) 0000175 ELSE 0000176 CALL RTEST(C,BETAO,Y,XXO,XO,PO,PY,N,G,IA,IB,LAMBDA,DESCRL,WORK) 0000177 ENDIF 0000178 PRINT*,' ' 0000179 PRINT*,' ' 0000180 PRINT 130,P,N-QT,G,LAMBDA 0000181 IF (DESCRL .GT. 0.05) THEN 0000182 PRINT 140,DESCRL 0000183 ELSE 0000184 PRINT 150,DESCRL 0000185 ENDIF 0000186 PRINT*,' ' 0000187 PRINT*,' ' 0000188 GOTO 157 0000189 158 CONTINUE 0000190 C OFFER USER THE CHANCE TO DELETE MORE COLUMNS 0000191 PRINT*,'DO YOU WISH TO MAKE ANOTHER RUN (1=YES) ?' 0000192 READ*,ANSWER 0000193 IF (ANSWER .EQ. 0) GOTO 190 0000194 C THE USER MAY SWITCH HYPOTHESES 0000195 PRINT*,'DO YOU WISH TO KEEP THE ALTERNATIVE HYPOTHESIS (1) OR' 0000196 PRINT*,'SUBSTITUTE THE NULL HYPOTHESIS (0)?' 0000197 READ*,ANSWER 0000198 IF (ANSWER .EQ. 0) THEN 0000199 DO 180 I = 1,N 0000200 DO 170 J = 1,PO 0000201 XA(I,J)=XO(I,J) 0000202 IF (I.LT.PO) THEN 0000203 XXA(I,J)=XXO(I,J) 0000204 IF (I.LT.PY) THEN 0000205 BETAA(J,I)=BETAO(J,I) 0000206 IF(J.LT.PY)THEN 0000207 SIGMAA(I,J)=SIGMAO(I,J) 0000208 ENDIF 0000209 ENDIF 0000210 ENDIF 0000211 170 CONTINUE 0000212 180 CONTINUE 0000213 PA=PO 0000214 ENDIF 0000215 GOTO 20 0000216 190 CONTINUE 0000217 END 0000218 ************************************************************************0000219 SUBROUTINE GTOPFN(UN,FN,NEW,OLD,TITLE) 0000220 C GET AND OPEN FILENAME 0000221 C REVISED VERSION FEB 2 1985 0000222 INTEGER UN 0000223 LOGICAL NEW,OLD 0000224 CHARACTER*60 FN 0000225 CHARACTER*80 TITLE 0000226 C UN (INPUT) UNIT TO OPEN FILE ON 0000227 C FN (OUTPUT) NAME OF FILE THAT HAS BEEN OPENED 0000228 C NEW (INPUT) 0000229 C .TRUE. IF THE FILE IS REQUIRED TO BE A NEW FILE 0000230 C .FALSE. OTHERWISE 0000231 C OLD (INPUT) 0000232 C .TRUE. IF THE FILE IS REQUIRED TO BE A OLD FILE 0000233 C .FALSE. OTHERWISE 0000234 C TITLE (INPUT) PROMPT TO BE PRINTED WHEN ASKING FOR FN 0000235 C LOCAL VARIABLES 0000236 LOGICAL EX 0000237 10 CONTINUE 0000238 PRINT* 0000239 PRINT 710,TITLE 0000240 710 FORMAT(A) 0000241 READ 710,FN 0000242 INQUIRE(FILE=FN,EXIST=EX) 0000243 IF(NEW.AND.EX) THEN 0000244 PRINT*,'FILE = ' 0000245 PRINT*,' ',FN 0000246 PRINT*,'EXISTS. PLEASE TRY AGAIN.' 0000247 GOTO 10 0000248 ENDIF 0000249 IF(OLD.AND.(.NOT.EX)) THEN 0000250 PRINT*,' FILE = ' 0000251 PRINT*,' ',FN 0000252 PRINT*,' DOES NOT EXIST. PLEASE TRY AGAIN.' 0000253 GOTO 10 0000254 ENDIF 0000255 OPEN(UNIT=UN,FILE=FN,ERR=20) 0000256 REWIND UN 0000257 RETURN 0000258 20 CONTINUE 0000259 PRINT*,' ERROR IN OPENING FILE =' 0000260 PRINT*,FN 0000261 PRINT*,' PLEASE TRY AGAIN.' 0000262 GOTO 10 0000263 END 0000264 ************************************************************************0000265 SUBROUTINE BPRINT(B,S,N,M,N1,IB,XX,XPX,CB) 0000266 C THIS SUBROUTINE PRINTS OUT BETA AND SIGMA MATRICES. 0000267 INTEGER IB,I,N,J,M,K,L,N1 0000268 REAL B(IB,IB), S(IB,IB), XX(IB,IB) 0000269 LOGICAL XPX,CB 0000270 PRINT*,' ' 0000271 PRINT*,' ' 0000272 PRINT*,'BETA ESTIMATES:' 0000273 PRINT*,' ' 0000274 DO 20 I=1,N 0000275 PRINT 10,(B(I,J),J=1,M) 0000276 10 FORMAT(15(1X,F8.3)) 0000277 PRINT*,' ' 0000278 20 CONTINUE 0000279 PRINT*,' ' 0000280 PRINT*,' ' 0000281 PRINT*,'SIGMA ESTIMATES (UNBIASED--N = ',N1,', P = ',N,' ):' 0000282 PRINT*,' ' 0000283 DO 30 I=1,M 0000284 PRINT 10, ( (S(I,J)*N1/(N1-N) ),J=1,M) 0000285 PRINT*,' ' 0000286 30 CONTINUE 0000287 PRINT*,' ' 0000288 PRINT*,' ' 0000289 IF (XPX) THEN 0000290 PRINT*,'INV(X"X):' 0000291 PRINT*,' ' 0000292 DO 40 I=1,N 0000293 PRINT 10,(XX(I,J),J=1,N) 0000294 PRINT*,' ' 0000295 40 CONTINUE 0000296 PRINT*,' ' 0000297 PRINT*,' ' 0000298 ENDIF 0000299 IF (CB) THEN 0000300 PRINT*,'COVARIANCE BETA:' 0000301 PRINT*,' ' 0000302 PRINT 50,(((10*I+J),J=1,M),I=1,N) 0000303 50 FORMAT (15(:4X,'B',I2,2X)) 0000304 PRINT*,' ' 0000305 DO 70 I=1,N 0000306 DO 60 J=1,M 0000307 PRINT 10, (((( S(J,L)*N1/(N1-N) )*XX(I,K)),L=1,M),K=1,N) 0000308 60 CONTINUE 0000309 70 CONTINUE 0000310 PRINT*,' ' 0000311 PRINT*,' ' 0000312 ENDIF 0000313 RETURN 0000314 END 0000315 ************************************************************************0000316 SUBROUTINE BETSIG(Y,X,N,PY,P,BETA,SIGMA,IA,IB,W,E) 0000317 C THIS SUBROUTINE CALCULATES THE BETA AND SIGMA ESTIMATES. 0000318 INTEGER IA,IB,I,J,K,N,IJOB,IER 0000319 REAL Y(IA,IB),X(IA,IB),BETA(IB,IB),SIGMA(IB,IB),W(IB,IB),E(IA,IB) 0000320 INTEGER IRB 0000321 PARAMETER(IRB=10) 0000322 REAL WKAREA(2*IRB),DUMMYB,Z(IRB,IRB) 0000323 REAL D1,D2 0000324 INTEGER PY,P 0000325 EXTERNAL LINV3F 0000326 C W = X'X 0000327 DO 3 I=1,P 0000328 DO 2 J=1,P 0000329 W(I,J)=0.0 0000330 DO 1 K=1,N 0000331 W(I,J)=W(I,J)+X(K,I)*X(K,J) 0000332 1 CONTINUE 0000333 2 CONTINUE 0000334 3 CONTINUE 0000335 C W = INV(X'X) 0000336 IJOB=1 0000337 D1=-1.0 0000338 CALL LINV3F(W,DUMMYB,IJOB,P,IB,D1,D2,WKAREA,IER) 0000339 C Z = X'Y 0000340 DO 13 I=1,P 0000341 DO 12 J=1,PY 0000342 Z(I,J)=0.0 0000343 DO 11 K=1,N 0000344 Z(I,J)=Z(I,J)+X(K,I)*Y(K,J) 0000345 11 CONTINUE 0000346 12 CONTINUE 0000347 13 CONTINUE 0000348 C BETA = INV(X'X)X'Y 0000349 DO 23 I=1,P 0000350 DO 22 J=1,PY 0000351 BETA(I,J)=0.0 0000352 DO 21 K=1,P 0000353 BETA(I,J)=BETA(I,J)+W(I,K)*Z(K,J) 0000354 21 CONTINUE 0000355 22 CONTINUE 0000356 23 CONTINUE 0000357 C CALCULATE E 0000358 DO 33 I=1,N 0000359 DO 32 J=1,PY 0000360 E(I,J)=Y(I,J) 0000361 DO 31 K=1,P 0000362 E(I,J)=E(I,J)-X(I,K)*BETA(K,J) 0000363 31 CONTINUE 0000364 32 CONTINUE 0000365 33 CONTINUE 0000366 C CALCULATE SIGMA=E'E/N 0000367 DO 43 I=1,PY 0000368 DO 42 J=1,PY 0000369 SIGMA(I,J)=0.0 0000370 DO 41 K=1,N 0000371 SIGMA(I,J)=SIGMA(I,J)+E(K,I)*E(K,J)/N 0000372 41 CONTINUE 0000373 42 CONTINUE 0000374 43 CONTINUE 0000375 RETURN 0000376 END 0000377 ************************************************************************0000378 SUBROUTINE WILKS(SIGMAA,SIGMAO,P,M,Q,VY,V1,LAMBDA,DESCRL) 0000379 C THIS SUBROUTINE CALCULATES WILKS LAMBDA AND THE ASSOCIATED 0000380 C DESCRIPTIVE LEVEL. 0000381 INTEGER IER,VY,V1,P,M,Q,N1,N2,IJOB 0000382 REAL SIGMAA(V1,V1),SIGMAO(V1,V1),DETO,DETA,D1,D2,LAMBDA 0000383 REAL L,CS,PR,LN,DESCRL,F,DUMMYB,S,DF 0000384 INTEGER IRB 0000385 PARAMETER(IRB=10) 0000386 REAL WKAREA(IRB) 0000387 EXTERNAL LINV3F,MDCH,MDFD 0000388 INTRINSIC ABS,LOG,REAL,SQRT 0000389 C CALCULATE THE DETERMINANTS OF SIGMAA AND SIGMAO (DETA, DETO) 0000390 D1=99.0 0000391 IJOB=4 0000392 CALL LINV3F(SIGMAA,DUMMYB,IJOB,VY,V1,D1,D2,WKAREA,IER) 0000393 DETA=D1*2**D2 0000394 PRINT*,' ' 0000395 PRINT*,'DETERMINANT FOR FULL MODEL: ',DETA 0000396 D1=99.0 0000397 CALL LINV3F(SIGMAO,DUMMYB,IJOB,VY,V1,D1,D2,WKAREA,IER) 0000398 DETO=D1*2**D2 0000399 PRINT*,' ' 0000400 PRINT*,'DETERMINANT FOR NULL MODEL: ',DETO 0000401 PRINT*,' ' 0000402 C CALCULATE WILKS LAMBDA 0000403 IF (ABS(DETO).LT.1.E-10) THEN 0000404 PRINT*,' ' 0000405 PRINT*,'ZERO DETERMINANT FOR NULL HYPOTHESIS.' 0000406 PRINT*,' ' 0000407 STOP 0000408 ENDIF 0000409 LAMBDA=DETA/DETO 0000410 C CALCULATE THE DESCRIPTIVE LEVEL 0000411 L=(1-LAMBDA)/LAMBDA 0000412 IF (P.EQ.1) GOTO 10 0000413 IF (P.EQ.2) GOTO 20 0000414 IF (Q.EQ.1) GOTO 30 0000415 IF (Q.EQ.2) GOTO 40 0000416 LN=LOG(LAMBDA) 0000417 CS=-(REAL(M)-REAL(P-Q+1)/2.0)*LN 0000418 DF=REAL(Q*P) 0000419 CALL MDCH(CS,DF,PR,IER) 0000420 DESCRL=1-PR 0000421 GOTO 60 0000422 10 F=L*(REAL(M)/REAL(Q)) 0000423 N1=Q 0000424 N2=M 0000425 GOTO 50 0000426 20 S=SQRT(LAMBDA) 0000427 F=((1.0-S)/S)*(REAL(M-1)/REAL(Q)) 0000428 N1=2*Q 0000429 N2=2*(M-1) 0000430 GOTO 50 0000431 30 F=L*REAL(M-P+1)/REAL(P) 0000432 N1=P 0000433 N2=M-P+1 0000434 GOTO 50 0000435 40 S=SQRT(LAMBDA) 0000436 F=((1.0-S)/S)*((REAL(M-P+1))/REAL(P)) 0000437 N1=2*P 0000438 N2=2*(M-P+1) 0000439 50 CONTINUE 0000440 CALL MDFD(F,N1,N2,PR,IER) 0000441 DESCRL=1-PR 0000442 60 RETURN 0000443 END 0000444 ************************************************************************0000445 SUBROUTINE RTEST(C,B,Y,XXI,X,Q,P,N,G,IA,IB,LAMBDA,DESCRL,T1) 0000446 C THIS SUBROUTINE CALCULATES LAMBDA AND DESCRL FOR A TEST 0000447 C OF THE FORM H: CB=0 0000448 INTEGER I,J,K,P,N,Q,G,IA,IB,IJOB,IER,M,IRB,K1,K2 0000449 REAL C(IB,IB),B(IB,IB),Y(IA,IB),XXI(IB,IB),X(IA,IB),T1(IB,IB) 0000450 REAL LAMBDA,DESCRL 0000451 REAL D1,D2 0000452 PARAMETER (IRB=10) 0000453 REAL E(IRB,IRB),HAE(IRB,IRB),T3(IRB,IRB),T4(IRB,IRB) 0000454 REAL WKAREA(2*IRB),CB(IRB,IRB),DUMMYB 0000455 EXTERNAL LINV3F,WILKS 0000456 INTRINSIC REAL 0000457 C CB=C*B 0000458 DO 3 I=1,G 0000459 DO 2 J=1,P 0000460 CB(I,J)=0.0 0000461 DO 1 K=1,Q 0000462 CB(I,J)=CB(I,J) + C(I,K)*B(K,J) 0000463 1 CONTINUE 0000464 2 CONTINUE 0000465 3 CONTINUE 0000466 C T1 = Y'*X 0000467 DO 13 I=1,P 0000468 DO 12 J=1,Q 0000469 T1(I,J)=0.0 0000470 DO 11 K=1,N 0000471 T1(I,J)=T1(I,J) + Y(K,I)*X(K,J) 0000472 11 CONTINUE 0000473 12 CONTINUE 0000474 13 CONTINUE 0000475 C T3=C*INV(X'X) 0000476 DO 17 I=1,G 0000477 DO 16 J=1,Q 0000478 T3(I,J)=0.0 0000479 DO 15 K=1,Q 0000480 T3(I,J)=T3(I,J)+C(I,K)*XXI(K,J) 0000481 15 CONTINUE 0000482 16 CONTINUE 0000483 17 CONTINUE 0000484 C T4 = T3*C' 0000485 DO 23 I=1,G 0000486 DO 22 J=1,G 0000487 T4(I,J)=0.0 0000488 DO 21 K=1,Q 0000489 T4(I,J)=T4(I,J)+T3(I,K)*C(J,K) 0000490 21 CONTINUE 0000491 22 CONTINUE 0000492 23 CONTINUE 0000493 C T4=INV(T4) 0000494 D1=-1.0 0000495 IJOB=1 0000496 CALL LINV3F(T4,DUMMYB,IJOB,G,IRB,D1,D2,WKAREA,IER) 0000497 C T3=(CB)'*T4 0000498 DO 33 I=1,P 0000499 DO 32 J=1,G 0000500 T3(I,J)=0.0 0000501 DO 31 K=1,G 0000502 T3(I,J)=T3(I,J)+CB(K,I)*T4(K,J) 0000503 31 CONTINUE 0000504 32 CONTINUE 0000505 33 CONTINUE 0000506 C E=Y'Y - [Y'*X]*[INV(X'*X)]*[X'*Y] HAE=E + T3*CB 0000507 DO 43 I=1,P 0000508 DO 42 J=1,P 0000509 E(I,J)=0.0 0000510 HAE(I,J)=0.0 0000511 DO 37 K1=1,N 0000512 E(I,J)=E(I,J)+Y(K1,I)*Y(K1,J) 0000513 37 CONTINUE 0000514 DO 38 K1=1,G 0000515 HAE(I,J)=HAE(I,J)+T3(I,K1)*CB(K1,J) 0000516 38 CONTINUE 0000517 DO 40 K1=1,Q 0000518 DO 39 K2=1,Q 0000519 E(I,J)=E(I,J)-T1(I,K1)*XXI(K1,K2)*T1(J,K2) 0000520 39 CONTINUE 0000521 40 CONTINUE 0000522 HAE(I,J)=HAE(I,J)+E(I,J) 0000523 E(I,J)=E(I,J)/REAL(N) 0000524 HAE(I,J)=HAE(I,J)/REAL(N) 0000525 42 CONTINUE 0000526 43 CONTINUE 0000527 C GET LAMBDA AND DESCRL 0000528 M=N-Q 0000529 CALL WILKS(E,HAE,P,M,G,P,IB,LAMBDA,DESCRL) 0000530 RETURN 0000531 END 0000532 ************************************************************************0000533 C SUBROUTINE LINV3F(A,B,IJOB,N,N,IA,D1,D2,WKAREA,IER) 0000534 C IMSL ROUTINE FOR COMPUTING DETERMINANTS 0000535 C TERMINAL ERROR: 0000536 C IER=130 INDICATES THAT MATRIX A IS ALGORITHMICALLY SINGULAR. 0000537 C END 0000538 ************************************************************************0000539 C SUBROUTINE MDCH(CS,DF,P,IER) 0000540 C IMSL ROUTINE TO DETERMINE P VALUE FOR CHI-SQUARED PROBABILITY 0000541 C DISTRIBUTION FUNCTION 0000542 C END 0000543 ************************************************************************0000544 C SUBROUTINE MDFH(F,N1,N2,P,IER) 0000545 C IMSL ROUTINE TO DETERMINE P VALUE FOR AN F 0000546 C DISTRIBUTION FUNCTION 0000547 C END 0000548 $$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ APPENDIX D: THE CALL STATEMENTS FROM GENLIN 1: PROGRAM GENLIN 39: CALL GTOPFN(UN,FN,NEW,OLD,TITLE) 125: CALL BETSIG(Y,XA,N,PY,PA,BETAA,SIGMAA,IA,IB,XXA,WORK) 127: CALL BETSIG(Y,XO,N,PY,PO,BETAO,SIGMAO,IA,IB,XXO,WORK) 133: CALL BPRINT(BETAO,SIGMAO,PO,PY,N,IB,XXO,XPXO,CBO) 138: CALL BPRINT(BETAA,SIGMAA,PA,PY,N,IB,XXA,XPXA,CBA) 144: CALL WILKS(SIGMAA,SIGMAO,P,M,Q,PY,IB,LAMBDA,DESCRL) 175: CALL RTEST(C,BETAA,Y,XXA,XA,PA,PY,N,G,IA,IB,LAMBDA,DESCRL,WORK) 177: CALL RTEST(C,BETAO,Y,XXO,XO,PO,PY,N,G,IA,IB,LAMBDA,DESCRL,WORK) 220: SUBROUTINE GTOPFN(UN,FN,NEW,OLD,TITLE) 266: SUBROUTINE BPRINT(B,S,N,M,N1,IB,XX,XPX,CB) 317: SUBROUTINE BETSIG(Y,X,N,PY,P,BETA,SIGMA,IA,IB,W,E) 339: CALL LINV3F(W,DUMMYB,IJOB,P,IB,D1,D2,WKAREA,IER) 379: SUBROUTINE WILKS(SIGMAA,SIGMAO,P,M,Q,VY,V1,LAMBDA,DESCRL) 393: CALL LINV3F(SIGMAA,DUMMYB,IJOB,VY,V1,D1,D2,WKAREA,IER) 398: CALL LINV3F(SIGMAO,DUMMYB,IJOB,VY,V1,D1,D2,WKAREA,IER) 420: CALL MDCH(CS,DF,PR,IER) 441: CALL MDFD(F,N1,N2,PR,IER) 446: SUBROUTINE RTEST(C,B,Y,XXI,X,Q,P,N,G,IA,IB,LAMBDA,DESCRL,T1) 497: CALL LINV3F(T4,DUMMYB,IJOB,G,IRB,D1,D2,WKAREA,IER) 530: CALL WILKS(E,HAE,P,M,G,P,IB,LAMBDA,DESCRL) THE FOLLOWING SUBROUTINE ARE FROM THE IMSL SUBROUTINE LIBRARY: LINV3F MDCH MDFD