Skip to content

Commit

Permalink
Fix the two piecewise-linear regression calculation
Browse files Browse the repository at this point in the history 
#66

The current implementation of the regression calculation has these flaws:

1. When processing (x[0], y[0]), L1 must be any line through (x0, y0) which meets
   L2 at a point (x’, y’) where x[0] < x' < x[1]. L1 has no error.
2. When processing (x[n - 2], y[n - 2]), L2 must be any line through
   (x[n - 1], y[n - 1]) which meets L1 at a point (x’, y’) where
   x[n - 2] < x' < x[n - 1]. L2 has no error.
3. The lambda calculation is incorrect. It includes a term called H which is equal
   to C - I. Looking at the algorithm of Kundu/Ubhaya, this should be just C.
4. lambda should to be used with calculating L1 and (1 - lambda) should to be used
   with calculating L2. Currently (1 - lambda) is used in calculating L1 and L2.
5. The current calculation has this condition if (t1 != t2) continue; This
   condition is almost always true even if t1 and t2 are essentiallyEqual.
  • Loading branch information
@shallawa
shallawa committed Oct 10, 2024
1 parent f1c7edb commit 40af19e
Showing 1 changed file with 187 additions and 147 deletions.
334 changes: 187 additions & 147 deletions MotionMark/resources/statistics.js
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
/*
* Copyright (C) 2015-2017 Apple Inc. All rights reserved.
* Copyright (C) 2015-2024 Apple Inc. All rights reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions
Expand Down Expand Up @@ -181,39 +181,39 @@ Regression = Utilities.createClass(
// All samples are analyzed. startIndex, endIndex are just stored for use by the caller.
function(samples, startIndex, endIndex, options)
{
const desiredFrameLength = options.desiredFrameLength;
var profile;
this.startIndex = Math.min(startIndex, endIndex);
this.endIndex = Math.max(startIndex, endIndex);

if (!options.preferredProfile || options.preferredProfile == Strings.json.profiles.slope) {
profile = this._calculateRegression(samples, {
shouldClip: true,
s1: desiredFrameLength,
t1: 0
});
this.profile = Strings.json.profiles.slope;
} else if (options.preferredProfile == Strings.json.profiles.flat) {
profile = this._calculateRegression(samples, {
shouldClip: true,
s1: desiredFrameLength,
this.s1 = 0;
this.t1 = 0;
this.n1 = 0;
this.e1 = Number.MAX_VALUE;

this.s2 = 0;
this.t2 = 0;
this.n2 = 0;
this.e2 = Number.MAX_VALUE;

this.complexity = 0;

if (options.preferredProfile == Strings.json.profiles.flat) {
this._calculateRegression(samples, {
s1: options.desiredFrameLength,
t1: 0,
t2: 0
});
this.profile = Strings.json.profiles.flat;
} else {
this._calculateRegression(samples, {
s1: options.desiredFrameLength,
t1: 0
});
this.profile = Strings.json.profiles.slope;
}

this.startIndex = Math.min(startIndex, endIndex);
this.endIndex = Math.max(startIndex, endIndex);

this.complexity = profile.complexity;
this.s1 = profile.s1;
this.t1 = profile.t1;
this.s2 = profile.s2;
this.t2 = profile.t2;
this.stdev1 = profile.stdev1;
this.stdev2 = profile.stdev2;
this.n1 = profile.n1;
this.n2 = profile.n2;
this.error = profile.error;
this.stdev1 = Math.sqrt(this.e1 / this.n1);
this.stdev2 = Math.sqrt(this.e2 / this.n2);
this.error = this._error();
}, {

valueAt: function(complexity)
Expand All @@ -223,6 +223,60 @@ Regression = Utilities.createClass(
return this.s1 + this.t1 * complexity;
},

_intersection: function(segment1, segment2)
{
return (segment1.s - segment2.s) / (segment2.t - segment1.t);
},

_error: function() {
return this.e1 + this.e2;
},

_areEssentiallyEqual: function(n1, n2) {
// Choose epsilon not too small to ensure the intersetion
// of the two segements is not too far from sampled data.
const epsilon = 0.0001;
return Math.abs(n1 - n2) < epsilon;
},

_setOptimal: function(segment1, segment2, x, xn, options) {
if (segment1.e + segment2.e > this.e1 + this.e2)
return false;

segment1.s = options.s1 !== undefined ? options.s1 : segment1.s;
segment1.t = options.t1 !== undefined ? options.t1 : segment1.t;
segment2.s = options.s2 !== undefined ? options.s2 : segment2.s;
segment2.t = options.t2 !== undefined ? options.t2 : segment2.t;

// The score is the x coordinate of the intersection of segment1 and segment2.
let complexity = this._intersection(segment1, segment2);

if (!this._areEssentiallyEqual(segment1.t, segment2.t)) {
// If segment1 and segment2 are not parallel, then they have to meet
// at complexity such that x <= complexity <= xn.
if (!(complexity >= x && complexity <= xn))
return false;
} else {
// If segment1 and segment2 are parallel, then they have to form one
// single line.
if (!this._areEssentiallyEqual(segment1.s, segment2.s))
return false;
}

this.s1 = segment1.s;
this.t1 = segment1.t;
this.n1 = segment1.n;
this.e1 = segment1.e;

this.s2 = segment2.s;
this.t2 = segment2.t;
this.n2 = segment2.n;
this.e2 = segment2.e;

this.complexity = complexity;
return true;
},

// A generic two-segment piecewise regression calculator. Based on Kundu/Ubhaya
//
// Minimize sum of (y - y')^2
Expand All @@ -239,153 +293,139 @@ Regression = Utilities.createClass(
const complexityIndex = 0;
const frameLengthIndex = 1;

if (samples.length == 1) {
// Only one sample point; we can't calculate any regression.
var x = samples[0][complexityIndex];
return {
complexity: x,
s1: x,
t1: 0,
s2: x,
t2: 0,
error1: 0,
error2: 0
};
}

// Sort by increasing complexity.
var sortedSamples = samples.slice().sort((a, b) => a[complexityIndex] - b[complexityIndex]);

// x is expected to increase in complexity
var lowComplexity = sortedSamples[0][complexityIndex];
var highComplexity = sortedSamples[samples.length - 1][complexityIndex];

var a1 = 0, b1 = 0, c1 = 0, d1 = 0, h1 = 0, k1 = 0;
var a2 = 0, b2 = 0, c2 = 0, d2 = 0, h2 = 0, k2 = 0;

for (var i = 0; i < sortedSamples.length; ++i) {
var x = sortedSamples[i][complexityIndex];
var y = sortedSamples[i][frameLengthIndex];
let sortedSamples = samples.slice().sort((a, b) => a[complexityIndex] - b[complexityIndex]);

let a1 = 0, b1 = 0, c1 = 0, d1 = 0, h1 = 0, k1 = 0;
let a2 = 0, b2 = 0, c2 = 0, d2 = 0, h2 = 0, k2 = 0;

for (let j = 0; j < sortedSamples.length; ++j) {
let x = sortedSamples[j][complexityIndex];
let y = sortedSamples[j][frameLengthIndex];

a2 += 1;
b2 += x;
c2 += x * x;
d2 += y;
h2 += y * x;
h2 += x * y;
k2 += y * y;
}

var s1_best, t1_best, s2_best, t2_best, n1_best, n2_best, error1_best, error2_best, x_best, x_prime;

function setBest(s1, t1, error1, s2, t2, error2, splitIndex, x_prime, x)
{
s1_best = s1;
t1_best = t1;
error1_best = error1;
s2_best = s2;
t2_best = t2;
error2_best = error2;
// Number of samples included in the first segment, inclusive of splitIndex
n1_best = splitIndex + 1;
// Number of samples included in the second segment
n2_best = samples.length - splitIndex - 1;
if (!options.shouldClip || (x_prime >= lowComplexity && x_prime <= highComplexity))
x_best = x_prime;
else {
// Discontinuous piecewise regression
x_best = x;
}
}
for (let j = 0; j < sortedSamples.length - 1; ++j) {
let x = sortedSamples[j][complexityIndex];
let y = sortedSamples[j][frameLengthIndex];
let xx = x * x;
let xy = x * y;
let yy = y * y;

// Iterate from 0 to n - 2, inclusive
for (var i = 0; i < sortedSamples.length - 1; ++i) {
var x = sortedSamples[i][complexityIndex];
var y = sortedSamples[i][frameLengthIndex];
var xx = x * x;
var yx = y * x;
var yy = y * y;
// a1, b1, etc. is sum from 0 to i, inclusive
a1 += 1;
b1 += x;
c1 += xx;
d1 += y;
h1 += yx;
h1 += xy;
k1 += yy;
// a2, b2, etc. is sum from i+1 to sortedSamples.length - 1, inclusive

a2 -= 1;
b2 -= x;
c2 -= xx;
d2 -= y;
h2 -= yx;
h2 -= xy;
k2 -= yy;

var A = c1*d1 - b1*h1;
var B = a1*h1 - b1*d1;
var C = a1*c1 - b1*b1;
var D = c2*d2 - b2*h2;
var E = a2*h2 - b2*d2;
var F = a2*c2 - b2*b2;
var s1 = options.s1 !== undefined ? options.s1 : (A / C);
var t1 = options.t1 !== undefined ? options.t1 : (B / C);
var s2 = options.s2 !== undefined ? options.s2 : (D / F);
var t2 = options.t2 !== undefined ? options.t2 : (E / F);
// Assumes that the two segments meet
var x_prime = (s1 - s2) / (t2 - t1);

var error1 = (k1 + a1*s1*s1 + c1*t1*t1 - 2*d1*s1 - 2*h1*t1 + 2*b1*s1*t1) || Number.MAX_VALUE;
var error2 = (k2 + a2*s2*s2 + c2*t2*t2 - 2*d2*s2 - 2*h2*t2 + 2*b2*s2*t2) || Number.MAX_VALUE;

if (i == 0) {
setBest(s1, t1, error1, s2, t2, error2, i, x_prime, x);
continue;
}
let A = (c1 * d1) - (b1 * h1);
let B = (a1 * h1) - (b1 * d1);
let C = (a1 * c1) - (b1 * b1);
let D = (c2 * d2) - (b2 * h2);
let E = (a2 * h2) - (b2 * d2);
let F = (a2 * c2) - (b2 * b2);

let s1 = A / C;
let t1 = B / C;
let s2 = D / F;
let t2 = E / F;

if (C == 0 || F == 0)
continue;

// Projected point is not between this and the next sample
var nextSampleComplexity = sortedSamples[i + 1][complexityIndex];
if (x_prime > nextSampleComplexity || x_prime < x) {
// Calculate lambda, which divides the weight of this sample between the two lines

// These values remove the influence of this sample
var I = c1 - 2*b1*x + a1*xx;
var H = C - I;
var G = A + B*x - C*y;

var J = D + E*x - F*y;
var K = c2 - 2*b2*x + a2*xx;

var lambda = (G*F + G*K - H*J) / (I*J + G*K);
if (lambda > 0 && lambda < 1) {
var lambda1 = 1 - lambda;
s1 = options.s1 !== undefined ? options.s1 : ((A - lambda1*(-h1*x + d1*xx + c1*y - b1*yx)) / (C - lambda1*I));
t1 = options.t1 !== undefined ? options.t1 : ((B - lambda1*(h1 - d1*x - b1*y + a1*yx)) / (C - lambda1*I));
s2 = options.s2 !== undefined ? options.s2 : ((D + lambda1*(-h2*x + d2*xx + c2*y - b2*yx)) / (F + lambda1*K));
t2 = options.t2 !== undefined ? options.t2 : ((E + lambda1*(h2 - d2*x - b2*y + a2*yx)) / (F + lambda1*K));
x_prime = (s1 - s2) / (t2 - t1);

error1 = ((k1 + a1*s1*s1 + c1*t1*t1 - 2*d1*s1 - 2*h1*t1 + 2*b1*s1*t1) - lambda1 * Math.pow(y - (s1 + t1*x), 2)) || Number.MAX_VALUE;
error2 = ((k2 + a2*s2*s2 + c2*t2*t2 - 2*d2*s2 - 2*h2*t2 + 2*b2*s2*t2) + lambda1 * Math.pow(y - (s2 + t2*x), 2)) || Number.MAX_VALUE;
} else if (t1 != t2)
continue;
let segment1;
let segment2;
let xp = (j == 0) ? 0 : sortedSamples[j - 1][complexityIndex];

if (j == 0) {
// Let segment1 be any line through (x[0], y[0]) which meets segment2 at
// a point (x’, y’) where x[0] < x' < x[1]. segment1 has no error.
let xMid = (x + sortedSamples[j + 1][complexityIndex]) / 2;
let yMid = s2 + t2 * xMid;
let tMid = (yMid - y) / (xMid - x);
segment1 = {
s: y - tMid * x,
t: tMid,
n: 1,
e: 0
};
} else {
segment1 = {
s: s1,
t: t1,
n: j + 1,
e: k1 + (a1 * s1 * s1) + (c1 * t1 * t1) - (2 * d1 * s1) - (2 * h1 * t1) + (2 * b1 * s1 * t1)
};
}

if (error1 + error2 < error1_best + error2_best)
setBest(s1, t1, error1, s2, t2, error2, i, x_prime, x);
}
if (j == sortedSamples.length - 2) {
// Let segment2 be any line through (x[n - 1], y[n - 1]) which meets segment1
// at a point (x’, y’) where x[n - 2] < x' < x[n - 1]. segment2 has no error.
let xMid = (x + sortedSamples[j + 1][complexityIndex]) / 2;
let yMid = s1 + t1 * xMid;
let tMid = (yMid - sortedSamples[j + 1][frameLengthIndex]) / (xMid - sortedSamples[j + 1][complexityIndex]);
segment2 = {
s: y - tMid * x,
t: tMid,
n: 1,
e: 0
};
} else {
segment2 = {
s: s2,
t: t2,
n: sortedSamples.length - (j + 1),
e: k2 + (a2 * s2 * s2) + (c2 * t2 * t2) - (2 * d2 * s2) - (2 * h2 * t2) + (2 * b2 * s2 * t2)
};
}

return {
complexity: x_best,
s1: s1_best,
t1: t1_best,
stdev1: Math.sqrt(error1_best / n1_best),
s2: s2_best,
t2: t2_best,
stdev2: Math.sqrt(error2_best / n2_best),
error: error1_best + error2_best,
n1: n1_best,
n2: n2_best
};
if (this._setOptimal(segment1, segment2, x, sortedSamples[j + 1][complexityIndex], options))
continue

// These values remove the influence of this sample
let G = A + B * x - C * y;
let J = D + E * x - F * y;

let I = c1 - 2 * b1 * x + a1 * xx;
let K = c2 - 2 * b2 * x + a2 * xx;

// Calculate lambda, which divides the weight of this sample between the two lines
let lambda = (G * F + G * K - J * C) / (I * J + G * K);
if (!(lambda > 0 && lambda < 1))
continue;

let lambda1 = 1 - lambda;

segment1 = {
s: (A + lambda * (-h1 * x + d1 * xx + c1 * y - b1 * xy)) / (C - lambda * I),
t: (B + lambda * (h1 - d1 * x - b1 * y + a1 * xy)) / (C - lambda * I),
n: j + 1,
e: (k1 + a1 * s1 * s1 + c1 * t1 * t1 - 2 * d1 * s1 - 2 * h1 * t1 + 2 * b1 * s1 * t1) - lambda * Math.pow(y - (s1 + t1 * x), 2)
};

segment2 = {
s: (D + lambda1 * (-h2 * x + d2 * xx + c2 * y - b2 * xy)) / (F + lambda1 * K),
t: (E + lambda1 * (h2 - d2 * x - b2 * y + a2 * xy)) / (F + lambda1 * K),
n: sortedSamples.length - (j + 1),
e: (k2 + a2 * s2 * s2 + c2 * t2 * t2 - 2 * d2 * s2 - 2 * h2 * t2 + 2 * b2 * s2 * t2) + lambda1 * Math.pow(y - (s2 + t2 * x), 2)
};

this._setOptimal(segment1, segment2, x, sortedSamples[j + 1][complexityIndex], options);
}
}
});

Expand Down

0 comments on commit 40af19e

Please sign in to comment.