Miscellaneous simulator improvements

impl/README.md

@@ -3,3 +3,6 @@ Usage:
 ```sh
 npm install && npm run pack && npm run serve-local
 ```
+
+A live version of the simulator can be found at
+https://w3c.github.io/ppa/simulator.html.

impl/package.json

@@ -17,6 +17,7 @@
     "html-webpack-plugin": "^5.6.3",
     "http-server": "^14.1.1",
     "prettier": "^3.6.2",
+    "simple-statistics": "^7.8.8",
     "ts-loader": "^9.5.2",
     "typescript": "^5.9.2",
     "typescript-eslint": "^8.39.0",

impl/src/backend.test.ts

@@ -0,0 +1,159 @@
+import { fairlyAllocateCredit } from "./backend";
+
+import { strict as assert } from "assert";
+import test from "node:test";
+import { inverseErrorFunction } from "simple-statistics";
+
+interface FairlyAllocateCreditTestCase {
+  name: string;
+  credit: number[];
+  value: number;
+  needsRand?: boolean;
+}
+
+function noRand(): number {
+  throw new Error("no rand expected");
+}
+
+type Interval = [min: number, max: number];
+
+// https://en.wikipedia.org/wiki/Probit
+function normalPpf(q: number, stdDev: number): number {
+  return stdDev * Math.sqrt(2) * inverseErrorFunction(2 * q - 1);
+}
+
+const minNForIntervalApprox = 1000;
+
+function getIntervalApprox(n: number, p: number, alpha: number): Interval {
+  if (n < minNForIntervalApprox) {
+    throw new RangeError(`n must be >= ${minNForIntervalApprox}`);
+  }
+
+  // Approximates a binomial distribution with a normal distribution which is a bit
+  // simpler as it is symmetric.
+  const mean = n * p;
+  const variance = mean * (1 - p);
+  const diff = normalPpf(1 - alpha / 2, Math.sqrt(variance));
+  return [mean - diff, mean + diff];
+}
+
+function getAllIntervals(
+  n: number,
+  creditFractions: readonly number[],
+  alphaTotal: number,
+): Interval[] {
+  // We are testing one hypothesis per dimension, so divide `alphaTotal` by
+  // the number of dimensions: https://en.wikipedia.org/wiki/Bonferroni_correction
+  const alpha = alphaTotal / creditFractions.length;
+  return creditFractions.map((cf) => getIntervalApprox(n, cf, alpha));
+}
+
+function runFairlyAllocateCreditTest(
+  tc: Readonly<FairlyAllocateCreditTestCase>,
+): void {
+  // TODO: replace with precise sum
+  const sumCredit = tc.credit.reduce((a, b) => a + b, 0);
+  const normalizedFloatCredit = tc.credit.map((item) => item / sumCredit);
+
+  const [rand, k] = tc.needsRand ? [Math.random, 1000] : [noRand, 1];
+
+  const totals = new Array<number>(tc.credit.length).fill(0);
+
+  for (let n = 0; n < k; ++n) {
+    const actualCredit = fairlyAllocateCredit(tc.credit, tc.value, rand);
+
+    assert.equal(actualCredit.length, tc.credit.length);
+
+    for (const [j, actual] of actualCredit.entries()) {
+      assert.ok(Number.isInteger(actual));
+
+      const normalized = normalizedFloatCredit[j]! * tc.value;
+      const diff = Math.abs(actual - normalized);
+      assert.ok(
+        diff < 1,
+        `credit error >= 1: actual=${actual}, normalized=${normalized}`,
+      );
+
+      totals[j]! += actual / tc.value;
+    }
+
+    assert.equal(
+      // TODO: replace with precise sum
+      actualCredit.reduce((a, b) => a + b, 0),
+      tc.value,
+      `actual credit does not sum to value: ${actualCredit.join(", ")}`,
+    );
+  }
+
+  const alpha = 0.00001; // Probability of test failing at random.
+
+  const intervals: Interval[] =
+    k > 1
+      ? getAllIntervals(
+          k,
+          normalizedFloatCredit.map((c) => c - Math.floor(c)),
+          alpha,
+        )
+      : normalizedFloatCredit.map((c) => [c, c]);
+
+  for (const [j, total] of totals.entries()) {
+    const [min, max] = intervals[j]!;
+    assert.ok(
+      total >= min && total <= max,
+      `total for credit[${j}] ${total} not in ${1 - alpha} confidence interval [${min}, ${max}]`,
+    );
+  }
+}
+
+const testCases: FairlyAllocateCreditTestCase[] = [
+  {
+    name: "credit-equal-to-value",
+    credit: [1],
+    value: 1,
+    needsRand: false,
+  },
+  {
+    name: "credit-less-than-value",
+    credit: [2],
+    value: 3,
+    needsRand: false,
+  },
+  {
+    name: "credit-less-than-1",
+    credit: [0.25],
+    value: 4,
+    needsRand: false,
+  },
+  {
+    name: "2-credit-divides-value-evenly",
+    credit: [3, 1],
+    value: 8,
+    needsRand: false,
+  },
+  {
+    name: "3-credit-divides-value-evenly",
+    credit: [2, 1, 1],
+    value: 8,
+    needsRand: false,
+  },
+  {
+    name: "2-credit-divides-value-unevenly",
+    credit: [1, 1],
+    value: 5,
+    needsRand: true,
+  },
+  {
+    name: "3-credit-divides-value-unevenly",
+    credit: [2, 1, 1],
+    value: 5,
+    needsRand: true,
+  },
+];
+
+void test("fairlyAllocateCredit", async (t) => {
+  await Promise.all(
+    testCases.map((tc) =>
+      t.test(tc.name, () => runFairlyAllocateCreditTest(tc)),
+    ),
+  );
+});

impl/src/backend.ts

@@ -50,7 +50,7 @@ interface ValidatedConversionOptions {
 }
 
 interface ValidatedLogicOptions {
-  credit: number[];
+  credit: readonly number[];
 }
 
 export function days(days: number): Temporal.Duration {
@@ -557,7 +557,7 @@ export class Backend {
     matchedImpressions: Set<Impression>,
     histogramSize: number,
     value: number,
-    credit: number[],
+    credit: readonly number[],
   ): number[] {
     if (matchedImpressions.size === 0) {
       throw new DOMException(
@@ -582,12 +582,13 @@ export class Backend {
 
     const lastNImpressions = sortedImpressions.slice(0, N);
 
-    const normalizedCredit = this.#fairlyAllocateCredit(credit, value);
+    const normalizedCredit = fairlyAllocateCredit(credit, value, () =>
+      this.#delegate.random(),
+    );
 
     const histogram = allZeroHistogram(histogramSize);
 
-    for (let i = 0; i < lastNImpressions.length; ++i) {
-      const impression = lastNImpressions[i]!;
+    for (const [i, impression] of lastNImpressions.entries()) {
       const value = normalizedCredit[i];
       const index = impression.histogramIndex;
       if (index < histogram.length) {
@@ -597,7 +598,7 @@ export class Backend {
     return histogram;
   }
 
-  #encryptReport(report: number[]): Uint8Array {
+  #encryptReport(report: readonly number[]): Uint8Array {
     void report;
     return new Uint8Array(0); // TODO
   }
@@ -606,10 +607,7 @@ export class Backend {
     const period = this.#delegate.privacyBudgetEpoch.total("seconds");
     let start = this.#epochStartStore.get(site);
     if (start === undefined) {
-      const p = this.#delegate.random();
-      if (!(p >= 0 && p < 1)) {
-        throw new RangeError("random must be in the range [0, 1)");
-      }
+      const p = checkRandom(this.#delegate.random());
       const dur = Temporal.Duration.from({
         seconds: p * period,
       });
@@ -667,42 +665,54 @@ export class Backend {
       );
     });
   }
+}
 
-  #fairlyAllocateCredit(credit: number[], value: number): number[] {
-    const sumCredit = credit.reduce((a, b) => a + b, 0);
+function checkRandom(p: number): number {
+  if (!(p >= 0 && p < 1)) {
+    throw new RangeError("random must be in the range [0, 1)");
+  }
+  return p;
+}
 
-    const roundedCredit = credit.map((item) => (value * item) / sumCredit);
+export function fairlyAllocateCredit(
+  credit: readonly number[],
+  value: number,
+  rand: () => number,
+): number[] {
+  // TODO: replace with precise sum
+  const sumCredit = credit.reduce((a, b) => a + b, 0);
 
-    let idx1 = 0;
+  const roundedCredit = credit.map((item) => (value * item) / sumCredit);
 
-    for (let n = 1; n < roundedCredit.length; ++n) {
-      let idx2 = n;
+  let idx1 = 0;
 
-      const frac1 = roundedCredit[idx1]! - Math.floor(roundedCredit[idx1]!);
-      const frac2 = roundedCredit[idx2]! - Math.floor(roundedCredit[idx2]!);
-      if (frac1 === 0 && frac2 === 0) {
-        continue;
-      }
+  for (let n = 1; n < roundedCredit.length; ++n) {
+    let idx2 = n;
 
-      const [incr1, incr2] =
-        frac1 + frac2 > 1 ? [1 - frac1, 1 - frac2] : [-frac1, -frac2];
+    const frac1 = roundedCredit[idx1]! - Math.floor(roundedCredit[idx1]!);
+    const frac2 = roundedCredit[idx2]! - Math.floor(roundedCredit[idx2]!);
+    if (frac1 === 0 && frac2 === 0) {
+      continue;
+    }
 
-      const p1 = incr2 / (incr1 + incr2);
+    const [incr1, incr2] =
+      frac1 + frac2 > 1 ? [1 - frac1, 1 - frac2] : [-frac1, -frac2];
 
-      const r = this.#delegate.random();
+    const p1 = incr2 / (incr1 + incr2);
 
-      let incr;
-      if (r < p1) {
-        incr = incr1;
-        [idx1, idx2] = [idx2, idx1];
-      } else {
-        incr = incr2;
-      }
+    const r = checkRandom(rand());
 
-      roundedCredit[idx2]! += incr;
-      roundedCredit[idx1]! -= incr;
+    let incr;
+    if (r < p1) {
+      incr = incr1;
+      [idx1, idx2] = [idx2, idx1];
+    } else {
+      incr = incr2;
     }
 
-    return roundedCredit.map((item) => Math.round(item));
+    roundedCredit[idx2]! += incr;
+    roundedCredit[idx1]! -= incr;
   }
+
+  return roundedCredit.map((item) => Math.round(item));
 }