TRANSFORM

' 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(R0 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",