Google
 

Trailing-Edge - PDP-10 Archives - decus_20tap2_198111 - decus/20-0026/polrg.smp
There are 2 other files named polrg.smp in the archive. Click here to see a list.
C                                                                       PLRG  10
C     ..................................................................PLRG  20
C                                                                       PLRG  30
C        SAMPLE MAIN PROGRAM FOR POLYNOMIAL REGRESSION - POLRG          PLRG  40
C                                                                       PLRG  50
C        PURPOSE                                                        PLRG  60
C           (1) READ THE PROBLEM PARAMETER CARD FOR A POLYNOMIAL REGRES-PLRG  70
C           SION, (2) CALL SUBROUTINES TO PERFORM THE ANALYSIS, (3)     PLRG  80
C           PRINT THE REGRESSION COEFFICIENTS AND ANALYSIS OF VARIANCE  PLRG  90
C           TABLE FOR POLYNOMIALS OF SUCCESSIVELY INCREASING DEGREES,   PLRG 100
C           AND (4) OPTIONALLY PRINT THE TABLE OF RESIDUALS AND A PLOT  PLRG 110
C           OF Y VALUES AND Y ESTIMATES.                                PLRG 120
C                                                                       PLRG 130
C        REMARKS                                                        PLRG 140
C           THE NUMBER OF OBSERVATIONS, N, MUST BE GREATER THAN M+1,    PLRG 150
C           WHERE M IS THE HIGHEST DEGREE POLYNOMIAL SPECIFIED.         PLRG 160
C           IF THERE IS NO REDUCTION IN THE RESIDUAL SUM OF SQUARES     PLRG 170
C           BETWEEN TWO SUCCESSIVE DEGREES OF THE POLYNOMIALS, THE      PLRG 180
C           PROGRAM TERMINATES THE PROBLEM BEFORE COMPLETING THE ANALY- PLRG 190
C           SIS FOR THE HIGHEST DEGREE POLYNOMIAL SPECIFIED.            PLRG 200
C                                                                       PLRG 210
C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED                  PLRG 220
C           GDATA                                                       PLRG 230
C           ORDER                                                       PLRG 240
C           MINV                                                        PLRG 250
C           MULTR                                                       PLRG 260
C           PLOT  (A SPECIAL PLOT SUBROUTINE PROVIDED FOR THE SAMPLE    PLRG 270
C                 PROGRAM.)                                             PLRG 280
C                                                                       PLRG 290
C        METHOD                                                         PLRG 300
C           REFER TO B. OSTLE, 'STATISTICS IN RESEARCH', THE IOWA STATE PLRG 310
C           COLLEGE PRESS', 1954, CHAPTER 6.                            PLRG 320
C                                                                       PLRG 330
C     ..................................................................PLRG 340
C                                                                       PLRG 350
C     THE FOLLOWING DIMENSION MUST BE GREATER THAN OR EQUAL TO THE      PLRG 360
C     PRODUCT OF N*(M+1), WHERE N IS THE NUMBER OF OBSERVATIONS AND M   PLRG 370
C     IS THE HIGHEST DEGREE POLYNOMIAL SPECIFIED..                      PLRG 380
C                                                                       PLRG 390
         DIMENSION X(1100)                                              PLRG 400
C                                                                       PLRG 410
C     THE FOLLOWING DIMENSION MUST BE GREATER THAN OR EQUAL TO THE      PLRG 420
C     PRODUCT OF M*M..                                                  PLRG 430
C                                                                       PLRG 440
         DIMENSION DI(100)                                              PLRG 450
C                                                                       PLRG 460
C     THE FOLLOWING DIMENSION MUST BE GREATER THAN OR EQUAL TO          PLRG 470
C     (M+2)*(M+1)/2..                                                   PLRG 480
C                                                                       PLRG 490
         DIMENSION D(66)                                                PLRG 500
C                                                                       PLRG 510
C     THE FOLLOWING DIMENSIONS MUST BE GREATER THAN OR EQUAL TO M..     PLRG 520
C                                                                       PLRG 530
         DIMENSION B(10),E(10),SB(10),T(10)                             PLRG 540
C                                                                       PLRG 550
C     THE FOLLOWING DIMENSIONS MUST BE GREATER THAN OR EQUAL TO (M+1).. PLRG 560
C                                                                       PLRG 570
         DIMENSION XBAR(11),STD(11),COE(11),SUMSQ(11),ISAVE(11)         PLRG 580
C                                                                       PLRG 590
C     THE FOLLOWING DIMENSION MUST BE GREATER THAN OR EQUAL TO 10..     PLRG 600
C                                                                       PLRG 610
         DIMENSION ANS(10)                                              PLRG 620
C                                                                       PLRG 630
C     THE FOLLOWING DIMENSION WILL BE USED IF THE PLOT OF OBSERVED DATA PLRG 640
C     AND ESTIMATES IS DESIRED.  THE SIZE OF THE DIMENSION, IN THIS     PLRG 650
C     CASE, MUST BE GREATER THAN OR EQUAL TO N*3.  OTHERWISE, THE SIZE  PLRG 660
C     OF DIMENSION MAY BE SET TO 1.                                     PLRG 670
C                                                                       PLRG 680
         DIMENSION P(300)                                               PLRG 690
C                                                                       PLRG 700
C     ..................................................................PLRG 710
C                                                                       PLRG 720
C        IF A DOUBLE PRECISION VERSION OF THIS ROUTINE IS DESIRED, THE  PLRG 730
C        C IN COLUMN 1 SHOULD BE REMOVED FROM THE DOUBLE PRECISION      PLRG 740
C        STATEMENT WHICH FOLLOWS.                                       PLRG 750
C                                                                       PLRG 760
C     DOUBLE PRECISION X,XBAR,STD,D,SUMSQ,DI,E,B,SB,T,ANS,DET,COE       PLRG 770
C                                                                       PLRG 780
C        THE C MUST ALSO BE REMOVED FROM DOUBLE PRECISION STATEMENTS    PLRG 790
C        APPEARING IN OTHER ROUTINES USED IN CONJUNCTION WITH THIS      PLRG 800
C        ROUTINE.                                                       PLRG 810
C                                                                       PLRG 820
C        ...............................................................PLRG 830
C                                                                       PLRG 840
    1 FORMAT(A4,A2,I5,I2,I1)                                            PLRG 850
    2 FORMAT(2F6.0)                                                     PLRG 860
    3 FORMAT(27H1POLYNOMIAL REGRESSION.....,A4,A2/)                     PLRG 870
    4 FORMAT(23H0NUMBER OF OBSERVATIONS,I6//)                           PLRG 880
    5 FORMAT(32H0POLYNOMIAL REGRESSION OF DEGREE,I3)                    PLRG 890
    6 FORMAT(12H0  INTERCEPT,E20.7)                                     PLRG 900
    7 FORMAT(26H0  REGRESSION COEFFICIENTS/(6E20.7))                    PLRG 910
    8 FORMAT(1H0/24X,24HANALYSIS OF VARIANCE FOR,I4,19H  DEGREE POLYNOMIPLRG 920
     1AL/)                                                              PLRG 930
    9 FORMAT(1H0,5X,19HSOURCE OF VARIATION,7X,9HDEGREE OF,7X,6HSUM OF,9XPLRG 940
     1,4HMEAN,10X,1HF,9X,20HIMPROVEMENT IN TERMS/33X,7HFREEDOM,8X,7HSQUAPLRG 950
     2RES,7X,6HSQUARE,7X,5HVALUE,8X,17HOF SUM OF SQUARES)               PLRG 960
   10 FORMAT(20H0  DUE TO REGRESSION,12X,I6,F17.5,F14.5,F13.5,F20.5)    PLRG 970
   11 FORMAT(32H   DEVIATION ABOUT REGRESSION   ,I6,F17.5,F14.5)        PLRG 980
   12 FORMAT(8X,5HTOTAL,19X,I6,F17.5///)                                PLRG 990
   13 FORMAT(17H0  NO IMPROVEMENT)                                      PLRG1000
   14 FORMAT(1H0//27X,18HTABLE OF RESIDUALS//16H OBSERVATION NO.,5X,7HX PLRG1010
     1VALUE,7X,7HY VALUE,7X,10HY ESTIMATE,7X,8HRESIDUAL/)               PLRG1020
   15 FORMAT(1H0,3X,I6,F18.5,F14.5,F17.5,F15.5)                         PLRG1030
C                                                                       PLRG1040
C     ..................................................................PLRG1050
C                                                                       PLRG1060
C     READ PROBLEM PARAMETER CARD                                       PLRG1070
C                                                                       PLRG1080
  100 READ (5,1,END=999) PR,PR1,N,M,NPLOT                               PLRG1090
C                                                                       PLRG1100
C        PR....PROBLEM NUMBER (MAY BE ALPHAMERIC)                       PLRG1110
C        PR1...PROBLEM NUMBER (CONTINUED)                               PLRG1120
C        N.....NUMBER OF OBSERVATIONS                                   PLRG1130
C        M.....HIGHEST DEGREE POLYNOMIAL SPECIFIED                      PLRG1140
C        NPLOT.OPTION CODE FOR PLOTTING                                 PLRG1150
C              0  IF PLOT IS NOT DESIRED.                               PLRG1160
C              1  IF PLOT IS DESIRED.                                   PLRG1170
C                                                                       PLRG1180
C     PRINT PROBLEM NUMBER AND N.                                       PLRG1190
C                                                                       PLRG1200
      WRITE (6,3) PR,PR1                                                PLRG1210
      WRITE (6,4) N                                                     PLRG1220
C                                                                       PLRG1230
C     READ INPUT DATA                                                   PLRG1240
C                                                                       PLRG1250
      L=N*M                                                             PLRG1260
      DO 110 I=1,N                                                      PLRG1270
      J=L+I                                                             PLRG1280
C                                                                       PLRG1290
C        X(I) IS THE INDEPENDENT VARIABLE, AND X(J) IS THE DEPENDENT    PLRG1300
C        VARIABLE.                                                      PLRG1310
C                                                                       PLRG1320
  110 READ (5,2) X(I),X(J)                                              PLRG1330
C                                                                       PLRG1340
      CALL GDATA (N,M,X,XBAR,STD,D,SUMSQ)                               PLRG1350
C                                                                       PLRG1360
      MM=M+1                                                            PLRG1370
      SUM=0.0                                                           PLRG1380
      NT=N-1                                                            PLRG1390
C                                                                       PLRG1400
      DO 200 I=1,M                                                      PLRG1410
      ISAVE(I)=I                                                        PLRG1420
C                                                                       PLRG1430
C     FORM SUBSET OF CORRELATION COEFFICIENT MATRIX                     PLRG1440
C                                                                       PLRG1450
      CALL ORDER (MM,D,MM,I,ISAVE,DI,E)                                 PLRG1460
C                                                                       PLRG1470
C     INVERT THE SUBMATRIX OF CORRELATION COEFFICIENTS                  PLRG1480
C                                                                       PLRG1490
      CALL MINV (DI,I,DET,B,T)                                          PLRG1500
C                                                                       PLRG1510
      CALL MULTR (N,I,XBAR,STD,SUMSQ,DI,E,ISAVE,B,SB,T,ANS)             PLRG1520
C                                                                       PLRG1530
C     PRINT THE RESULT OF CALCULATION                                   PLRG1540
C                                                                       PLRG1550
      WRITE (6,5) I                                                     PLRG1560
      IF(ANS(7)) 140,130,130                                            PLRG1570
  130 SUMIP=ANS(4)-SUM                                                  PLRG1580
      IF(SUMIP) 140, 140, 150                                           PLRG1590
  140 WRITE (6,13)                                                      PLRG1600
      GO TO 210                                                         PLRG1610
  150 WRITE (6,6) ANS(1)                                                PLRG1620
      WRITE (6,7) (B(J),J=1,I)                                          PLRG1630
      WRITE (6,8) I                                                     PLRG1640
      WRITE (6,9)                                                       PLRG1650
      SUM=ANS(4)                                                        PLRG1660
      WRITE (6,10) I,ANS(4),ANS(6),ANS(10),SUMIP                        PLRG1670
      NI=ANS(8)                                                         PLRG1680
      WRITE (6,11) NI,ANS(7),ANS(9)                                     PLRG1690
      WRITE (6,12) NT,SUMSQ(MM)                                         PLRG1700
C                                                                       PLRG1710
C     SAVE COEFFICIENTS FOR CALCULATION OF Y ESTIMATES                  PLRG1720
C                                                                       PLRG1730
      COE(1)=ANS(1)                                                     PLRG1740
      DO 160 J=1,I                                                      PLRG1750
  160 COE(J+1)=B(J)                                                     PLRG1760
      LA=I                                                              PLRG1770
  200 CONTINUE                                                          PLRG1780
C                                                                       PLRG1790
C     TEST WHETHER PLOT IS DESIRED                                      PLRG1800
C                                                                       PLRG1810
  210 IF(NPLOT) 100, 100, 220                                           PLRG1820
C                                                                       PLRG1830
C        CALCULATE ESTIMATES                                            PLRG1840
C                                                                       PLRG1850
  220 NP3=N+N                                                           PLRG1860
      DO 230 I=1,N                                                      PLRG1870
      NP3=NP3+1                                                         PLRG1880
      P(NP3)=COE(1)                                                     PLRG1890
      L=I                                                               PLRG1900
      DO 230 J=1,LA                                                     PLRG1910
      P(NP3)=P(NP3)+X(L)*COE(J+1)                                       PLRG1920
  230 L=L+N                                                             PLRG1930
C                                                                       PLRG1940
C        COPY OBSERVED DATA                                             PLRG1950
C                                                                       PLRG1960
      N2=N                                                              PLRG1970
      L=N*M                                                             PLRG1980
      DO 240 I=1,N                                                      PLRG1990
      P(I)=X(I)                                                         PLRG2000
      N2=N2+1                                                           PLRG2010
      L=L+1                                                             PLRG2020
  240 P(N2)=X(L)                                                        PLRG2030
C                                                                       PLRG2040
C     PRINT TABLE OF RESIDUALS                                          PLRG2050
C                                                                       PLRG2060
      WRITE (6,3) PR,PR1                                                PLRG2070
      WRITE (6,5) LA                                                    PLRG2080
      WRITE (6,14)                                                      PLRG2090
      NP2=N                                                             PLRG2100
      NP3=N+N                                                           PLRG2110
      DO 250 I=1,N                                                      PLRG2120
      NP2=NP2+1                                                         PLRG2130
      NP3=NP3+1                                                         PLRG2140
      RESID=P(NP2)-P(NP3)                                               PLRG2150
  250 WRITE (6,15) I,P(I),P(NP2),P(NP3),RESID                           PLRG2160
C                                                                       PLRG2170
      CALL PLOT (LA,P,N,3,0,1)                                          PLRG2180
C                                                                       PLRG2190
      GO TO 100                                                         PLRG2200
999	STOP
      END                                                               PLRG2210