Moments of Estimator of Cp

 /****************************************************************/
 /*       S A S   S A M P L E   L I B R A R Y                    */
 /*                                                              */
 /*    NAME: CPMOM                                               */
 /*   TITLE: Moments of Estimator of Cp                          */
 /* PRODUCT: QC                                                  */
 /*  SYSTEM: ALL                                                 */
 /*    KEYS: Capability Analysis, Capability Indices,            */
 /*   PROCS: TABULATE                                            */
 /*    DATA:                                                     */
 /*     REF: W. L. Pearn, S. Kotz, N. L. Johnson (1992).         */
 /*          "Distributional and Inferential Properties of       */
 /*          Process Capability Indices".  Journal of Quality    */
 /*          Technology 24, pp. 216-231.                         */
 /*                                                              */
 /*   NOTES: This program computes the expected value and        */
 /*          standard deviation for the estimator of Cp.         */
 /*                                                              */
 /*    MISC:                                                     */
 /*                                                              */
 /****************************************************************/

   options ps=60 ls=100;

   data cp;
     keep ds f ex va st;
     label f  = 'f'
           ds = 'd/sigma'
           ex = '  Exp  '
           st = '  Std  '
           va = '  Var  ';
     do ds = 2 to 6;
       do f = 4 to 59 by 5;
         ex = (gamma((f-1)/2)/gamma(f/2))*sqrt((f/2));
         ex = (ds/3)*ex;
         va = (gamma((f-1)/2)/gamma(f/2))**2;
         va = (f/(f-2))-((f/2)*va);
         va = ((ds/3)**2)*va;
         st = sqrt(va);
         output;
       end;
     end;
   run;

   proc tabulate data = cp noseps;
      class f ds;
      var ex va st;
      table f, ds*(ex st)*sum=' '*f = 7.4 / rts = 9;
   run;