/**********************************************************************/ /* Appendix A: Listing of Macro SERIALT */ /* Macro: SERIALT */ /* Description: Executes serial test on RANUNI */ /* Parameters: */ /* DIM # of dimensions (DIM-tuples) */ /* INT # of intervals along each dimension */ /* SAMPLES # of samples of random numbers */ /* NOFRV # of random vectors per sample */ /* SEED Seed of function RANUNI */ %macro serialt(dim=,int=,samples=,nofrv=,seed=0); %do sample=1 %to &samples; /* Generate the random vectors. File work has */ /* variables xi's (i=1,2,...,&dim), where xi */ /* is the interval number in which the random */ /* vector falls along the ith dimension. */ data work(drop=i); do i=1 to &nofrv; %do i=1 %to &dim; x&i=1+int(&int*ranuni(&seed)); %end; output; end; /* Count number of vectors in each hypercube */ proc freq data=work; table x1 %if &dim>1 %then %do; %do i=2 %to &dim; * x&i %end; %end; / noprint out=work; /* Calculate the serial test */ data work(keep=prob); set work end=end; retain chisq 0 expfreq; if _n_=1 then expfreq=&nofrv/(&int**&dim); chisq=chisq+(count-expfreq)**2/expfreq; if end then do; prob=1-probchi(chisq,&int**&dim-1); output; end; run; /* Save the results of the serial test */ proc append base=results data=work; %end; proc format; value signif 0.0000-0.0001='P<0.0001' 0.0001-0.001 ='P<0.001' 0.001 -0.01 ='P<0.01' 0.01 -0.05 ='P<0.05' 0.05 -0.10 ='P<0.10' 0.10 -0.15 ='P<0.15' 0.15 -1 ='P>0.15'; /* Summarize the results */ proc freq data=results; table prob; format prob signif.; title1 "Result of the Serial Test"; title2 "-------------------------"; title3 "# of Samples: &samples"; title4 "# of Random Vectors/Sample: &nofrv"; title5 "# of Dimensions: &dim"; title6 "# of Intervals/Dimension: &int"; run; %mend; /**********************************************************************/ /* Appendix B: Listing of Macro TCOIN */ /* Macro: TCOIN */ /* Description: The macro simulates the */ /* coin toss and */ /* executes a special run test. */ /* Parameters: */ /* NTOSS Number of tosses */ /* P Probability of head */ /* SEED Seed of random number generator */ %macro tcoin(ntoss=,p=0.5,seed=0); data outcome(keep=outcome); array side(2) $1 sideh sidet; retain sideh 'H' sidet 'T'; do i=1 to &ntoss; outcome=side(rantbl(&seed,&p)); output; end; proc freq data=outcome; table outcome; /* Determine the runs */ data runs(keep=run); set outcome; by notsorted outcome; retain run; if first.outcome then run=0; run=run+1; if last.outcome then output; /* Produce statistics of the runs */ proc univariate data=runs freq; /* Produce run frequencies for the test */ proc freq data=runs; table run / out=runfreq noprint; /* Cut the run-distribution at 10 */ data runfreq; set runfreq end=end; retain rest 0; if run<=10 then do; output; return; end; rest=rest+count; if end then do; count=rest; run=11; output; end; /* Calculate the expected run frequencies */ data expfreq; retain cexpprob 0; do run=1 to 10; expprob=1/(2**run); expfreq=(&ntoss/2)*expprob; cexpprob=cexpprob+expprob; output; end; run=11; expfreq=(&ntoss/2)*(1-cexpprob); output; /* Calculate and print the goodness-of-fit test */ data _null_; merge runfreq expfreq; by run; file print; retain chisq 0; if _n_=1 then put 'RUN FREQ EXP.FREQ' / 24*'-'; if count=. then count=0; chisq=chisq+(count-expfreq)**2/expfreq; if run<=10 then put run 2. count 11. expfreq 11.; else put '11+' count 10. expfreq 11.; if run=11 then do; prob=1-probchi(chisq,10); put 24*'-' / 'RUN TEST:' chisq 15.4 / 'PROBABILITY:' prob 12.4; end; %mend; /**********************************************************************/ /* Appendix C: Listing of Macro TBCOIN */ /* Macro: TBCOIN */ /* Description: The macro simulates the */ /* "tossing of a biased coin" */ /* where the tosses are defined */ /* as two-coin flips */ /* Parameters: */ /* NTOSS Number of tosses */ /* P Probability of head */ /* SEED Seed of random number generator */ %macro tbcoin(ntoss=,p=,seed=0); data outcome(keep=toss1 toss2); array side(2) $1 sideh sidet; retain sideh 'H' sidet 'T'; do i=1 to &ntoss; toss1=side(rantbl(&seed,&p)); toss2=side(rantbl(&seed,&p)); if toss1^=toss2 then output; end; proc freq; table toss1*toss2/list; %mend; /**********************************************************************/ /* Appendix D: Listing of Macro MHS1 */ /* Macro: MHS1 */ /* Description: Simulation of the Monty Hall */ /* problem */ /* Parameters: */ /* HOST Behavior of host. Possible values: */ /* SO - switch is always offered */ /* NSO - switch is never offered */ /* M - host is malevolent (switch is */ /* offered only when initial */ /* choice is a winning one) */ /* SPROB Probability of player to switch when */ /* switch is offered (0<=SPROB<=1) */ /* OUT output data set of the simulation */ /* NPLAYS Number of plays */ /* SEED Seed of function RANUNI */ %macro mhs1(host=so,sprob=1,out=out,nplays=, seed=0); data &out(drop=i); array env(3) $1 env1-env3; length kenv showenv initch finalch 3; do i=1 to &nplays; env(1)='E'; env(2)='E'; env(3)='E'; /* Place the key into an envelope */ kenv=rantbl(&seed,1/3,1/3); env(kenv)='K'; /* Make initial choice */ initch=rantbl(&seed,1/3,1/3); initprz =env(initch); /* Host offers no switch */ showenv=.; %if &host=NSO %then %do; finalch=initch; %end; /* Host offers switch */ %if &host=SO %then %do; link show; link switch; %end; /* Host is malevolent */ %if &host=M %then %do; if env(initch)='E' then finalch=initch; else do; link show; link switch; end; %end; finalprz=env(finalch); output; end; return; %if &host=SO | &host=M %then %do; /* Choose an empty envelope to show */ show: if initch^=kenv then showenv=6-(initch+kenv); else do; showenv=rantbl(&seed,0.5); if kenv=1 | (kenv=2 & showenv=2) then showenv=showenv+1; end; return; /* Make the switch */ switch: if ranuni(&seed)<&sprob then finalch=6-(initch+showenv); else finalch=initch; return; %end; run; /* Present simulation results */ proc freq data=out; table finalprz initprz*finalprz env1*env2*env3 env1*env2*env3*initch*initprz*showenv* finalch*finalprz/list missing; %mend; /**********************************************************************/ /* Appendix E: /* Macro MHS2 */ /* Description: Simulates a generalized Monty */ /* Hall problem */ /* Parameters */ /* NEMPTY number of empty envelopes */ /* NKEY number of envelops with key */ /* restriction: NEMPTY+NKEY<=200 */ /* NREMOVE number of empty envelopes to be removed */ /* NPLAYS number of games to simulate */ /* OFFER probability of host offering the */ /* switch */ /* SWITCH probability of player taking the */ /* switch */ /* SEED seed of the random number generator */ /* OUT name of the output dataset which */ /* has one variable (PROB) */ %macro mhs2(nempty=,nkey=,nremove=,nplays=, offer=,switch=,seed=0,out=); %let nenv=%eval(&nempty+&nkey); data &out(keep=prob); array env(&nenv) $1 env1-env&nenv; length pack0 pack $&nenv; /* Make the empty envelopes */ do i=1 to &nempty; env(i)='E'; end; /* Make the envelopes with key */ do i=&nempty+1 to &nenv; env(i)='K'; end; /* Shuffle the envelopes */ do i=1 to &nenv-1; r=int((&nenv-i+1)*ranuni(&seed))+1; temp=env(&nenv-i+1); env(&nenv-i+1)=env(r); env(r)=temp; end; /* pack0 saves the contents of all envelopes */ pack0=env(1); do i=2 to &nenv; pack0=trim(pack0)||env(i); end; /* Play the game &nplays times */ nwins=0; do p=1 to &nplays; pack=pack0; /* Make initial selection */ sel=int(&nenv*ranuni(&seed))+1; prize=substr(pack,sel,1); if ranuni(&seed)<&offer and ranuni(&seed)<&switch then do; /* Remove empty envelopes */ i=0; j=1; do while (i<&nremove); if j^=sel and env(j)='E' then do; i=i+1; substr(pack,j,1)=' '; end; j=j+1; end; substr(pack,sel,1)=' '; pack=compress(pack); /* Make final selection */ sel=int((&nenv-&nremove-1)*ranuni(&seed))+1; prize=substr(pack,sel,1); end; if prize='K' then nwins=nwins+1; end; /* Output the probability of winning */ prob=nwins/&nplays; output; run; %mend; /**********************************************************************/