SigmaPlot Product Uses - Transform to Compute Shelf Life

SigmaPlot Transform to Compute Shelf Life

' TRANSFORM TO COMPUTE SHELF LIFE
     'X DATA IS PLACED IN COLUMN X_COL (3 OR MORE DATA POINTS)
     'Y DATA IS PLACED IN COLUMN Y_COL
     'RESULTS ARE PLACED IN COLUMNS RES THROUGH RES+7
        'COLS RES & RES+1 CONTAIN THE REGRESSION LINE
        'COL RES+2 CONTAINS THE LOWER CONFIDENCE LINE
        'COLS RES+3 & RES+4 CONTAIN THE (T90,90) POINT FOR THE DROP LINES
        'COLS RES+5 & RES+6 CONTAIN THE T90 VALUE
        'COL RES+7 IS A WORKING COLUMN
     'COPY A PAGE TEMPLATE ON YOUR DATA TO CREATE THE GRAPH
  ' INPUT
  X_COL=1           'COLUMN NUMBER FOR X DATA
  Y_COL=2           'COLUMN NUMBER FOR Y DATA
  RES=4              'FIRST RESULTS COLUMN
  ' PROGRAM
  Z=1.96             'Z FOR 95% CONFIDENCE
  Y0 = 90            'Y SPECIFICATION LIMIT (%)
  X1=COL(X_COL)      'DEFINE X VALUES
  Y1=COL(Y_COL)      'DEFINE Y VALUES
     'ROWWISE DELETE MISSING VALUES
  FOR ROW = 1 TO SIZE(X1) DO
     CELL(RES+7,ROW)=MISSING(BLOCK(X_COL, ROW, Y_COL, ROW))
  END FOR
  X=IF(COL(RES+7)=0,X1)
  Y=IF(COL(RES+7)=0,Y1)
  N=SIZE(X)           'NUMBER OF DATA POINTS
  V=N-2                 'N MUST BE > 2
  XBAR=MEAN(X)                     'MEAN OF X
  DENOM=TOTAL((X-XBAR)^2)          'SUM OF SQS ABOUT MEAN
  ALPHA=TOTAL(X^2)/(N*DENOM)       '1,1 COEFF OF (X'X)^-1
  BETA=-XBAR/DENOM                 '1,2 COEFF OF (X'X)^-1
  DELTA=1/DENOM                    '2,2 COEFF OF (X'X)^-1
  R1=TOTAL(Y)                      '1ST ROW OF X'Y
  R2=TOTAL(X*Y)                    '2ND ROW OF X'Y
  B0=ALPHA*R1+BETA*R2              'INTERCEPT PARAMETER
  B1=BETA*R1+DELTA*R2              'SLOPE PARAMETER     'COMPUTE T VALUE
  T123=Z+(Z^3+Z)/(4*V)+(5*Z^5+16*Z^3+3*Z)/(96*V^2)
  T4=(3*Z^7+19*Z^5+17*Z^3-15*Z)/(384*V^3)
  T5=79*Z^9+776*Z^7+1482*Z^5-1920*Z^3-945*Z
  T6=27*Z^11+339*Z^9+930*Z^7-1782*Z^5-765*Z^3+17955*Z
  T1=T123+T4+T5/(92160*V^4)+T6/(368640*V^5)
  T=IF(V=1, 12.706, IF(V=2, 4.303,T1))
     'ESTIMATE OF S
  S=SQRT(TOTAL(((Y-(B0+B1*X))^2))/V)
     'QUADRATIC EQUATION FOR EXPIRATION TIME
  DELTA0 = B0-Y0
  A = DELTA - (B1/(T*S))^2
  B = 2*BETA - 2*B1*DELTA0/(T*S)^2
  C = ALPHA - (DELTA0/(T*S))^2
  B24AC=B^2-4*A*C
  ROOT1=(-B + SQRT(B24AC))/(2*A)
  ROOT2=(-B - SQRT(B24AC))/(2*A)
     'FIND APPROPRIATE ROOT
  MINX=0  R={ROOT1,ROOT2}
  ROOT=IF(S=0, (Y0-B0)/B1, MAX(IF(R<(Y0-B0)/B1,R)))  
MAXX=IF(S=0,1.1*ROOT, IF(B1>0 OR B24AC <0, MAX(X), 1.1*ROOT))
  XREG=DATA(MINX,MAXX,(MAXX-MINX)/20)
  'XREG=DATA(MINX,MAXX,1)
  'USED FOR COMPARISON WITH SIGMASTAT
  'XREG=X1                 ' USED FOR COMPARISON WITH SAS
  YREG=B0+B1*XREG
  ' CONFIDENCE LIMITS
  TERM=ALPHA+2*BETA*XREG+DELTA*XREG^2
  CONF_LIM=SQRT(TERM)
  LOW_CONF=YREG-T*S*CONF_LIM        ;LOWER LIMIT
  ' OUTPUT
     'REGRESSION
  COL(RES)=XREG           ' X VALUES OF REGRESSION LINE
  COL(RES+1)=YREG         ' Y VALUES OF REGRESSION LINE
     'LOWER CONFIDENCE INTERVAL
  COL(RES+2)=LOW_CONF     ' LOWER CONFIDENCE LIMIT
  'COL(RES+9)=LOW_CONF    ' LOWER CONFIDENCE LIMIT
     'DEGENERATE LINE PLOT FOR SPECIFICATION LIMIT DROP LINES
  CELL(RES+3,1)=IF(N<3, "N MUST > 2", IF(B1>0 OR B24AC<0, "NO SOLUTION", ROOT))
  CELL(RES+4,1)=Y0
     'SHELF LIFE
  CELL(RES+5,1)= "   T90  =  "
  CELL(RES+6,1) = IF(B1>0, " + INFINITY",
  IF(B24AC<0, " NO SOLUTION", ROOT))