macro
twosamplepsize p1 p2;
   mfg alpha;
   cons beta.

#This macro finds the smallest values of n and dx for the two sample test
#to compare fractions defective. The user must specify the fractions defective
#p1 and p2 and the corresponding type 1 and 2 error probabilities. This
#macro provides the exact solution analogous to the normal approximation method.
#The hypotheses tested are Ho:p1=p2 vs. Ha:p1<p2. p1 must always be 
#less than p2. If you want to test Ha:p1>p2 enter the smaller fraction
#as p1, the larger as p2, and swap the alpha and beta values. Follow this macro
#with the occurve.mac macro to create the OC curve.

#DO NOT USE THIS MACRO UNLESS YOU UNDERSTAND IT AND KNOW WHAT YOU ARE DOING!
#FISHERS METHOD IS CONSERVATIVE AND DOES NOT HAVE THE POWER FOR THE CONDITIONS
#RETURNED BY THIS MACRO. USE FISHERSPOWER AND FISHERS TEST TO DESIGN AND
#ANALYZE YOUR EXPERIMENT INSTEAD.

#If the null hypothesis is true then both p1 and p2 are equal to p1.
#If the null hypothesis if false then p2 has shifted to a larger value than p1.

#By: Paul G. Mathews, 217 Third Street, Fairport Harbor, OH 44077
#pmathews@apk.net
#440-350-0911
#Copyright 3 April 2001 All rights reserved
#Rev. 1.1 18July2001 PGM Added clarifying comments

#Warning: This macro has known bugs and does not always give correct answers.
#Use at your own risk. PGMathews takes no responsibility for use or misuse.

#Warning: This macro runs very very slowly, however, it might be worth it.

#Example call:
#   mtb> %twosamplepsize 0.01 0.10

#%twosamplepsize 0.01 0.10;
#mfg 0.05;
#cons 0.20.

mconstant n p1 p2 dx alpha beta
mconstant notalpha notbeta Pp1eqp2 Pp1nep2 zalpha zbeta
default alpha=0.05
default beta=0.10
notitle

let notalpha=1-alpha
let notbeta=1-beta

note
note The macro is running. This may take a while ...
note

#Find the normal approximation sample size as a starting point
invcdf alpha zalpha
invcdf beta zbeta
let n=round(((zalpha+zbeta)/(p2-p1))**2*(p1*(1-p1)+p2*(1-p2)))
#This n is a little smaller than Minitab's because Minitab uses pbar=(p1+p2)/2 to estimate sigma-p.
#To duplicate Minitab's result use: let n=round(((zalpha+zbeta)/(p2-p1))**2*2*pbar*(1-pbar)
#where let pbar=(p1+p2)/2
note
note The normal approximation sample size is:
print n
note

if p1>p2
   note p1 must be less than or equal to p2. Please reformulate your hypotheses.
   exit
endif

let dx=2				#dx will be initialized to 0, this is for dx memory/speed
let Pp1nep2=0				#(p1,notalpha) condition is met already
while (Pp1nep2<beta)	
   let Pp1eqp2=notalpha
   let n=n-10
   let dx=dx-2					#This is dx memory from the preceeding iteration
   while (Pp1eqp2<=notalpha) and (dx<=n)	#Find the smallest dx that meets alpha condition
      let dx=dx+1
      call Ponedx n dx p1 p1 Pp1eqp2
   endwhile
   call Ponedx n dx p1 p2 Pp1nep2		#Now find the corresponding beta value
endwhile					#When this loop ends alpha condition is met, beta is not.
			
while (Pp1nep2>beta)
   let n=n+1					#n is too small so come up from below
   call Ponedx n dx p1 p1 Pp1eqp2
   call Ponedx n dx p1 p2 Pp1nep2
endwhile

note The input parameters were:
print p1 p2 alpha beta
note
let alpha=1-Pp1eqp2			#Warning: Overwriting alpha and beta for output!
let beta=Pp1nep2
note The sample size, acceptance number, and actual alpha and beta are:
print n dx alpha beta			#Print dx, 1-alpha, and beta

title
brief 1
endmacro


#*********************************************
macro
   Ponedx nn DxA p1 p2 PDxA

#This macro finds P(dx<=DxA;n,p1,p2)

#Example call:
#   mtb> %Ponedx 40 2 0.01 0.10 k1

mconstant nn dx DxA p1 p2 Px1 PDxA bottom
mcolumn x1 x2 bxnp1 bxnp2 bxnp 

brief 0
let bottom=-nn
let PDxA=0
set x1
   0:nn
end
do dx=bottom:DxA
   erase bxnp1 bxnp2
   pdf x1 bxnp1;
      binomial nn p1.
   let x2=x1+dx
   pdf x2 bxnp2;
      binomial nn p2.
   let bxnp=bxnp1*bxnp2
   sum bxnp Px1
   let PDxA=PDxA+Px1
enddo

brief 1
endmacro