Trailing-Edge
-
PDP-10 Archives
-
decuslib10-06
-
43,50422/cmodt.bas
There are 2 other files named cmodt.bas in the archive. Click here to see a list.
00020REM **************************************************************
00030REM CMODT CMODT CMODT CMODT CMDOT CMODT
00040REM **************************************************************
00050 DIM L(10)
00060 FILES RFILE1,RFILE2,RFILE3
00100RESTORE#1
00101 INPUT# 1,I1,I2,I3
00110SCRATCH#1
00111 PRINT # 1,29,I2,I3
00120X=0
00130V7=0
00140V6=0
00150V4=0
00160 DIM S(99)
00170 DIM Q(99)
00180GOSUB 5160
00190V8=0
00200PRINT L$
00210PRINT " COMPARISON OF TWO PROPORTIONS-INDEPENDENT BETAS"
00220PRINT
00230PRINT "THIS MODULE ALLOWS YOU TO COMPARE TWO INDEPENDENTLY BETA"
00240PRINT "DISTRIBUTED PROPORTIONS."
00250PRINT
00260PRINT "THE MODULE WILL COMPUTE AND PRINT THE PROBABILITY THAT THE "
00270PRINT "DIFFERENCE OF THE TWO PROPORTIONS IS GREATER THAN K WHICH"
00280PRINT "IS TO BE SPECIFIED BY YOU. THE MODULE ALWAYS USES AS THE "
00290PRINT "THE DIFFERENCE THE BETA DISTRIBUTION WITH THE LARGER MEAN"
00300PRINT "MINUS THE ONE WITH THE SMALLER MEAN."
00310REM**********************************************************************
00320PRINT "INPUT THE PARAMETERS OF THE DISTRIBUTION ON PI-ONE (PI-1)."
00330REM********************************************************************
00340V7=0
00350P=0
00360PRINT
00370PRINT "INPUT PARAMETER A.";
00380GOSUB 400
00390GOTO 460
00400GOSUB 9000
00410IF O1>2000 THEN 430
00420IF O1 >= 1.15 THEN 450
00430PRINT "REENTER. MINIMUM IS 1.15. MAXIMUM IS 2000."
00440GOTO 400
00450RETURN
00460A1=O1
00470PRINT "INPUT PARAMETER B.";
00480GOSUB 400
00490PRINT
00500B1=O1
00510A=A1
00520B=B1
00530GOSUB 5220
00540K7=.01
00550K8=.03
00560GOSUB 620
00570GOTO 790
00580K7=K8+.01
00590K8=K7+.02
00600GOSUB 620
00610RETURN
00620IF K7<1 THEN 640
00630RETURN
00640IF K8 <= .995 THEN 660
00650K8=.995
00660IF V4=1 THEN 780
00670FOR P1=K7 TO K8+.001 STEP .01
00680K5=INT(P1*100+.55)
00690IF P>.9999 THEN 730
00700GOSUB 2500
00710IF P<1.E-17 THEN 750
00720IF P<1 THEN 760
00730P=1
00740GOTO 760
00750P=0
00760Q(K5)=P
00770NEXT P1
00780RETURN
00790REM*********************************************************************
00800PRINT "INPUT THE PARAMETERS OF THE DISTRIBUTION ON PI-TWO (PI-2)."
00810REM*********************************************************************
00820GOSUB 580
00830PRINT
00840PRINT "INPUT PARAMETER A.";
00850GOSUB 400
00860A2=O1
00870GOSUB 580
00880PRINT "INPUT PARAMETER B.";
00890GOSUB 400
00900B2=O1
00910GOSUB 580
00920O7=6
00930PRINT L$
00940M1=A1/(A1+B1)
00950M2=A2/(A2+B2)
00960IF V6=1 THEN 1390
00970V6=1
00980IF M1>M2 THEN 1000
00990GOTO 1070
01000PRINT "THE MODULE WILL COMPUTE AND PRINT THE PROBABILITIES FOR "
01010GOSUB 580
01020PRINT "THE DIFFERENCE P1-1 MINUS PI-2."
01030GOSUB 580
01040GOSUB 580
01050V7=1
01060GOTO 1120
01070PRINT "THE MODULE WILL COMPUTE AND PRINT THE PROBABILITIES FOR "
01080GOSUB 580
01090PRINT "THE DIFFERENCE P1-2 MINUS PI-1."
01100GOSUB 580
01110PRINT
01120PRINT
01130GOSUB 1240
01140PRINT
01150IF V7=0 THEN 1200
01160PRINT "IF, FOR EXAMPLE, YOU WANTED THE PROBABILITY THAT P1-1 IS"
01170GOSUB 580
01180PRINT "GREATER THAN PI-2 YOU WOULD SPECIFY A VALUE OF 0 FOR K."
01190GOTO 1290
01200PRINT "IF, FOR EXAMPLE, YOU WANTED THE PROBABILITY THAT PI-2 IS"
01210GOSUB 580
01220PRINT "GREATER THAN PI-ONE YOU WOULD SPECIFY A VALUE OF 0 FOR K."
01230GOTO 1290
01240PRINT "YOU CAN SPECIFY UP TO 5 K VALUES AT A TIME. THE MODULE WILL"
01250PRINT "COMPUTE AND PRINT THE PROBABILITIES FOR THESE VALUES AND THEN"
01260PRINT "ALLOW YOU TO SPECIFY MORE VALUES IF YOU WANT."
01270GOSUB 580
01280RETURN
01290PRINT
01300GOSUB 580
01310GOSUB 580
01320GOSUB 580
01330PRINT "WHEN YOU ARE READY TO CONTINUE TYPE '1'.";
01340GOSUB 9000
01350GOSUB 580
01360GOSUB 580
01370PRINT L$
01380GOTO 1460
01390PRINT "THE PROBABILITIES WILL BE COMPUTED FOR THE DIFFERENCE"
01400IF M1>M2 THEN 1430
01410PRINT " PI-2 MINUS PI-1."
01420GOTO 1450
01430PRINT " PI-1 MINUS PI-2."
01440V7=1
01450PRINT
01460PRINT "INPUT THE NUMBER OF K VALUES YOU WANT TO SPECIFY.";
01470GOSUB 9000
01480IF O1 >= 1 THEN 1510
01490PRINT "REENTER. N MUST BE AT LEAST 1 BUT NOT MORE THAN 5."
01500GOTO 1470
01510IF O1>5 THEN 1490
01520GOSUB 580
01530N=O1
01540FOR K0=1 TO N
01550REM
01560 PRINT"INPUT VALUE ";K0;" ";
01570GOSUB 9000
01580GOSUB 580
01590IF O1 <= -1 THEN 1740
01600IF O1 >= 1 THEN 1740
01610L(K0)=O1
01620NEXT K0
01630GOSUB 580
01640GOSUB 580
01650GOSUB 2580
01670IF V7=1 THEN 1700
01680K9=1
01690GOTO 1710
01700K9=-1
01710REM
01720V4=1
01730GOTO 1760
01740PRINT "REENTER. K MUST BE BETWEEN -1 AND 1 (EXCLUSIVE)."
01750GOTO 1570
01760FOR K5=1 TO N
01770K=L(K5)
01780A0=0
01790Q0=0
01800FOR P1=.01 TO .995 STEP .01
01810P2=P1+K9*K
01820F4=INT(P2*100)
01830IF F4>99 THEN 1870
01840:P2##.## F4######.##
01850IF P2<.01 THEN 1980
01860IF P2<1 THEN 1890
01870A0=A0+1-Q0
01880GOTO 2000
01890A0=A0+(Q(INT(100*P1+.55))-Q0)*S(F4)
01900IF F4+1 <= 99 THEN 1930
01910A7=1-S(F4)
01920GOTO 1940
01930A7=S(F4+1)-S(F4)
01940A7=A7*(P2*100-INT(P2*100))
01950A7=A7*Q(INT(100*P1+.55))-A7*Q0
01960A0=A0+A7
01970GOTO 1980
01980Q0=Q(INT(P1*100+.55))
01990NEXT P1
02000J2=.995+K*K9
02010IF J2<.01 THEN 2060
02020IF J2>.9999 THEN 2060
02030J1=0
02040GOSUB 5000
02050A0=(1-Q0)*P+A0
02060IF K9=1 THEN 2100
02070PRINT USING 2090,K,A0
02080GOTO 2110
02090: PROB( DIFF >###.### )=##.##
02100PRINT USING 2090,K,1-A0
02110NEXT K5
02120GOTO 3130
02130K7=INT(100*A/(A+B))/100
02140P=1
02150K0=K7
02160K1=K0-.03
02170IF K1>.01 THEN 2190
02180K1=.01
02190FOR P2=K0 TO K1-.003 STEP -.01
02200IF V4=1 THEN 2380
02210IF P<.00001 THEN 2260
02220GOSUB 2420
02230IF P<1 THEN 2270
02240P=1
02250GOTO 2270
02260P=0
02270S(INT(P2*100+.55))=P
02280NEXT P2
02290RETURN
02300FOR P2=K7+.01 TO .995 STEP .01
02310IF V4=1 THEN 2370
02320IF P>.9999 THEN 2350
02330GOSUB 2420
02340GOTO 2360
02350P=1
02360S(INT(P2*100+.55))=P
02370NEXT P2
02380RETURN
02390REM******************************************************************
02400REM CDF FOR PI-2
02410REM******************************************************************
02420REM
02430:######.## STANDARD DEVIATION ######.##
02440J2=P2-.005
02450GOSUB 5000
02460RETURN
02470REM******************************************************************
02480REM CDF FOR PI-1
02490REM******************************************************************
02500J1=0
02510J2=P1
02520GOSUB 5000
02530RETURN
02540REM**********************************************************
02550REM**********************************************************
02560REM**********************************************************
02570REM************************************************************
02580PRINT L$
02590GOSUB 580
02600PRINT " PI-1 BETA DISTRIBUTIONS PI-2"
02610GOSUB 580
02620PRINT "-------------------------------------------------------"
02630GOSUB 580
02640:######.## PARAMETER A ######.##
02650PRINT USING 2640,A1,A2
02660GOSUB 580
02670:######.## PARAMETER B ######.##
02680PRINT USING 2670,B1,B2
02690GOSUB 580
02700:######.## MEAN ######.##
02710PRINT USING 2700,M1,M2
02720GOSUB 580
02730GOTO 2770
02740:######.## MODE ######.##
02750PRINT USING 2740,(A1-1)/(A1+B1-2),(A2-1)/(A2+B2-2)
02760GOSUB 580
02770S5=(A1*B1)/(((A1+B1)**2)*(A1+B1+1))
02780S6=(A2*B2)/(((A2+B2)**2)*(A2+B2+1))
02790: STANDARD DEVIATION ##.##
02800GOSUB 620
02810K7=K8+.01
02820K8=.995
02830GOSUB 620
02840PRINT "------------------------------------------------------"
02850IF V7=1 THEN 2970
02860PRINT " DIFFERENCE (PI-2 MINUS PI-1)"
02870A=A2
02880B=B2
02890GOSUB 5220
02900: MEAN ##.##
02910PRINT USING 2900,M2-M1
02920GOSUB 2130
02930PRINT USING 2790,SQR(S5+S6)
02940K0=K1-.01
02950GOSUB 2160
02960GOTO 3070
02970PRINT " DIFFERENCE (PI-1 MINUS PI-2)"
02980A=A2
02990B=B2
03000GOSUB 5220
03010: MEAN ##.##
03020PRINT USING 3010,M1-M2
03030GOSUB 2130
03040PRINT USING 2790,SQR(S5+S6)
03050K0=K1-.01
03060GOSUB 2160
03070PRINT "-------------------------------------------------------"
03080K0=K1-.01
03090GOSUB 2180
03100GOSUB 2300
03110V4=1
03120RETURN
03130PRINT "-------------------------------------------------------"
03140PRINT "WHEN YOU ARE READY TO CONTINUE TYPE '1'.";
03150GOSUB 9000
03160PRINT L$
03170PRINT "TYPE THE NUMBER OF THE OPTION YOU WANT."
03180PRINT " 1. FURTHER COMPARE THESE TWO PROPORTIONS"
03190PRINT " 2. COMPARE TWO OTHER PROPORTIONS"
03200PRINT " 3. EXIT MODULE"
03210GOSUB 9000
03220IF O1=3 THEN 3320
03230IF O1=2 THEN 3270
03240IF O1=1 THEN 3300
03250PRINT "REENTER. INPUT MUST BE NUMBER OF OPTION."
03260GOTO 3210
03270PRINT L$
03280V4=0
03290GOTO 320
03300PRINT L$
03310GOTO 1370
03320CHAIN "RSTRT"
05000REM ***************************************************
05005REM BETA CDF ROUTINE
05010REM INPUT A B J2
05015REM OUTPUT P
05020REM GOSUB'S TO BE CALLED PRIOR 5160 AND 5220
05025 DIM W(16),O(16)
05030GOSUB 5340
05035IF A+B>85 THEN 5280
05040P=0
05045C6=0
05050IF A>1 THEN 5080
05055C6=A
05060C7=B
05065A=C7
05070B=C6
05075J2=1-J2
05080D0=(J2-J1)*.5
05085D1=(J1+J2)*.5
05090FOR I1=1 TO 16
05095D9=D0*O(I1)+D1
05100IF D9=0 THEN 5115
05105IF D9=1 THEN 5115
05107D9=LOG(D9)*(A-1)+LOG(1-D9)*(B-1)
05108IF D9<-80 THEN 5115
05110P=P+W(I1)*EXP(D9)
05115NEXT I1
05120P=P*F0
05125P=P*D0
05130IF C6=0 THEN 5155
05135A=C6
05140B=C7
05145P=1-P
05150J2=1-J2
05155RETURN
05160FOR I1=1 TO 16
05165READ W(I1),O(I1)
05170NEXT I1
05175DATA 2.71525E-02,-.989401
05180DATA 6.22535E-02,-.944575,9.51585E-02,-.865631
05185DATA .124629,-.755404,.149596,-.617876
05190DATA .169156,-.458017,.182603,-.281604,.189451,-9.50125E-02
05195DATA .189451,9.50125E-02,.182603,.281604,.169156,.458017
05200DATA .149596,.617876,.124629,.755404
05205DATA 9.51585E-02,.865631,6.22535E-02,.944575,2.71525E-02
05210DATA .989401
05215RETURN
05220G9=A+B
05225GOSUB 5850
05230F0=G0
05235G9=A
05240GOSUB 5850
05245F0=F0-G0
05250G9=B
05255GOSUB 5850
05260F0=F0-G0
05265IF A+B>85 THEN 5275
05270F0=EXP(F0)
05275RETURN
05280W1=(B*J2)**(1/3)
05285W2=(A*(1-J2))**(1/3)
05290GOSUB 5325
05295I1=P
05300W1=(B*J1)**(1/3)
05305W2=(A*(1-J1))**(1/3)
05310GOSUB 5325
05315P=I1-P
05320RETURN
05325Y3=3*(W1*(1-1/9/B)-W2*(1-1/9/A))/SQR(W1*W1/B+W2*W2/A)
05330GOSUB 8000
05335RETURN
05340REM 2/16/76 CHANGED TO ALL LOG
05345D2=F0+(A-1)*LOG(J2)+(B-1)*LOG(1-J2)
05350IF D2<-80 THEN 5370
05355IF D2>85 THEN 5380
05360D2=EXP(D2)
05365RETURN
05370D2=1.E-37
05375RETURN
05380D2=1.E+37
05385RETURN
05390REM
05395REM END OF BETA CDF ROUTINE
05400REM *******************************************************
05850REM ****************************************************
05852REM LOG GAMMA ROUTINE
05853REM INPUT G9
05854REM OUTPUT G0
05860G5=G9
05863IF G9 <= 1.E+30 THEN 5872
05866G0=1.E+38
05869RETURN
05872IF G9>1.E-09 THEN 5881
05875G0=0
05878RETURN
05881IF G9<1.E+10 THEN 5890
05884G0=G9*(LOG(G9)-1)
05887RETURN
05890G6=1
05893IF 18<G5 THEN 5905
05896G6=G6*G5
05899G5=G5+1
05902GOTO 5893
05905R8=1/G5**2
05908G0=(G5-.5)*LOG(G5)-G5+.918939-LOG(G6)
05911C1=8.33333E-02
05914C2=2.77778E-03
05917C3=7.93651E-04
05920C4=5.95238E-04
05923G0=G0+1/G5*(C1-(R8*(C2+(R8*(C3-(R8*(C4)))))))
05926RETURN
05927REM END OF LOG GAMMA ROUTINE
05928REM ****************************************************
08000REM **********************************************************
08001REM ROUTINE CALCULATES THE CDF FOR NORMAL DISTRIBUTION
08002REM INPUT Y3
08003REM OUTPUT P
08004REM
08005Y4=ABS(Y3)
08010X1=X
08015X=Y3
08020T=1/(1+.231642*Y4)
08021IF X*X/2<80 THEN 8025
08022D=0
08023GOTO 8030
08025D=.398942*EXP(-X*X/2)
08030C1=1.33027
08035C2=1.82126
08040C3=1.78148
08045C4=.356564
08050C5=.319382
08055P=1-D*T*((((C1*T-C2)*T+C3)*T-C4)*T+C5)
08060IF X >= 0 THEN 8070
08065P=1-P
08070X=X1
08075RETURN
08076REM
08077REM END OF NORMAL CDF ROUTINE
08078REM **********************************************************
09000REM--SUBROUTINE THAT DETERMINES IF RESTART HAS BEEN REQUESTED.
09005INPUT O1
09015IF O1=-9999 THEN 9025
09020RETURN
09025CHAIN "RSTRT"
09035REM*************END ROUTINE
09999 END