**FREE
CTL-OPT OPTION(*NOSRCSTMT:*NODEBUGIO:*NOSHOWCPY) NOMAIN BNDDIR('STATS');

/COPY *LIBL/STATS,STATSCOPY

//******************************************************************************
// standalone fields & constants
//******************************************************************************

DCL-C @MaxElements 500;
DCL-S Real# FLOAT(8);
DCL-S Integer INT(10);
DCL-S Index INT(5);

//******************************************************************************
// procedure: least squares line equation
//******************************************************************************

DCL-PROC $LeastSquaresLineEq EXPORT;
DCL-PI *N LIKEDS(ds$LeastSquaresLineEq);
##DataPts# LIKE(Index) CONST;
##Y LIKE(Integer) CONST DIM(@MaxElements);
##X LIKE(Integer) CONST DIM(@MaxElements) OPTIONS(*NOPASS);
END-PI;

DCL-S i Like(Index);
DCL-S SumX LIKE(##X);
DCL-S AvgX LIKE(##X);
DCL-S SumY LIKE(##Y);
DCL-S AvgY LIKE(##Y);
DCL-S Numerator LIKE(Real#);
DCL-S Denominator LIKE(Real#);

CLEAR ds$LeastSquaresLineEq;

MONITOR;
   FOR i = 1 TO ##DataPts#;
      SumX += ##X(i);
      SumY += ##Y(i);
   ENDFOR;

   AvgX = SumX / ##DataPts#;
   AvgY = SumY / ##DataPts#;

   FOR i = 1 TO ##DataPts#;
      Numerator += (##X(i) - AvgX) * (##Y(i) - AvgY);
      Denominator += (##X(i) - AvgX) ** 2;
   ENDFOR;

   ds$LeastSquaresLineEq.Slope = Numerator / Denominator;
   ds$LeastSquaresLineEq.YIntercept =
      AvgY - ds$LeastSquaresLineEq.Slope * AvgX;
ON-ERROR;
   ds$LeastSquaresLineEq.Error = *ON;
ENDMON;

RETURN ds$LeastSquaresLineEq;

END-PROC;

//******************************************************************************
// procedure: least squares mean average deviation
//******************************************************************************

DCL-PROC $LeastSquaresMAD EXPORT;
DCL-PI *N LIKEDS(ds$LeastSquaresMAD);
##DataPts# LIKE(Index) CONST;
##Duration LIKE(Index) CONST;
##Y LIKE(Integer) CONST DIM(@MaxElements);
##X LIKE(Integer) CONST DIM(@MaxElements);
END-PI;

DCL-S @Error IND INZ;
DCL-S Slope Like(Real#);
DCL-S YIntercept Like(Real#);
DCL-S p Like(Index);

CLEAR ds$LeastSquaresMAD;

MONITOR;
   FOR p = ##Duration + 1 TO ##DataPts#;
      ds$LeastSquaresLineEq =
         $LeastSquaresLineEq(
            ##Duration :
            %SubArr(##Y : p - ##Duration : ##Duration) :
            %SubArr(##X : p - ##Duration : ##Duration));

     ds$LeastSquaresMAD.MeanAvgDev +=
        %ABS((ds$LeastSquaresLineEq.YIntercept +
              ds$LeastSquaresLineEq.Slope * p) - ##Y(p));
   ENDFOR;

   ds$LeastSquaresMAD.MeanAvgDev =
      ds$LeastSquaresMAD.MeanAvgDev / (##DataPts# - 2);

ON-ERROR;
   ds$LeastSquaresLineEq.Error = *ON;
ENDMON;

RETURN ds$LeastSquaresMAD;

END-PROC;

//******************************************************************************
// procedure: winters forecast cycle
//******************************************************************************

DCL-PROC $WintersFcstCycle EXPORT;
DCL-PI *N LIKEDS(ds$WintersFcstCycle);
##CycleSize LIKE(Index) CONST;
##Alpha PACKED(4:3) CONST;
##Beta PACKED(4:3) CONST;
##Gamma PACKED(4:3) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S Accumulator Like(Real#) INZ;
DCL-S i Like(Index);
DCL-S Seasonality Like(Real#) INZ DIM(@MaxElements);
DCL-S SmoothedEst Like(Real#) INZ DIM(@MaxElements);
DCL-S TrendQty Like(Real#) INZ DIM(@MaxElements);

CLEAR ds$WintersFcstCycle;

MONITOR;
   FOR i = 1 TO ##CycleSize;
      Accumulator += ##Y(i);
   ENDFOR;

   FOR i = 1 TO ##CycleSize;
      Seasonality(i) = ##CycleSize * (##Y(i) / Accumulator);
   ENDFOR;

   SmoothedEst(##CycleSize) = Accumulator / ##CycleSize;

   CLEAR TrendQty(##CycleSize);

   FOR i = ##CycleSize + 1 TO ##CycleSize * 2;

      SmoothedEst(i) =
         (##Alpha * ##Y(i) / Seasonality(i-##CycleSize)) +
         (1 - ##Alpha) * (SmoothedEst(i-1) + TrendQty(i-1));

      TrendQty(i) =
         ##Beta * (SmoothedEst(i) - SmoothedEst(i-1)) +
         (1 - ##Beta) * TrendQty(i-1);

      Seasonality(i) =
         ##Gamma * ##Y(i) / SmoothedEst(i) +
         (1 - ##Alpha) * Seasonality(i - ##CycleSize);

      ds$WintersFcstCycle.Forecast(i-##CycleSize) =
         (SmoothedEst(i-1) + TrendQty(i-1)) *
          Seasonality(i-##CycleSize);
   ENDFOR;
ON-ERROR;
   ds$WintersFcstCycle.Error = *ON;
ENDMON;

RETURN ds$WintersFcstCycle;

END-PROC;

//******************************************************************************
// procedure: winters mean average deviation
//******************************************************************************

DCL-PROC $WintersMAD EXPORT;
DCL-PI *N LIKEDS(ds$WintersMAD);
##CyclePts# LIKE(Index) CONST;
##Alpha PACKED(4:3) CONST;
##Beta PACKED(4:3) CONST;
##Gamma PACKED(4:3) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S Forecast Like(Real#) DIM(@MaxElements);
DCL-S i Like(Index);

CLEAR ds$WintersMAD;

MONITOR;
   ds$WintersMAD.Error =
      $WintersFcstCycle(##CyclePts# :
                        ##Alpha :
                        ##Beta :
                        ##Gamma :
                        ##Y);
   IF ds$WintersMAD.Error = *OFF;
      FOR i = 1 TO ##CyclePts#;
         ds$WintersMAD.MeanAvgDev +=
            %ABS(Forecast(i) - ##Y(i + ##CyclePts#));
      ENDFOR;
      ds$WintersMAD.MeanAvgDev = ds$WintersMAD.MeanAvgDev / ##CyclePts#;
   ENDIF;
ON-ERROR;
   ds$WintersMAD.Error = *ON;
ENDMON;

RETURN ds$WintersMAD;

END-PROC;

//******************************************************************************
// procedure: moving average forecast
//******************************************************************************

DCL-PROC $MovingAvgFcst EXPORT;
DCL-PI *N LIKEDS(ds$MovingAvgFcst);
##DataPts# LIKE(Index) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S i Like(Index);

CLEAR ds$MovingAvgFcst;

MONITOR;
   FOR i = 1 TO ##DataPts#;
      ds$MovingAvgFcst.Forecast += ##Y(i);
   ENDFOR;
   ds$MovingAvgFcst.Forecast = ds$MovingAvgFcst.Forecast / ##DataPts#;
ON-ERROR;
   ds$MovingAvgFcst.Error = *ON;
ENDMON;

RETURN ds$MovingAvgFcst;

END-PROC;

//******************************************************************************
// procedure: moving average mean average deviation
//******************************************************************************

DCL-PROC $MovingAvgMAD EXPORT;
DCL-PI *N LIKEDS(ds$MovingAvgMAD);
##DataPts# LIKE(Index) CONST;
##Duration LIKE(INTEGER) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S i Like(Index);

CLEAR ds$MovingAvgMAD;

MONITOR;
   FOR i = 1 TO ##DataPts# - ##Duration;
      ds$MovingAvgFcst =
         $MovingAvgFcst(##Duration :
                        %SubArr(##Y : i : ##Duration));
     IF ds$MovingAvgFcst.Error;
        ds$MovingAvgMAD.Error = *ON;
        LEAVE;
     ENDIF;

     ds$MovingAvgMAD.MeanAvgDev +=
        %ABS(ds$MovingAvgFcst.Forecast - ##Y(i+##Duration));
   ENDFOR;
   ds$MovingAvgMAD.MeanAvgDev =
      ds$MovingAvgMAD.MeanAvgDev / (##DataPts# - ##Duration);
ON-ERROR;
   ds$MovingAvgMAD.Error = *ON;
ENDMON;

RETURN ds$MovingAvgMAD;

END-PROC;

//******************************************************************************
// procedure: exponential smoothing linear trend forecast
//******************************************************************************

DCL-PROC $ExpSmthngLinTrndFcst EXPORT;
DCL-PI *N LIKEDS(ds$ExpSmthngLinTrndFcst);
##DataPts# LIKE(Index) CONST;
##Alpha PACKED(4:3) CONST;
##Beta PACKED(4:3) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S SmoothedEst Like(Real#) INZ DIM(@MaxElements);
DCL-S TrendQty Like(Real#) INZ DIM(@MaxElements);

DCL-S i Like(Index);

CLEAR ds$ExpSmthngLinTrndFcst;

MONITOR;
   SmoothedEst(1) = ##Y(1);
   FOR i = 2 TO ##DataPts#;
      SmoothedEst(i) =
         (##Alpha * ##Y(i)) +
         (1 - ##Alpha) * (SmoothedEst(i-1) + TrendQty(i-1));

      TrendQty(i) =
         ##Beta * (SmoothedEst(i) - SmoothedEst(i-1)) +
         (1 - ##Beta) * (TrendQty(i-1));

      ds$ExpSmthngLinTrndFcst.Forecast(i) =
         SmoothedEst(i-1) + TrendQty(i-1);
   ENDFOR;
ON-ERROR;
   ds$ExpSmthngLinTrndFcst.Error = *ON;
ENDMON;

RETURN ds$ExpSmthngLinTrndFcst;

END-PROC;

//******************************************************************************
// procedure: exponential smoothing linear trend mean avg deviation
//******************************************************************************

DCL-PROC $ExpSmthngLinTrndMAD EXPORT;
DCL-PI *N LIKEDS(ds$ExpSmthngLinTrndMAD);
##DataPts# LIKE(Index) CONST;
##Alpha LIKE(Real#) CONST;
##Beta LIKE(Real#) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S i Like(Index);

CLEAR ds$ExpSmthngLinTrndMAD;

MONITOR;
   ds$ExpSmthngLinTrndFcst =
      $ExpSmthngLinTrndFcst(##DataPts# :
                            ##Alpha :
                            ##Beta :
                            ##Y);

   IF ds$ExpSmthngLinTrndFcst.Error;
      ds$ExpSmthngLinTrndMAD.Error = *ON;
   ELSE;
      FOR i = 2 TO ##DataPts#;
         ds$ExpSmthngLinTrndMAD.MeanAvgDev +=
            %ABS(ds$ExpSmthngLinTrndFcst.Forecast(i) - ##Y(i));
      ENDFOR;
      ds$ExpSmthngLinTrndMAD.MeanAvgDev =
         ds$ExpSmthngLinTrndMAD.MeanAvgDev / (##DataPts# - 1);
   ENDIF;
ON-ERROR;
   ds$ExpSmthngLinTrndMAD.Error = *ON;
ENDMON;

RETURN ds$ExpSmthngLinTrndMAD;

END-PROC;

//******************************************************************************
// procedure: exponential smoothing forecast
//******************************************************************************

DCL-PROC $ExpSmthngFcst EXPORT;
DCL-PI *N LIKEDS(ds$ExpSmthngFcst);
##DataPts# LIKE(Index) CONST;
##Alpha PACKED(4:3) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S i Like(Index);

CLEAR ds$ExpSmthngFcst;

MONITOR;
   ds$ExpSmthngFcst.Forecast(2) = ##Y(1);
   FOR i = 3 TO ##DataPts#;
      ds$ExpSmthngFcst.Forecast(i) =
         (##Alpha * ##Y(i-1)) + ((1 - ##Alpha) *
          ds$ExpSmthngFcst.Forecast(i-1));
   ENDFOR;
ON-ERROR;
   ds$ExpSmthngFcst.Error = *ON;
ENDMON;

RETURN ds$ExpSmthngFcst;

END-PROC;

//******************************************************************************
// procedure: exponential smoothing mean avg deviation
//******************************************************************************

DCL-PROC $ExpSmthngMAD EXPORT;
DCL-PI *N LIKEDS(ds$ExpSmthngMAD);
##DataPts# LIKE(Index) CONST;
##Alpha LIKE(Real#) CONST;
##Y LIKE(INTEGER) DIM(@MaxElements) CONST;
END-PI;

DCL-S i Like(Index);
DCL-S Forecast Like(Real#) DIM(@MaxElements);

CLEAR ds$ExpSmthngMAD;

MONITOR;
   ds$ExpSmthngFcst =
      $ExpSmthngFcst(##DataPts# :
                     ##Alpha :
                     ##Y);

   IF ds$ExpSmthngFcst.Error;
      ds$ExpSmthngMAD.Error = *ON;
   ELSE;
      FOR i = 2 TO ##DataPts#;
         ds$ExpSmthngMAD.MeanAvgDev += %ABS(Forecast(i) - ##Y(i));
      ENDFOR;
      ds$ExpSmthngMAD.MeanAvgDev =
         ds$ExpSmthngMAD.MeanAvgDev / (##DataPts# - 1);
   ENDIF;
ON-ERROR;
   ds$ExpSmthngMAD.Error = *ON;
ENDMON;

RETURN ds$ExpSmthngMAD;

END-PROC;