Google
 

Trailing-Edge - PDP-10 Archives - decuslib20-04 - decus/20-0113/cmodl.bas
There are 2 other files named cmodl.bas in the archive. Click here to see a list.
00020REM **************************************************************
00030REM       CMODL   CMODL   CMODL    CMODL   CMODL   CMODL   CMODL
00040REM **************************************************************
00050REM
00060REM           NORMAL FIT OF LEAST SQUARE UTILITY
00070REM
00080REM   THIS MODULES ATTEMPTS TO FIT A NORMAL OGIVE FUNCTION
00090REM   TO THE UTILITIES PASSED FROM CMODK (ASSESSMENT OF UTILITY).
00100REM   IF ONLY A VERY POOR FIT IS POSSIBLE THEN THE MODULES PRINTS
00110REM   THAT IT IS NOT POSSIBLE TO FIT A NORMAL OGIVE.
00120REM
00130REM   IF A NORMAL OGIVE DOES MAKE A REASONABLE FIT THE THE LEAST
00140REM   SQUARES, NORMAL OGIVE AND CORRECTED NORMAL FITS ARE PRINTED
00150REM
00160REM *****************************************************************
00170DIM X(84),Y(84),Z(84),P(9),U(9),Q(9)
00180FILES RFILE1,RFILE2,RFILE3
00230RESTORE#1
00231  INPUT#  1,I1,I2,I3
00240SCRATCH#1
00241  PRINT #  1,20,I2,I3
00250RESTORE#2
00260FOR I=1 TO 84
00271  INPUT#  2,Y(I)
00280NEXT I
00290FOR I=1 TO 9
00301  INPUT#  2,P(I)
00310NEXT I
00320FOR I=1 TO 9
00331  INPUT#  2,U(I)
00340NEXT I
00350E8=.005
00360N0=7
00370N1=N0+1
00380N2=N0+2
00390N3=N0-1
00400GOTO 1100
00410REM--***** normal approximation--*****
00420F9=0
00430REM initial guesses
00440U0=0
00450FOR I=2 TO 9
00460U0=U0+(U(I)-U(I-1))*((P(I)+P(I-1))/2)
00470NEXT I
00480U1=0
00490FOR I=2 TO 9
00500U1=U1+(U(I)-U(I-1))*((P(I)+P(I-1))/2)^2
00510NEXT I
00520U1=SQR(U1-U0*U0)
00530U8=U1
00540U9=U0
00550REM--u0 is mu  u1 is sigma
00560FOR I=1 TO 9
00570X=((P(I)-U0)/U1)
00580Y3=X
00590GOSUB 1580
00600W(I)=P
00610NEXT I
00620F9=F9+1
00630FOR I=1 TO 9
00640Q0=(P(I)-U0)/U1
00650Q(I)=.398942*EXP(-.5*Q0*Q0)
00660NEXT I
00670S0=0
00680G3=0
00690G2=0
00700FOR J=2 TO N1
00710J0=J-1
00720J1=J+1
00730FOR I=1 TO J0
00740FOR K=J1 TO N2
00750S0=S0+1
00760X(S0)=-1
00770Z(S0)=-1
00780IF Y(S0)<0 THEN 870
00790X(S0)=(Q(J)-Q(K))/(W(K)-W(J))-(Q(I)-Q(J))/(W(J)-W(I))
00800X(S0)=X(S0)/U1
00810Z(S0)=((P(J)-U0)*Q(J)-(P(K)-U0)*Q(K))/(W(K)-W(J))
00820Z(S0)=Z(S0)-((P(I)-U0)*Q(I)-(P(J)-U0)*Q(J))/(W(J)-W(I))
00830G0=(LOG(Y(S0)/(1-Y(S0)))-LOG((W(J)-W(I))/(W(K)-W(J))))
00840Z(S0)=Z(S0)/U1/U1
00850G2=G2+G0*X(S0)
00860G3=G3+G0*Z(S0)
00870NEXT K
00880NEXT I
00890NEXT J
00900REM--g2 is gmu   g3 is gsigma
00910A1=0
00920A2=0
00930A4=0
00940FOR I=1 TO 84
00950IF Y(I)<0 THEN 990
00960A1=A1+X(I)*X(I)
00970A2=A2+X(I)*Z(I)
00980A4=A4+Z(I)*Z(I)
00990NEXT I
01000A3=A2
01010X2=(G2/A1-G3/A3)/(A2/A1-A4/A3)
01020X1=G2/A1-A2/A1*X2
01030U0=U0-X1
01040U1=U1-X2
01050IF U0>50 THEN 1450
01060IF U1>50 THEN 1450
01070IF ABS(X1)>E8 THEN 560
01080IF ABS(X2)>E8 THEN 560
01090GOTO 1160
01100PRINT L$
01110PRINT
01120PRINT "                     UTILITIES"
01130PRINT "-------------------------------------------------------------"
01140PRINT " OUTCOME        LEAST SQUARES   CORRECTED NORMAL    NORMAL"
01150GOTO 410
01160:######.##         ##.##            ##.##          ##.##
01170Y3=(P(1)-U0)/U1
01180GOSUB 1580
01190U7=P
01200Y3=(P(9)-U0)/U1
01210GOSUB 1580
01220U6=P
01230FOR I=1 TO N2
01240Y3=(P(I)-U0)/U1
01250GOSUB 1580
01260PRINT  USING 1160,P(I),U(I),(P-U7)/(U6-U7),P
01270NEXT I
01280PRINT "--------------------------------------------------------------"
01290PRINT "THIS COMPLETES THE FITTING OF A NORMAL DISTRIBUTION."
01300PRINT "TYPE THE NUMBER OF THE OPTION YOU WANT."
01310PRINT "     1. FIND THE EXPECTED UTILITY USING LEAST SQUARES FIT"
01320PRINT "     2. FIND THE EXPECTED UTILITY USING NORMAL FIT (ASSUMES"
01330PRINT "        THE DISTRIBUTION IS ALSO NORMAL)"
01340PRINT "     3. EXIT MODULE"
01350GOSUB 1810
01360IF O1=3 THEN 1390
01370IF O1=1 THEN 1400
01380IF O1=2 THEN 1420
01390CHAIN "RSTRT"
01400SCRATCH#2
01402FOR I=1 TO 9
01403 PRINT #2,P(I)
01404NEXT I
01405FOR I=1 TO 9
01406PRINT#2,U(I)
01407NEXT I
01410CHAIN "CMODP"
01420SCRATCH#3
01421  PRINT #  3,U0,U1
01430SCRATCH#2
01432FOR I=1 TO 9
01433PRINT#2,P(I)
01434NEXT I
01435FOR I=1 TO 9
01436PRINT#2,U(9)
01437NEXT I
01440CHAIN "CMODS"
01450PRINT "***********************************************************"
01460PRINT "YOUR LEAST SQUARES UTILITIES CAN NOT BE FITTED BY THIS "
01470PRINT "NORMAL ROUTINE."
01480PRINT "***********************************************************"
01490PRINT
01500PRINT "TYPE THE NUMBER OF THE OPTION YOU WANT."
01510PRINT "   1. FIND EXPECTED UTILITY USING THE LSQ UTILITIES"
01520PRINT "   2. EXIT MODULE"
01530GOSUB 1810
01540IF O1=2 THEN 1390
01550IF O1=1 THEN 1400
01560PRINT "REENTER.  INPUT MUST BE 1 OR 2."
01570GOTO 1530
01580REM **********************************************************
01590REM      ROUTINE CALCULATES THE CDF FOR NORMAL DISTRIBUTION
01600REM               INPUT       Y3
01610REM               OUTPUT      P
01620REM
01630Y4=ABS(Y3)
01640X1=X
01650X=Y3
01660T=1/(1+.231642*Y4)
01670D=.398942*EXP(-X*X/2)
01680C1=1.33027
01690C2=1.82126
01700C3=1.78148
01710C4=.356564
01720C5=.319382
01730P=1-D*T*((((C1*T-C2)*T+C3)*T-C4)*T+C5)
01740IF X >= 0 THEN 1760
01750P=1-P
01760X=X1
01770RETURN
01780REM
01790REM        END OF NORMAL CDF ROUTINE
01800REM **********************************************************
01810REM--SUBROUTINE THAT DETERMINES IF RESTART HAS BEEN REQUESTED.
01820INPUT O1
01830IF O1=-9999 THEN 1850
01840RETURN
01850CHAIN "RSTRT"
01860REM*************END ROUTINE
09999END