LCL for Cpk (Guirguis, Rodriguez)
/****************************************************************/
/* S A S S A M P L E L I B R A R Y */
/* */
/* NAME: CPKCON2 */
/* TITLE: LCL for Cpk (Guirguis, Rodriguez) */
/* PRODUCT: QC */
/* SYSTEM: ALL */
/* KEYS: Capability Analysis, Capability Indices, */
/* PROCS: CAPABILITY IML */
/* DATA: */
/* */
/* REF: Youn-Min Chou, D. B. Owen and S. A. Borrego (1990). */
/* "Lower Confidence Limits on Process Capability */
/* Indices". Journal of Quality Technology 22, */
/* pp. 223-229. Corrigenda, JQT 24, p. 251. */
/* */
/* Kushler, R. H. and Hurley, P. (1992). "Confidence */
/* Bounds for Capability Indices". Journal of Quality */
/* Technology 24, pp. 188-195. */
/* */
/* Guirguis, G. H. and Rodriguez, R. N. (1992). */
/* "Computation of Owen's Q Function Applied to */
/* Process Capability Analysis". Journal of Quality */
/* Technology 24, pp. 236-246. */
/* */
/* NOTES: This program computes lower confidence limits for */
/* Cpk using the method of COB (1990). The program */
/* is more general than Table 5 of COB in that it */
/* computes the confidence limit as a function of */
/* CPL, CPU, and n, following the approach of Guirguis */
/* and Rodriguez (1992). Note that this method can */
/* yield conservative results as pointed out by */
/* Kushler and Hurley (1992) and Guirguis and */
/* Rodriguez (1992). */
/* */
/* MISC: */
/* */
/****************************************************************/
options ps=60 ls=80;
data a;
input x @@;
cards;
1.38 1.49 1.43 1.60 1.59
1.34 1.44 1.64 1.83 1.57
1.45 1.74 1.61 1.39 1.63
1.73 1.61 1.35 1.51 1.47
1.46 1.41 1.56 1.40 1.58
1.43 1.53 1.53 1.58 1.62
1.58 1.46 1.26 1.57 1.41
1.53 1.36 1.63 1.36 1.66
1.49 1.55 1.67 1.41 1.39
1.75 1.37 1.36 1.86 1.49
;
proc capability data=a noprint;
spec lsl=0.8 usl=2.4;
var x;
output out = summary
n = n
mean = mean
std = std
cpu = cpu
cpl = cpl
cpk = cpk
;
proc iml;
/****************************************************************/
/* */
/* IML MODULES */
/* */
/****************************************************************/
/* Integrand */
start fun(x) global(a1,a2,b,nn,mlog);
fx = probnorm(a1*x-b) - probnorm(-a2*x+b);
fx = fx*exp((nn-2)*log(x)-(0.5*x*x)-mlog);
return(fx);
finish;
/* Integral */
start probd(cpl,cpu,n,ck) global(a1,a2,b,nn,mlog);
nn = n;
a1 = (3*cpl*sqrt(n))/sqrt(n-1);
a2 = (3*cpu*sqrt(n))/sqrt(n-1);
b = 3*ck*sqrt(n);
int = 2*sqrt(n-1)*ck/(cpl+cpu) ;
peak = int + 1.;
int = int || .P;
mlog = lgamma(0.5*(n-1))+0.5*(n-3)*log(2);
call quad( z1, "fun", int ) peak=peak;
return(z1);
finish;
/* Initial value for the root search */
start ini(cpl,cpu,n,level);
cko=0.0;
paux = 1.0;
aa:
if(paux < 0.99 | paux < level) then goto bb;
cko = cko+(3.0/sqrt(n));
paux=probd(cpl,cpu,n,cko);
goto aa;
bb:
return(cko);
finish;
/* Root search */
start findck(cpl,cpu,n,level);
ckr = ini(cpl,cpu,n,level);
count = 1;
del = 1;
do while( abs(del) > 1.e-5 & count < 100 );
p1 = probd(cpl,cpu,n,ckr);
p2 = probd(cpl,cpu,n,ckr*1.001);
if (p1 = . | p2 = .)
then
do;
ckr = .;
goto out;
end;
der = (p2-p1)*1000/ckr;
del = p1/der*log(level/p1);
ckr = ckr + del;
count = count+1;
end;
if ( count >= 100 )
then
cr = .;
out:
return(ckr);
finish;
/****************************************************************/
/* */
/* MAIN PROGRAM */
/* */
/****************************************************************/
use summary;
read all var {n,mean,std,cpu,cpl,cpk};
level = 0.95;
reset noname;
ck = findck(cpl,cpu,n,level);
print ,;
rowa = {" Sample Size ",
" Mean ",
" Standard Deviation "};
rowb = {" Lower Confidence Limit ",
" Cpk "};
auxa = n//mean//std;
auxb = ck//cpk;
mattrib auxa rowname=(rowa) format=6.2;
mattrib auxb rowname=(rowb) format=6.2;
print " Summary Statistics";
print auxa;
print " 95% Lower Confidence Limit for Cpk";
print auxb;
quit;
run;
title;