Use bit search for max epsilon (#1429)

ts/src/flexible-event/privacy.test.ts

@@ -9,6 +9,7 @@ import {
   binaryEntropy,
   flipProbabilityDp,
   maxInformationGain,
+  epsilonToBoundInfoGainAndDp,
 } from './privacy'
 
 const flipProbabilityTests = [
@@ -101,6 +102,94 @@ void test('maxInformationGain', async (t) => {
   )
 })
 
+void test('epsilonToBoundInfoGainAndDp', async (t) => {
+  const infoGainUppers = [11.5, 6.5]
+  const epsilonUpper = 14
+  const numStatesRange = 100000
+  const numTests = 500
+
+  await Promise.all(
+    Array(numTests)
+      .fill(0)
+      .map((_, i) =>
+        t.test(`${i}`, () => {
+          const numStates = Math.ceil(Math.random() * numStatesRange)
+          const infoGainUpper = infoGainUppers[Math.round(Math.random())]!
+
+          const epsilon = epsilonToBoundInfoGainAndDp(
+            numStates,
+            infoGainUpper,
+            epsilonUpper
+          )
+
+          assert(maxInformationGain(numStates, epsilon) <= infoGainUpper)
+
+          if (epsilon < epsilonUpper) {
+            assert(
+              maxInformationGain(numStates, epsilon + 1e-15) > infoGainUpper
+            )
+          }
+        })
+      )
+  )
+})
+
+const epsilonSearchTests = [
+  {
+    numStates: 5545,
+    infoGainUpper: 6.5,
+    epsilonUpper: 14,
+    expected: 9.028709123768687,
+  },
+  {
+    numStates: 2106,
+    infoGainUpper: 6.5,
+    epsilonUpper: 14,
+    expected: 8.366900276574821,
+  },
+  {
+    numStates: 16036,
+    infoGainUpper: 6.5,
+    epsilonUpper: 14,
+    expected: 9.829279343808693,
+  },
+  {
+    numStates: 84121,
+    infoGainUpper: 11.5,
+    epsilonUpper: 14,
+    expected: 12.45087042924698,
+  },
+  {
+    numStates: 24895,
+    infoGainUpper: 11.5,
+    epsilonUpper: 14,
+    expected: 11.723490703852157,
+  },
+  {
+    numStates: 3648,
+    infoGainUpper: 11.5,
+    epsilonUpper: 14,
+    expected: 12.233993328032184,
+  },
+]
+
+void test('epsilonSearch', async (t) => {
+  await Promise.all(
+    epsilonSearchTests.map((tc) =>
+      t.test(`${tc.numStates}, ${tc.infoGainUpper}, ${tc.epsilonUpper}`, () => {
+        assert.deepStrictEqual(
+          tc.expected,
+          epsilonToBoundInfoGainAndDp(
+            tc.numStates,
+            tc.infoGainUpper,
+            tc.epsilonUpper
+          )
+        )
+      })
+    )
+  )
+})
+
 const binaryEntropyTests = [
   { x: 0, expected: 0 },
   { x: 0.5, expected: 1 },

ts/src/flexible-event/privacy.ts

@@ -184,33 +184,49 @@ export function maxInformationGain(numStates: number, epsilon: number): number {
 
 // Returns the effective epsilon needed to satisfy an information gain bound
 // given a number of output states in the q-ary symmetric channel.
-function epsilonToBoundInfoGainAndDp(
+//
+// The exponent section of the double is used as a power with base 2, which is
+// multiplied by the significand, which is at least 1 and less than 2. In our
+// case, the sign bit remains unset since we use a positive epsilon. The double
+// is 2^exponent * significand. Since the significand is at least 1, changing it
+// can only lower 2^exponent to at least its own value. And since the
+// significand is less than 2, changing it can only raise 2^exponent to a value
+// less than 2^(exponent + 1). We can therefore first find the highest
+// 2^exponent less than or equal to max-settable-event-level-epsilon that
+// produces information gain within limit, and then search for a significand
+// that raises the overall value as much as possible.
+//
+// In an additive binary representation, the value represented by one bit is
+// higher than all the values added by lower bits together. This inequality
+// holds when multiplying by the constant, 2^exponent. This means that if we
+// search the significand from high bits to low, each choice to set a bit either
+// is already too high, or provides the opportunity to get closer to a higher
+// target using some combination of the lower bits.
+export function epsilonToBoundInfoGainAndDp(
   numStates: number,
   infoGainUpperBound: number,
-  epsilonUpperBound: number,
-  tolerance: number = 0.00001
+  epsilonUpperBound: number
 ): number {
-  // Just perform a simple binary search over values of epsilon.
-  let epsLow = 0
-  let epsHigh = epsilonUpperBound
+  const buffer = new ArrayBuffer(8)
+  const dataView = new DataView(buffer)
 
-  for (;;) {
-    const epsilon = (epsHigh + epsLow) / 2
-    const infoGain = maxInformationGain(numStates, epsilon)
+  for (let bit = 1n << 62n; bit > 0n; bit >>= 1n) {
+    dataView.setBigUint64(0, dataView.getBigUint64(0) | bit)
 
-    if (infoGain > infoGainUpperBound) {
-      epsHigh = epsilon
+    const epsilon = dataView.getFloat64(0)
+    if (epsilon > epsilonUpperBound) {
+      dataView.setBigUint64(0, dataView.getBigUint64(0) & ~bit)
       continue
     }
 
-    // Allow slack by returning something slightly non-optimal (governed by the tolerance)
-    // that still meets the privacy bar. If epsHigh == epsLow we're now governed by the epsilon
-    // bound and can return.
-    if (infoGain < infoGainUpperBound - tolerance && epsHigh !== epsLow) {
-      epsLow = epsilon
-      continue
-    }
+    const infoGain = maxInformationGain(numStates, epsilon)
 
-    return epsilon
+    if (infoGain > infoGainUpperBound) {
+      dataView.setBigUint64(0, dataView.getBigUint64(0) & ~bit)
+    } else if (epsilon === epsilonUpperBound) {
+      return epsilon
+    }
   }
+
+  return dataView.getFloat64(0)
 }