1 10 REM BAYES
  2 12REM WRITTEN BY R.C.HEWITT, JAN 7,1969
  3 17 DIM Y(200)
  4 100 GO TO 675
  5 101 GO TO 440
  6 199 GO TO 300
  7 200 DATA .50,.75,.90,.95,.999,.9999,.99999
  8 299 DATA 10
  9 300 READ N1
 10 305 LET Y1=0
 11 310 LET N3=N1-1
 12 315 FOR I=1 TO N1
 13 320 READ Y(I)
 14 325 LET Y1=Y1+Y(I)
 15 330 NEXT I
 16 335 LET Y3=Y1/N1
 17 340 LET Y1=0
 18 345 FOR I=1 TO N1
 19 350 LET Y1=Y1+(Y(I)-Y3)*(Y(I)-Y3)
 20 355 NEXT I
 21 360 LET S1=Y1/N1
 22 365 READ N2
 23 370 LET Y1=0
 24 375 LET N4=N2-1
 25 380 FOR I=1 TO N2
 26 385 READ Y(I)
 27 390 LET Y1=Y1+Y(I)
 28 395 NEXT I
 29 400 LET Y4=Y1/N2
 30 405 LET Y1=0
 31 410 FOR I=1 TO N2
 32 415 LET Y1=Y1+(Y(I)-Y4)*(Y(I)-Y4)
 33 420 NEXT I
 34 425 LET S2=Y1/N2
 35 430 LET Y5=Y3-Y4
 36 435 GO TO 470
 37 440 READ Y,N1,Y3,S1,Z,N2,Y4,S2
 38 445 LET N3=N1-1
 39 450 LET N4=N2-1
 40 455 LET Y5=Y3-Y4
 41 460 LET S1=S1*S1
 42 465 LET S2=S2*S2
 43 470 PRINT TAB(8);"STATISTIC";TAB(30);"SAMPLE 1";TAB(45);"SAMPLE 2"
 44 475 PRINT
 45 480 PRINT "SAMPLE MEAN";TAB(30);Y3,Y4
 46 485 PRINT "SAMPLE STD DEVIATION";TAB(30);SQR(S1),SQR(S2)
 47 490 PRINT "SAMPLE SIZE";TAB(30);N1,N2
 48 495 PRINT "SAMPLE VARIANCE";TAB(30);S1,S2
 49 500 PRINT
 50 505 PRINT "DIFF OF MEANS";TAB(37);Y5
 51 510 PRINT
 52 515 PRINT
 53 520 PRINT "HIGHEST POSTERIOR DENSITY INTERVALS"
 54 521 PRINT "DIFFERENCE BETWEEN MEANS:"
 55 525 PRINT
 56 530 PRINT" HPD LEVEL","LOWER LIM","UPPER LIM"
 57 535 PRINT
 58 540 READ P
 59 545 IF P=10 THEN 885
 60 550 LET P=(1-P)/2
 61 555 LET T=SQR(LOG(1/P^2))
 62 560 LET X=2.515517+0.802853*T+0.010328*T*T
 63 565 LET X=T-X/(1+1.432788*T+0.189269*T*T+0.001308*T*T*T)
 64 570 LET G1=X*(1+X*X)/4
 65 575 LET G2=X*(1+16*X*X/3*(1+5*X*X/16))/32
 66 580 LET G3=-5*X*(1-17*X*X/15*(1+19*X*X/17*(1+3*X*X/19)))/128
 67 585 LET G4=1-1482*X*X/1920*(1+776*X*X/1482*(+79*X*X/776))
 68 590 LET G4=-945*X*(1+1920*X*X/245*G4)/92160
 69 595 LET S4=S2/N2
 70 600 LET S3=S1/N1
 71 605 LET C1=S4/(S3+S4)
 72 610 LET C2=S3/(S3+S4)
 73 615 LET N6=N4/(N4-2)
 74 620 LET N5=N3/(N3-2)
 75 625 LET F=(N6*C1+N5*C2)*(N6*C1+N5*C2)
 76 630 LET F=F/((N6*N6*C1*C1)/(N4-4)+(N5*N5*C2*C2)/(N3-4))
 77 635 LET B=F+4
 78 640 LET A=(B-2)*(N6*C1+N5*C2)/B
 79 645 LET T5=X+G1/B+G2/B^2+G3/B^3+G4/B^4
 80 650 LET T6=SQR(A*(S3+S4))*T5
 81 655 LET Y6=Y5-T6
 82 660 LET Y9=Y5+T6
 83 665 PRINT (1-2*P)*100,Y8,Y9
 84 670 GO TO 540
 85 675 PRINT " THIS PROGRAM COMPUTES THE HIGHEST POSTERIOR DENSITY "
 86 680 PRINT "INTERVALS FOR THE DIFFERENCE BETWEEN TWO POPULATION MEANS"
 87 685 PRINT "USING BAYESIAN STATISTICAL INFERENCE.  THE DATA MAY BE"
 88 690 PRINT "ENTERED IN EITHER OF TWO WAYS:"
 89 695 PRINT
 90 700 PRINT "     A. TYPE THE DATA ON ONE LINE AS:"
 91 705 PRINT
 92 710 PRINT "      100 DATA H1,N1,M1,S1,H2,N2,M2,S2"
 93 715 PRINT "      RUN"
 94 720 PRINT
 95 725 PRINT "     WHERE,"
 96 730 PRINT "       H1=DUMMY TO AGREE WITH CONDIF*** INPUT"
 97 735 PRINT "          (ENTER '0' AS AN APPROPRIATE NUMBER)"
 98 740 PRINT "       N1=SIZE OF SAMPLE 1"
 99 745 PRINT "       M1=ARITHMETIC MEAN OF SAMPLE 1"
100 750 PRINT "       S1=STANDARD DEVIATION OF SAMPLE 1(BASED ON"
101 755 PRINT "          DIVISOR OF N1)"
102 760 PRINT
103 765 PRINT "      H2,N2,M2,S2, ARE THE SAME FOR SAMPLE 2"
104 770 PRINT
105 775 PRINT "     B. TYPE THE DATA ON MORE THAN ONE LINE AS:"
106 780 PRINT
107 785 PRINT "      100 DATA N1,Y(1),Y(2),Y(3),Y(4)"
108 790 PRINT "      101 DATA Y(5),   . . .    ,Y(N1)"
109 795 PRINT "      102 DATA N2,Z(1),Z(2),Z(3),Z(4)"
110 800 PRINT "      103 DATA Z(5),   . . .    ,Z(N2)"
111 805 PRINT "      RUN"
112 810 PRINT
113 815 PRINT "     WHERE,"
114 820 PRINT "       N1=THE SIZE OF SAMPLE 1"
115 825 PRINT "       N2=THE SIZE OF SAMPLE 2"
116 830 PRINT "       Y(I),I=1 TO N1  ARE THE ELEMENTS OF SAMPLE 1"
117 835 PRINT "       Z(I),I=1 TO N2  ARE THE ELEMENTS OF SAMPLE 2"
118 840 PRINT
119 845 PRINT
120 850 PRINT " IF DIFFERENT HIGHEST POSTERIOR DENSITY INTERVALS ARE"
121 855 PRINT "DESIRED, THEN THEY MAY BE ENTERED AS:"
122 860 PRINT
123 865 PRINT "       200 DATA L1,L2,L3,L4,L5,..."
124 870 PRINT
125 875 PRINT "       WHERE ALL L=L1,L2,L3,... ARE SUCH THAT"
126 880 PRINT "            0.5< L < 1.0"
127 885 STOP
128 890 END