Holt forecasting

holt.juttle

// Holt forecaster, which models an unobserved
// level and trend component in a noisy signal. The result is a smoothed
// version of the signal which may be used to forecast future values.
//
// Y is the point field to smooth.
// slevel and strend are smoothing factors, numbers ranging [0..1]. 
// They determine how quickly the smoother adjusts to new values as they arrive.
// Setting a factor to 0 freezes the feature at its initial estimate,
// while setting a factor to 1 causes that feature to depend solely
// on the most recent point. Setting strend to null removes that feature
// from the model.
//
// smooth assumes arriving points are equally spaced in time.
// smooth returns the next-step forecast of the model, an array containing:
// [level, trend]
// level is the estimate of the next series value, trend the estimate of the
//  next series slope.
//
export reducer smooth(Y, slevel, strend) {
    var level, trend, level_1, trend_1;
    var initial = true;
     
    function update() {
        if (initial) {
            level = *Y;
            trend = 0;
            level_1 = *Y;
            trend_1 = 0;
            initial = false;
        } else {
            level = slevel * *Y + (1 - slevel) * (level_1 + trend_1);
            if (strend != null) {
                trend = strend * (level - level_1) + (1 - strend) * trend_1;
            }
            level_1 = level;
            trend_1 = trend;
        }
    }
    function result() {
        return [level +  trend, trend];
    }
    function reset() {
        initial = true;
    }
}

// smoother/estimator for standard deviation of a time series
// alpha [0..1] is the smoothing factor.
export reducer smooth_sd(Y, alpha) {
    var sigma2, last;
    var initial = true;
 
    function update() {
        if (initial) {
            sigma2 = 0;
            last = *Y;
            initial = false;
        } else {
            sigma2 = alpha * (*Y - last) * (*Y - last) + (1 - alpha) * sigma2;
            last = *Y;
        }
    }
    function result() {
        return Math.sqrt(sigma2);
    }
    function reset() {
        initial = true;
    }
}

// Perform a Holt forecast along with prediction intervals based
// on a smoothed stderr estimate. Y is the field to forecast. z
// is the number of stdandard deviations from mean to include in
// the prediction interval (z = 1.96 is the 97.5 percentile and thus
// a 95% confidence bound, if errors were idependent and normally
// distributed (often untrue)).
// slevel and strend are smoothing factors for Holt.
//
// For a field name Y, this populates the current point with:
// Y_pred: predicted next y value
// Y_trend: predicted next slope
// Y_err: estimated stderr (RMS prediction error over time)
// Y_hi, Y_low: prediction interval bounds
// 
export sub forecastPI(Y, z, alpha, slevel, strend) {
    put state = smooth(Y, slevel, strend)
    | put *(Y+"_pred") = state[0], *(Y+"_trend")=state[1]
    | put err2 = (*Y - state[0]) * (*Y - state[0])
    | put stderr = Math.sqrt(smooth(err2, alpha, null)[0])
    | put *(Y+"_err") = stderr
    | put *(Y+"_hi") = *(Y+"_pred") + z * stderr
    | put *(Y+"_low") = *(Y+"_pred") - z * stderr
    | remove state, err2, stderr
}

// Data generation utilities for demos:

sub hourly_cpu(from, to) { 
  demo cdn metrics 'cpu' -from from -to to -every :m:
    -nhosts 4 -dos .5 -dos_dur :5m: -ripple .3 -cpu_alpha 0.8
  | filter host ~ 'sjc*'
  | put cpu=value
  | remove type, host, pop, value
}
 
sub daily_cpu(from, to) { 
  demo cdn metrics 'cpu' -from from -to to -every :10m:
    -nhosts 4 -dos .3 -dos_dur :15m: -cpu_alpha 0.8 -daily 4
  | filter host ~ 'sjc*'
  | put cpu=value
  | remove type, host, pop, value
}

sub bumpy_cpu(from, to) { 
  demo cdn metrics 'cpu' -from from -to to -every :5m:
    -nhosts 4 -dos .3 -dos_dur :15s: -daily 3 
  | filter host ~ 'sjc*' 
  | put cpu=value
  | remove type, host, pop, value
}

// Demo 1: exponential forecasting for simulated CPU
//
// This demo applies a simple exponential smoother (the "level" portion of
// smooth() without trend estimation). It displays timeseries of
// the original input, the predicted/smoothed input, a 95% confidence
// interval around this prediction, and the stderr. 
// 
bumpy_cpu -from :2014-01-01T13:00:00:  -to :2014-01-01T21:00:00:
| forecastPI -Y "cpu" -z 1.96 -alpha .8 -slevel .5 -strend null
| remove cpu_trend
| split cpu, cpu_pred, cpu_low, cpu_hi, cpu_err
| @timechart;

//cpu -from :2014-01-01: -to :2014-01-01T01:20:00:
//daily_cpu -from :2014-01-01: -to :2014-01-02T01:20:00:

// Demo 2: Holt forecasting for simulated CPU
//
// This demo performs a Holt forecast. It displays timeseries of
// the original input, the predicted/smoothed input, a 95% confidence
// interval around this prediction, and the stderr. It also displays
// the predicted trend, which turns out to be useful at small time
// horizons for CPU anomaly detection (its upswings react to a cpu
// on the move immediately, whereas stderr lags as the graph shows).
// At longer horizons, CPUs do not have trends, just variation
// over the 0-100% range.
//
bumpy_cpu -from :2014-01-01T13:00:00:  -to :2014-01-01T21:00:00:
| forecastPI -Y "cpu" -z 1.96 -alpha .6 -slevel .6 -strend .6
| split cpu, cpu_pred, cpu_trend, cpu_low, cpu_hi, cpu_err
| @timechart;