Graff-Roeloffs' Modification of Dorfman
/****************************************************************/
/* S A S S A M P L E L I B R A R Y */
/* */
/* NAME: IEGRAFF */
/* TITLE: Graff-Roeloffs' Modification of Dorfman */
/* PRODUCT: QC */
/* SYSTEM: ALL */
/* KEYS: Inspection Sampling, */
/* PROCS: TABULATE */
/* DATA: */
/* */
/* MISC: */
/* */
/* NOTES: This program tabulates the probabilities of correct */
/* classification PC(NC) and PC(C) of nonconforming */
/* and conforming items respectively, and the percent */
/* reduction in expected overall number of tests for */
/* the modified Dorfman procedure proposed by Graff */
/* and Roeloffs. This consists of repeating inspec- */
/* tion until either r1 NC decisions or r2 C decisions */
/* have been obtained, whichever comes first. Here a */
/* hierarchal version is considered. */
/* */
/* Notation: */
/* */
/* k = stage */
/* */
/* r1j = r1 for jth stage */
/* r2j = r2 for jth stage */
/* */
/* n0 = number of items in group */
/* n1 = number of items in primary subgroup */
/* n2 = number of items in secondary subgroup */
/* */
/* omega = proportion of NC items in group */
/* */
/* p = Pr[ nonconforming item is classified */
/* as nonconforming ] */
/* pprime = Pr[ conforming item is classified as */
/* nonconforming ] */
/* */
/* */
/* REF: Johnson, N. L., Kotz, S., and Rodriguez, R. N. */
/* (1990), Statistical Effects of Imperfect Inspection */
/* Sampling: IV. Modified Dorfman Screening Pro- */
/* cedures, Journal of Quality Technology 22, 128-137. */
/* See table on page 133. */
/* */
/* Johnson, N. L., Kotz, S., and Wu, X. (1991). */
/* Inspection Errors for Attributes in Quality */
/* Control. London: Chapman & Hall. See Chapter 9. */
/* */
/****************************************************************/
data table;
keep pcc pcnc omega k n0 n1 n2 n3 r10 r20 p pp outersum pct;
format p 3.2 pp 3.2 pct 5.2 omega 4.2 pcc 5.4 pcnc 5.4 ;
label outersum = 'Expected/Number/of Tests'
pcc = 'PC(C)'
pcnc = 'PC(NC)'
k = 'k'
omega = 'omega'
p = 'p'
pp = 'p'''
pct = 'Percent/Reduction'
n0 = 'n0'
n1 = 'n1'
n2 = 'n2'
n3 = 'n3'
r10 = 'r10'
r20 = 'r20'
;
/*---set basic parameters---*/
do omega=0.05 to 0.10 by 0.05;
do k=0 to 2 by 1;
do r10 = 1 to 2 by 1;
do r20 = 1 to 2 by 1;
r11 = 1; r21 = 1;
r12 = 1; r22 = 1;
r13 = 1; r23 = 1;
/*---set parameters for group sizes---*/
n0 = 12; h0 = 1;
if k=0 then do;
n1 = 1; h1 = 12;
n2 = .; h2 = .;
n3 = .; h3 = .;
end;
else
if k=1 then do;
n1 = 4; h1 = 3;
n2 = 1; h2 = 4;
n3 = .; h3 = .;
end;
else
if k=2 then do;
n1 = 4; h1 = 3;
n2 = 2; h2 = 2;
n3 = 1; h3 = 2;
end;
/*---loop over p error probabilities---*/
do p=0.90 to 0.95 by 0.05;
/*---compute p0* ---*/
if (p<1.0) then do;
p0s = probbnml(1-p,r10+r20-1,r20-1);
p1s = probbnml(1-p,r11+r21-1,r21-1);
p2s = probbnml(1-p,r12+r22-1,r22-1);
p3s = probbnml(1-p,r13+r23-1,r23-1);
end;
else do;
p0s = 1;
p1s = 1;
p2s = 1;
p3s = 1;
end;
/*---loop over p' error probabilities---*/
do pp=0.05 to 0.10 by 0.05;
if pp>0.0 then do;
p0ps = probbnml(1-pp,r10+r20-1,r20-1);
p1ps = probbnml(1-pp,r11+r21-1,r21-1);
p2ps = probbnml(1-pp,r12+r22-1,r22-1);
p3ps = probbnml(1-pp,r13+r23-1,r23-1);
end;
else do;
p0ps = 0;
p1ps = 0;
p2ps = 0;
p3ps = 0;
end;
/*---sum over i ---*/
outersum = 0.0;
do i=0 to k+1;
/*---evaluate product---*/
if i=0 then iproduct = 1;
else if i=1 then iproduct = h1;
else if i=2 then iproduct = h1*h2;
else if i=3 then iproduct = h1*h2*h3;
/*---evaluate first inner term---*/
iterm1 = 0.0;
do u=0 to i;
/* evaluate diff = p0(nu) - p0(nu-1) */
if u=0 then diff=(1-omega)**n0;
else if u=1 then diff=(1-omega)**n1 - (1-omega)**n0;
else if u=2 then diff=(1-omega)**n2 - (1-omega)**n1;
else if u=3 then diff=(1-omega)**n3 - (1-omega)**n2;
else put 'error in diff for u > 3';
/* evaluate prod1 = product of pv* */
if u=0 then prod1=1;
else if u=1 then prod1=p0s;
else if u=2 then prod1=p0s*p1s;
else if u=3 then prod1=p0s*p1s*p2s;
else put 'error in prod1 for u > 3';
/*---evaluate prod2 = product of pw'* ---*/
if i=0 then do;
if u=0 then prod2 = 1.0 ;
else put 'error in prod2 for i=0';
end;
else if i=1 then do;
if u=0 then prod2 = p0ps;
else if u=1 then prod2 = 1;
else put 'error in prod2 for i=1';
end;
else if i=2 then do;
if u=0 then prod2 = p0ps*p1ps;
else if u=1 then prod2 = p1ps;
else if u=2 then prod2 = 1;
else put 'error in prod2 for i=2';
end;
else if i=3 then do;
if u=0 then prod2 = p0ps*p1ps*p2ps;
else if u=1 then prod2 = p1ps*p2ps;
else if u=2 then prod2 = p2ps;
else if u=3 then prod2 = 1;
else put 'error in prod2 for i=3';
end;
else put 'error for i > 4';
iterm1 = iterm1 + diff*prod1*prod2;
end; /* finish summing over u */
/*---find eip = ei' ---*/
if i=0 then do;
if pp>0.0 then
eip = r20*probbnml(pp ,r10+r20,r10-1)/(1-pp) +
r10*probbnml(1-pp,r10+r20,r20-1)/pp;
else
eip = r20 ;
end;
else
if i=1 then do;
if pp>0.0 then
eip = r21*probbnml(pp ,r11+r21,r11-1)/(1-pp) +
r11*probbnml(1-pp,r11+r21,r21-1)/pp;
else
eip = r21 ;
end;
else if i=2 then do;
if (pp>0.0) then
eip = r22*probbnml(pp ,r12+r22,r12-1)/(1-pp) +
r12*probbnml(1-pp,r12+r22,r22-1)/pp;
else
eip = r22 ;
end;
else
if i=3 then do;
if pp>0.0 then
eip = r23*probbnml(pp ,r13+r23,r13-1)/(1-pp) +
r13*probbnml(1-pp,r13+r23,r23-1)/pp;
else
eip = r23 ;
end;
else put 'error in eip for i > 3';
iterm1 = iterm1*eip;
/*---evaluate comp = 1-p0(ni) ---*/
if i=0 then comp = 1 - (1 - omega)**n0;
else if i=1 then comp = 1 - (1 - omega)**n1;
else if i=2 then comp = 1 - (1 - omega)**n2;
else if i=3 then comp = 1 - (1 - omega)**n3;
else put 'error in comp for i > 3';
/*---evaluate prod3 = product of pv'* ---*/
if i=0 then prod3 = 1;
else if i=1 then prod3 = p0s;
else if i=2 then prod3 = p0s*p1s;
else if i=3 then prod3 = p0s*p1s*p2s;
else put 'error in prod3 for i > 3';
/*---find ei = ei---*/
if i=0 then do;
if (p<1.0) then
ei = r20*probbnml(p ,r10+r20,r10-1)/(1-p) +
r10*probbnml(1-p,r10+r20,r20-1)/p;
else
ei = r10;
end;
else
if i=1 then do;
if (p<1.0) then
ei = r21*probbnml(p ,r11+r21,r11-1)/(1-p) +
r11*probbnml(1-p,r11+r21,r21-1)/p;
else
ei = r11;
end;
else
if i=2 then do;
if (p<1.0) then
ei = r22*probbnml(p ,r12+r22,r12-1)/(1-p) +
r12*probbnml(1-p,r12+r22,r22-1)/p;
else
ei = r12;
end;
else
if i=3 then do;
if (p<1.0) then
ei = r23*probbnml(p ,r13+r23,r13-1)/(1-p) +
r13*probbnml(1-p,r13+r23,r23-1)/p;
else
ei = r13;
end;
else put 'error in prod3 for i > 3';
iterm2 = comp*prod3*ei;
/*---update outer sum---*/
outersum = outersum + iproduct*(iterm1 + iterm2);
end; /* over i */
/*---percent savings---*/
pct = 100.0*(1 - outersum/n0);
/* pcnc = pc( nc ) */
if k = 0 then pcnc = p0s * p1s;
else
if k = 1 then pcnc = p0s * p1s * p2s;
else
if k = 2 then pcnc = p0s * p1s * p2s * p3s;
else put 'error for pcnc with k > 2';
/*---pcc = pc( c )---*/
sumpcc = 0.0;
do s = 0 to k + 1;
/*---compute diff = p0*(ns) - p0*(ns-1)---*/
if s = 0 then diff = (1-omega)**(n0-1);
else
if s = 1 then do;
if s=k+1 then
diff = 1 -(1-omega)**(n0-1);
else
diff = (1-omega)**(n1-1)-(1-omega)**(n0-1);
end;
else
if s = 2 then do;
if s=k+1 then
diff = 1 -(1-omega)**(n1-1);
else
diff = (1-omega)**(n2-1)-(1-omega)**(n1-1);
end;
else
if s = 3 then do;
if s=k+1 then
diff = 1 -(1-omega)**(n2-1);
else
diff = (1-omega)**(n3-1)-(1-omega)**(n2-1);
end;
else put 'error for pcc diff with s > 3';
/*---prod1 = product of pj* ---*/
if s = 0 then prod1 = 1;
else if s = 1 then prod1 = p0s ;
else if s = 2 then prod1 = p0s * p1s ;
else if s = 3 then prod1 = p0s * p1s * p2s ;
else put 'error for pcc prod1 with s > 3 ';
/*---prod2 = product of pj* ---*/
if s = 0 then do;
if k = 0 then prod2 = p0ps * p1ps;
else if k = 1 then prod2 = p0ps * p1ps * p2ps;
else if k = 2 then prod2 = p0ps * p1ps * p2ps *p3ps;
else put 'error with s = 0 and k > 2';
end;
else if s = 1 then do;
if k = 0 then prod2 = p1ps;
else if k = 1 then prod2 = p1ps * p2ps;
else if k = 2 then prod2 = p1ps * p2ps * p3ps;
else put 'error with s = 1 and k > 2';
end;
else if s = 2 then do;
if k = 0 then prod2 = 1;
else if k = 1 then prod2 = p2ps;
else if k = 2 then prod2 = p2ps * p3ps;
else put 'error with s = 2 and k > 2';
end;
else if s = 3 then do;
if k = 0 then prod2 = 1;
else if k = 1 then prod2 = 1;
else if k = 2 then prod2 = p3ps;
else put 'error with s = 3 and k > 2';
end;
else put 'error for pcc prod2 with s > 3 ';
sumpcc = sumpcc + diff * prod1 * prod2;
end; /* over s */
pcc = 1 - sumpcc;
output;
end; /* over pp */
end; /* over p */
end; /* over r10 */
end; /* over r20 */
end; /* over k */
end; /* over omega */
run;
proc sort; by n0 omega k;
proc print label noobs split='/' ;
var n1 n2 n3 r10 r20 p pp pct pcc pcnc;
by n0 omega k;
run;