ChangeLog

2015-04-03  Ryosuke Niwa  <rniwa@webkit.org>

        Add time series segmentation algorithms as moving averages

        https://bugs.webkit.org/show_bug.cgi?id=143362

        Reviewed by NOBODY (OOPS!).

        This patch implements two preliminary time series segmentation algorithms as moving averages.

        Recursive t-test: Compute Welch's t-statistic at each point in a given segment of the time series.

        If Welch's t-test implicates a statistically significance difference, then split the segment into two

        sub segments with the maximum t-statistic (i.e. the point at which if split would yield the highest

        probability that two segments do not share the same "underlying" mean in classical / frequentist sense).

        We repeat this process recursively. See [1] for the evaluation of this particular algorithm.

        Schwarz criterion: Use Schwarz or Bayesian information criterion to heuristically find the optimal

        segmentation. Intuitively, the problem of finding the best segmentation comes down to minimizing the

        residual sum of squares in each segment as in linear regressions. That is, for a given segment with

        values y_1 through y_n with mean y_avg, we want to minimize the sum of (y_i - y_avg)^2 over i = 1

        through i = n. However, we also don't want to split every data point into a separate segment so we need

        to account the "cost" of introducing new segments. We use a cost function that's loosely based on two

        models discussed in [2] for simplicity. We will tune this cost function further in the future.

        The problem of finding the best segmentation then reduces to a search problem. Unfortunately, our problem

        space is exponential with respect to the size of the time series since we could split at each data point.

        We workaround this problem by first splitting the time series into a manageable smaller grids, and only

        considering segmentation of a fixed size (i.e. the number of segments is constant). Since time series

        tend to contain a lot more data points than segments, this strategy finds the optimal solution without

        exploring much of the problem space.

        Finding the optimal segmentation of a fixed size is, itself, another search problem that is equivalent to

        finding the shortest path of a fixed length in DAG. Here, we use dynamic programming with a matrix of size

        n by n where n is the length of the time series (grid). Each entry in this matrix at (i, k) stores

        the minimum cost of segmenting data points 1 through i using k segments. We start our search at i = 0.

        Clearly C(1, 0) = 0 (note the actual code uses 0-based index). In i-th iteration, we compute the cost

        S(i, j) of each segment starting at i and ending at another point j after i and update C(j, k + 1) by

        min( C(j, k + 1), C(i, k) + S(i, j) ) for all values of j above i.

        [1] Kensuke Fukuda, H. Eugene Stanley, and Luis A. Nunes Amaral, "Heuristic segmentation of

        a nonstationary time series", Physical Review E 69, 021108 (2004)

        [2] Marc Lavielle, Gilles Teyssi`ere, "Detection of Multiple Change–Points in Multivariate Time Series"

        Lithuanian Mathematical Journal, vol 46, 2006

        * public/v2/index.html: Show the optional description for the chosen moving average strategy.

        * public/v2/js/statistics.js:

        (Statistics.testWelchsT):

        (Statistics.computeWelchsT): Extracted from testWelchsT. Generalized to take the offset and the length

        of each value array between which Welch's t-statistic is computed. This generalization helps the

        Schwarz criterion segmentation algorithm avoid splitting values array O(n^2) times.

        (.sampleMeanAndVarianceForValues): Ditto for the generalization.

        (.recursivelySplitIntoTwoSegmentsAtMaxTIfSignificantlyDifferent): Added. Implements recursive t-test.

        (.splitIntoSegmentsUntilGoodEnough): Added. Implements Schwarz criterion.

        (.findOptimalSegmentation): Added. Implements the algorithm to find the optimal segmentation of a fixed

        segment count.

        (.SampleVarianceUpperTriangularMatrix): Added. Stores S(i, j) used by findOptimalSegmentation.

        (.SampleVarianceUpperTriangularMatrix.prototype.costBetween): Added.

2015-04-02  Ryosuke Niwa  <rniwa@webkit.org>

        Perf dashboard should have UI to test out anomaly detection strategies

182308

Websites/perf.webkit.org/public/v2/js/statistics.js

@@var Statistics = new (function () {

    // Welch's t-test (http://en.wikipedia.org/wiki/Welch%27s_t_test)

    this.testWelchsT = function (values1, values2, probability) {

61 var stat1 = sampleMeanAndVarianceForValues(values1);

62 var stat2 = sampleMeanAndVarianceForValues(values2);

        return this.computeWelchsT(values1, 0, values1.length, values2, 0, values2.length, probability).significantlyDifferent;

    }

    this.computeWelchsT = function (values1, startIndex1, length1, values2, startIndex2, length2, probability) {

        var stat1 = sampleMeanAndVarianceForValues(values1, startIndex1, length1);

        var stat2 = sampleMeanAndVarianceForValues(values2, startIndex2, length2);

6367 var sumOfSampleVarianceOverSampleSize = stat1.variance / stat1.size + stat2.variance / stat2.size;

64 var t = (stat1.mean - stat2.mean) / Math.sqrt(sumOfSampleVarianceOverSampleSize);

68 var t = Math.abs((stat1.mean - stat2.mean) / Math.sqrt(sumOfSampleVarianceOverSampleSize));

        // http://en.wikipedia.org/wiki/Welch–Satterthwaite_equation

        var degreesOfFreedom = sumOfSampleVarianceOverSampleSize * sumOfSampleVarianceOverSampleSize

            / (stat1.variance * stat1.variance / stat1.size / stat1.size / stat1.degreesOfFreedom

                + stat2.variance * stat2.variance / stat2.size / stat2.size / stat2.degreesOfFreedom);

        // They're different beyond the confidence interval of the specified probability.

        return Math.abs(t) > tDistributionQuantiles[probability || 0.9][Math.round(degreesOfFreedom - 1)];

        return {

            t: t,

            degreesOfFreedom: degreesOfFreedom,

            significantlyDifferent: t > tDistributionQuantiles[probability || 0.9][Math.round(degreesOfFreedom - 1)],

        };

7379 }

7480

    function sampleMeanAndVarianceForValues(values) {

        var sum = Statistics.sum(values);

        var squareSum = Statistics.squareSum(values);

        var sampleMean = sum / values.length;

    function sampleMeanAndVarianceForValues(values, startIndex, length) {

        var sum = 0;

        for (var i = 0; i < length; i++)

            sum += values[startIndex + i];

        var squareSum = 0;

        for (var i = 0; i < length; i++)

            squareSum += values[startIndex + i] * values[startIndex + i];

        var sampleMean = sum / length;

7989 // FIXME: Maybe we should be using the biased sample variance.

80 var unbiasedSampleVariance = (squareSum - sum * sum / ~~values.~~length) / (~~values.~~length - 1);

90 var unbiasedSampleVariance = (squareSum - sum * sum / length) / (length - 1);

        return {

            mean: sampleMean,

            variance: unbiasedSampleVariance,

84 size: values.length,

85 degreesOfFreedom: values.length - 1,

            size: length,

            degreesOfFreedom: length - 1,

        }

    }

    function recursivelySplitIntoTwoSegmentsAtMaxTIfSignificantlyDifferent(values, startIndex, length, minLength, segments) {

        var tMax = 0;

        var argTMax = null;

        for (var i = 1; i < length - 1; i++) {

            var firstLength = i;

            var secondLength = length - i;

            if (firstLength < minLength || secondLength < minLength)

                continue;

            var result = Statistics.computeWelchsT(values, startIndex, firstLength, values, startIndex + i, secondLength, 0.9);

            if (result.significantlyDifferent && result.t > tMax) {

                tMax = result.t;

                argTMax = i;

            }

        }

        if (!tMax) {

            segments.push(values.slice(startIndex, startIndex + length));

            return;

86116 }

117 recursivelySplitIntoTwoSegmentsAtMaxTIfSignificantlyDifferent(values, startIndex, argTMax, minLength, segments);

118 recursivelySplitIntoTwoSegmentsAtMaxTIfSignificantlyDifferent(values, startIndex + argTMax, length - argTMax, minLength, segments);

    }

    // One-sided t-distribution.

@@var Statistics = new (function () {

                return averages;

            }

        },

        {

            id: 4,

            label: 'Segmentation: Recursive t-test',

            description: "Recursively split values into two segments if Welch's t-test detects a statistically significant difference.",

            parameterList: [{label: "Minimum segment length", value: 20, min: 5}],

            execute: function (minLength, values) {

                if (values.length < 2)

                    return null;

                var averages = new Array(values.length);

                var segments = new Array;

                recursivelySplitIntoTwoSegmentsAtMaxTIfSignificantlyDifferent(values, 0, values.length, minLength, segments);

                var averageIndex = 0;

                for (var j = 0; j < segments.length; j++) {

                    var values = segments[j];

                    var mean = Statistics.sum(values) / values.length;

                    for (var i = 0; i < values.length; i++)

                        averages[averageIndex++] = mean;

                }

                return averages;

            }

        },

        {

            id: 5,

            label: 'Segmentation: Schwarz criterion',

            description: 'Adoptive algorithnm that maximizes the Schwarz criterion (BIC).',

            // Based on Detection of Multiple Change–Points in Multivariate Time Series by Marc Lavielle (July 2006).

            execute: function (values) {

                if (values.length < 2)

                    return null;

                var averages = new Array(values.length);

                var averageIndex = 0;

                // Split the time series into grids since splitIntoSegmentsUntilGoodEnough is O(n^2).

                var gridLength = 500;

                var totalSegmentation = [0];

                for (var gridCount = 0; gridCount < Math.ceil(values.length / gridLength); gridCount++) {

                    var gridValues = values.slice(gridCount * gridLength, (gridCount + 1) * gridLength);

                    var segmentation = splitIntoSegmentsUntilGoodEnough(gridValues);

                    if (Statistics.debuggingSegmentation)

                        console.log('grid=' + gridCount, segmentation);

                    for (var i = 1; i < segmentation.length - 1; i++)

                        totalSegmentation.push(gridCount * gridLength + segmentation[i]);

                }

                if (Statistics.debuggingSegmentation)

                    console.log('Final Segmentation', totalSegmentation);

                totalSegmentation.push(values.length);

                for (var i = 1; i < totalSegmentation.length; i++) {

                    var segment = values.slice(totalSegmentation[i - 1], totalSegmentation[i]);

                    var average = Statistics.sum(segment) / segment.length;

                    for (var j = 0; j < segment.length; j++)

                        averages[averageIndex++] = average;

                }

                return averages;

            }

        },

196292 ];

197293

    this.debuggingSegmentation = false;

    function splitIntoSegmentsUntilGoodEnough(values) {

        if (values.length < 2)

            return [0, values.length];

        var matrix = new SampleVarianceUpperTriangularMatrix(values);

        var SchwarzCriterionBeta = Math.log1p(values.length - 1) / values.length;

        var BirgeAndMassartC = 2.5; // Suggested by the authors.

        var BirgeAndMassartPenalization = function (segmentCount) {

            return segmentCount * (1 + BirgeAndMassartC * Math.log1p(values.length / segmentCount - 1));

        }

        var segmentation;

        var minTotalCost = Infinity;

        var maxK = 50;

        for (var k = 1; k < maxK; k++) {

            var start = Date.now();

            var result = findOptimalSegmentation(values, matrix, k);

            var cost = result.cost / values.length;

            var penalty = SchwarzCriterionBeta * BirgeAndMassartPenalization(k);

            if (cost + penalty < minTotalCost) {

                minTotalCost = cost + penalty;

                segmentation = result.segmentation;

            } else

                maxK = Math.min(maxK, k + 3);

            if (Statistics.debuggingSegmentation)

                console.log('splitIntoSegmentsUntilGoodEnough', k, Date.now() - start, cost + penalty);

        }

        return segmentation;

    }

    function findOptimalSegmentation(values, costMatrix, segmentCount) {

        // Dynamic programming. cost[i][k] = The cost to segmenting values up to i into k segments.

        var cost = new Array(values.length);

        for (var i = 0; i < values.length; i++) {

            cost[i] = new Float32Array(segmentCount + 1);

        }

        var previousNode = new Array(values.length);

        for (var i = 0; i < values.length; i++)

            previousNode[i] = new Array(segmentCount + 1);

        cost[0] = [0]; // The cost of segmenting single value is always 0.

        previousNode[0] = [-1];

        for (var segmentStart = 0; segmentStart < values.length; segmentStart++) {

            var costBySegment = cost[segmentStart];

            for (var count = 0; count < segmentCount; count++) {

                if (previousNode[segmentStart][count] === undefined)

                    continue;

                for (var segmentEnd = segmentStart + 1; segmentEnd < values.length; segmentEnd++) {

                    var newCost = costBySegment[count] + costMatrix.costBetween(segmentStart, segmentEnd);

                    if (previousNode[segmentEnd][count + 1] === undefined || newCost < cost[segmentEnd][count + 1]) {

                        cost[segmentEnd][count + 1] = newCost;

                        previousNode[segmentEnd][count + 1] = segmentStart;

                    }

                }

            }

        }

        if (Statistics.debuggingSegmentation) {

            console.log('findOptimalSegmentation with k=', segmentCount);

            for (var i = 0; i < cost.length; i++) {

                var t = cost[i];

                var s = '';

                for (var j = 0; j < t.length; j++) {

                    var p = previousNode[i][j];

                    s += '(k=' + j;

                    if (p !== undefined)

                        s += ' c=' + t[j] + ' p=' + p

                    s += ')';

                }

                console.log(i, values[i], s);

            }

        }

        var currentIndex = values.length - 1;

        var segmentation = new Array(segmentCount);

        segmentation[0] = values.length;

        for (var i = 0; i < segmentCount; i++) {

            currentIndex = previousNode[currentIndex][segmentCount - i];

            segmentation[i + 1] = currentIndex;

        }

        return {segmentation: segmentation.reverse(), cost: cost[values.length - 1][segmentCount]};

    }

    function SampleVarianceUpperTriangularMatrix(values) {

        // The cost of segment (i, j].

        var costMatrix = new Array(values.length - 1);

        for (var i = 0; i < values.length - 1; i++) {

            var remainingValueCount = values.length - i - 1;

            costMatrix[i] = new Float32Array(remainingValueCount);

            var sum = values[i];

            var squareSum = sum * sum;

            costMatrix[i][0] = 0;

            for (var j = i + 1; j < values.length; j++) {

                var currentValue = values[j];

                sum += currentValue;

                squareSum += currentValue * currentValue;

                var sampleSize = j - i + 1;

                var stdev = Statistics.sampleStandardDeviation(sampleSize, sum, squareSum);

                costMatrix[i][j - i - 1] = sampleSize * Math.log1p(stdev * stdev - 1);

            }

        }

        this.costMatrix = costMatrix;

    }

    SampleVarianceUpperTriangularMatrix.prototype.costBetween = function(from, to) {

        if (from >= this.costMatrix.length || from == to)

            return 0; // The cost of the segment that starts at the last data point is 0.

        return this.costMatrix[from][to - from - 1];

    }

    this.EnvelopingStrategies = [

        {

            id: 100,

182308

Websites/perf.webkit.org/ChangeLog

Websites/perf.webkit.org/public/v2/index.html

Websites/perf.webkit.org/public/v2/js/statistics.js