Add noise with zCDP in Java

Java:
  * Add noise with zCDP in Java Gaussian Noise

Python dp_accounting:
  * Create DpEvent for the Randomized Response mechanism
  * Parameterize privacy accountant test
  * Add tests for composed NonPrivate event
  * Format RDP accountant files

Change-Id: I45707a23650a9c9be2649c49b035d8de96cfeca2
GitOrigin-RevId: 752cba93d6620ccbbf907fddfcf658ccd816db8b

java/main/com/google/privacy/differentialprivacy/DpPreconditions.java

@@ -38,6 +38,15 @@ static void checkEpsilon(double epsilon) {
         epsilon);
   }
 
+  static void checkRho(double rho) {
+    double rhoLowerBound = 1.0 / (1L << 50);
+    checkArgument(
+        Double.isFinite(rho) && rho >= rhoLowerBound,
+        "rho must be >= %s and < infinity. Provided value: %s",
+        rhoLowerBound,
+        rho);
+  }
+
   static void checkNoiseDelta(Double delta, Noise noise) {
     if (noise.getMechanismType() == MechanismType.LAPLACE
         || noise.getMechanismType() == MechanismType.DISCRETE_LAPLACE) {
@@ -80,6 +89,13 @@ static void checkL1Sensitivity(double l1Sensitivity) {
         l1Sensitivity);
   }
 
+  static void checkL2Sensitivity(double l2Sensitivity) {
+    checkArgument(
+        Double.isFinite(l2Sensitivity) && l2Sensitivity > 0,
+        "l2Sensitivity must be > 0 and finite. Provided value: %s",
+        l2Sensitivity);
+  }
+
   static void checkMaxPartitionsContributed(int maxPartitionsContributed) {
     // maxPartitionsContributed is the user-facing parameter, which is technically the same as
     // L0 sensitivity used by the noise internally.

java/main/com/google/privacy/differentialprivacy/GaussianNoise.java

@@ -99,17 +99,7 @@ public double addNoise(
     checkParameters(l0Sensitivity, lInfSensitivity, epsilon, delta);
 
     double l2Sensitivity = Noise.getL2Sensitivity(l0Sensitivity, lInfSensitivity);
-    double sigma = getSigma(l2Sensitivity, epsilon, delta);
-
-    double granularity = getGranularity(sigma);
-
-    // The square root of n is chosen in a way that places it in the interval between BINOMIAL_BOUND
-    // and BINOMIAL_BOUND / 2. This ensures that the respective binomial distribution consists of
-    // enough Bernoulli samples to closely approximate a Gaussian distribution.
-    double sqrtN = 2.0 * sigma / granularity;
-    long binomialSample = sampleSymmetricBinomial(sqrtN);
-    return SecureNoiseMath.roundToMultipleOfPowerOfTwo(x, granularity)
-        + binomialSample * granularity;
+    return addNoiseDefinedBySigma(x, getSigma(l2Sensitivity, epsilon, delta));
   }
 
   /**
@@ -122,23 +112,56 @@ public long addNoise(
     checkParameters(l0Sensitivity, lInfSensitivity, epsilon, delta);
 
     double l2Sensitivity = Noise.getL2Sensitivity(l0Sensitivity, lInfSensitivity);
-    double sigma = getSigma(l2Sensitivity, epsilon, delta);
+    return addNoiseDefinedBySigma(x, getSigma(l2Sensitivity, epsilon, delta));
+  }
 
-    double granularity = getGranularity(sigma);
+  /**
+   * Adds Gaussian noise to {@code x} such that the output is {@code rho}-zero Concentrated DP
+   * (zCDP) with respect to the specified L_2 sensitivity. For more details on rho-zCDP see
+   * https://eprint.iacr.org/2016/816.pdf
+   */
+  public double addNoiseDefinedByRho(double x, double l2Sensitivity, double rho) {
+    checkParametersForRhozCDP(l2Sensitivity, rho);
+    return addNoiseDefinedBySigma(x, getSigmaForRho(l2Sensitivity, rho));
+  }
+
+  /**
+   * Adds Gaussian noise to {@code x} such that the output is {@code rho}-zCDP with respect to the
+   * specified L_2 sensitivity.
+   */
+  public long addNoiseDefinedByRho(long x, double l2Sensitivity, double rho) {
+    checkParametersForRhozCDP(l2Sensitivity, rho);
+    return addNoiseDefinedBySigma(x, getSigmaForRho(l2Sensitivity, rho));
+  }
+
+  private double addNoiseDefinedBySigma(double x, double noiseSigma) {
+    double granularity = getGranularity(noiseSigma);
 
     // The square root of n is chosen in a way that places it in the interval between BINOMIAL_BOUND
     // and BINOMIAL_BOUND / 2. This ensures that the respective binomial distribution consists of
     // enough Bernoulli samples to closely approximate a Gaussian distribution.
-    double sqrtN = 2.0 * sigma / granularity;
+    double sqrtN = 2.0 * noiseSigma / granularity;
+    long binomialSample = sampleSymmetricBinomial(sqrtN);
+    return SecureNoiseMath.roundToMultipleOfPowerOfTwo(x, granularity)
+        + binomialSample * granularity;
+  }
+
+  private long addNoiseDefinedBySigma(long x, double noiseSigma) {
+    double granularity = getGranularity(noiseSigma);
+
+    // The square root of n is chosen in a way that places it in the interval between BINOMIAL_BOUND
+    // and BINOMIAL_BOUND / 2. This ensures that the respective binomial distribution consists of
+    // enough Bernoulli samples to closely approximate a Gaussian distribution.
+    double sqrtN = 2.0 * noiseSigma / granularity;
     long binomialSample = sampleSymmetricBinomial(sqrtN);
     if (granularity <= 1.0) {
       return x + Math.round(binomialSample * granularity);
-    } else {
-      return SecureNoiseMath.roundToMultiple(x, (long) granularity)
-          + binomialSample * (long) granularity;
     }
+    return SecureNoiseMath.roundToMultiple(x, (long) granularity)
+        + binomialSample * (long) granularity;
   }
 
+
   @Override
   public MechanismType getMechanismType() {
     return MechanismType.GAUSSIAN;
@@ -236,6 +259,11 @@ private void checkParameters(
         Double.isFinite(twoLInf), "2 * lInfSensitivity must be finite but is %s", twoLInf);
   }
 
+  private void checkParametersForRhozCDP(double l2Sensitivity, double rho) {
+    DpPreconditions.checkL2Sensitivity(l2Sensitivity);
+    DpPreconditions.checkRho(rho);
+  }
+
   private void checkConfidenceIntervalParameters(
       int l0Sensitivity, double lInfSensitivity, double epsilon, Double delta, double alpha) {
     DpPreconditions.checkAlpha(alpha);
@@ -283,6 +311,15 @@ public static double getSigma(double l2Sensitivity, double epsilon, double delta
     return upperBound;
   }
 
+  /*
+   * Computes sigma of Gaissian noise to satisify rho Zero Concentrated DP (rho-zCDP).
+   * For more details on rho-zCDP see https://eprint.iacr.org/2016/816.pdf
+   */
+  public static double getSigmaForRho(double l2Sensitivity, double rho) {
+    // From https://eprint.iacr.org/2016/816.pdf Propositon 6.
+    return l2Sensitivity / Math.sqrt(2 * rho);
+  }
+
   /**
    * Returns the smallest delta such that the Gaussian mechanism with standard deviation {@code
    * sigma} obtains {@code (epsilon, delta)}-differential privacy with respect to the provided L_2

java/tests/com/google/privacy/differentialprivacy/GaussianNoiseTest.java

@@ -266,7 +266,7 @@ public void addNoise_deltaNull_throwsException() {
                 DEFAULT_L_0_SENSITIVITY,
                 DEFAULT_L_INF_SENSITIVITY,
                 DEFAULT_EPSILON,
-                 /* delta= */ 0.0));
+                /* delta= */ 0.0));
   }
 
   @Test
@@ -348,6 +348,26 @@ public void addNoise_lInfSensitivityNan_throwsException() {
     assertThrows(IllegalArgumentException.class, () -> NOISE.addNoise(0, 1, NaN, 1, DEFAULT_DELTA));
   }
 
+  @Test
+  public void addNoiseDefinedByRho_Nan_rho_throwsException() {
+    assertThrows(IllegalArgumentException.class, () -> NOISE.addNoiseDefinedByRho(0, 1, NaN));
+  }
+
+  @Test
+  public void addNoiseDefinedByRho_zero_rho_throwsException() {
+    assertThrows(IllegalArgumentException.class, () -> NOISE.addNoiseDefinedByRho(0, 1, 0));
+  }
+
+  @Test
+  public void addNoiseDefinedByRho_negativeSensitivity_throwsException() {
+    assertThrows(IllegalArgumentException.class, () -> NOISE.addNoiseDefinedByRho(0, -1, 1));
+  }
+
+  @Test
+  public void addNoiseDefinedByRho_zeroSensitivity_throwsException() {
+    assertThrows(IllegalArgumentException.class, () -> NOISE.addNoiseDefinedByRho(0, 0, 1));
+  }
+
   @Test
   public void addNoise_returnsMultipleOfGranularity() {
     SecureRandom random = new SecureRandom();
@@ -395,8 +415,8 @@ public void addNoise_returnsMultipleOfGranularity() {
   public void addNoise_integralX_returnsMultipleOfGranularity() {
     SecureRandom random = new SecureRandom();
     for (int i = 0; i < NUM_SAMPLES; i++) {
-      // The rounding pricess should be independent of the value of x. Set x to a value between
-      // -1*10^6 and 10^6 at random should covere a broad range of congruence classes.
+      // The rounding process should be independent of the value of x. Set x to a value between
+      // -1*10^6 and 10^6 at random should cover a broad range of congruence classes.
       long x = (long) random.nextInt(2000000) - 1000000;
 
       // The following choice of epsilon, delta, l0 sensitivity and linf sensitivity should result
@@ -450,4 +470,39 @@ public void sampleSymmetricBinomial_hasAccurateStatisticalProperties() {
   public void getMechanismType_returnsGaussian() {
     assertThat(NOISE.getMechanismType()).isEqualTo(GAUSSIAN);
   }
+
+  @Test
+  public void getSigmaForRho_returnsCorrectly() {
+    // sigma = l2Sensitivity / sqrt(2*rho)
+    assertThat(GaussianNoise.getSigmaForRho(/* l2Sensitivity= */ 1.0, /* rho= */ 1))
+        .isWithin(1e-12)
+        .of(0.7071067811865475);
+    assertThat(GaussianNoise.getSigmaForRho(/* l2Sensitivity= */ 3.0, /* rho= */ 1))
+        .isWithin(1e-12)
+        .of(2.1213203435596424);
+    assertThat(GaussianNoise.getSigmaForRho(/* l2Sensitivity= */ 1.0, /* rho= */ 2))
+        .isWithin(1e-12)
+        .of(0.5);
+    assertThat(GaussianNoise.getSigmaForRho(/* l2Sensitivity= */ 10.0, /* rho= */ 8))
+        .isWithin(1e-12)
+        .of(2.5);
+  }
+
+  @Test
+  public void addNoiseDefinedByRho_hasAccurateStatisticalProperties() {
+    ImmutableList.Builder<Double> samples = ImmutableList.builder();
+    for (int i = 0; i < NUM_SAMPLES; i++) {
+      samples.add(NOISE.addNoiseDefinedByRho(/* x */ 1.0, /* l2Sensitivity */ 10.0, /* rho */ 2));
+    }
+    Stats stats = Stats.of(samples.build());
+
+    double variance = 25; // std_dev = 10/sqrt(2*2) = 5
+    // The tolerance is chosen according to the 99.9995% quantile of the anticipated distributions
+    // of the sample mean and variance. Thus, the test falsely rejects with a probability of 10^-5.
+    final double stdNormQuantile = 4.41717; // 99.9995% quantile of the standard normal distribution
+    double sampleMeanTolerance = stdNormQuantile * Math.sqrt(variance / NUM_SAMPLES);
+    double sampleVarianceTolerance = stdNormQuantile * variance * Math.sqrt(2.0 / NUM_SAMPLES);
+    assertThat(stats.mean()).isWithin(sampleMeanTolerance).of(1.0);
+    assertThat(stats.populationVariance()).isWithin(sampleVarianceTolerance).of(variance);
+  }
 }

python/dp_accounting/dp_accounting/__init__.py

@@ -27,6 +27,7 @@
 from dp_accounting.dp_event import NonPrivateDpEvent
 from dp_accounting.dp_event import NoOpDpEvent
 from dp_accounting.dp_event import PoissonSampledDpEvent
+from dp_accounting.dp_event import RandomizedResponseDpEvent
 from dp_accounting.dp_event import SampledWithoutReplacementDpEvent
 from dp_accounting.dp_event import SampledWithReplacementDpEvent
 from dp_accounting.dp_event import SelfComposedDpEvent

python/dp_accounting/dp_accounting/dp_event.py

@@ -242,6 +242,24 @@ class UnsupportedDpEvent(DpEvent):
   """
 
 
+@attr.s(frozen=True, slots=True, auto_attribs=True)
+class RandomizedResponseDpEvent(DpEvent):
+  """Represents application of the randomized response mechanism.
+
+  The Randomized Response over k buckets with noise parameter p takes in an
+  input which is one of the k buckets. With probability 1 - p, it simply
+  outputs the input bucket. Otherwise, with probability p, it outputs a bucket
+  drawn uniformly at random from the k buckets.
+
+  The noise parameter p can be any value in [0, 1], with p=0 corresponding to
+  the case where the mechanism always outputs the input bucket, and p=1
+  corresponding to the case where the mechanism outputs a bucket drawn
+  uniformly at random from the k buckets regardless of the input bucket.
+  """
+  noise_parameter: float
+  num_buckets: int
+
+
 @attr.s(frozen=True, slots=True, auto_attribs=True)
 class GaussianDpEvent(DpEvent):
   """Represents an application of the Gaussian mechanism.

python/dp_accounting/dp_accounting/privacy_accountant_test.py

@@ -21,13 +21,18 @@
 from typing import Collection
 
 from absl.testing import absltest
+from absl.testing import parameterized
 
 from dp_accounting import dp_event
 from dp_accounting import privacy_accountant
 
 
+class UnknownDpEvent(dp_event.DpEvent):
+  pass
+
+
 @absltest.skipThisClass('only intended to be run by subclasses')
-class PrivacyAccountantTest(absltest.TestCase):
+class PrivacyAccountantTest(parameterized.TestCase):
 
   def _make_test_accountants(
       self) -> Collection[privacy_accountant.PrivacyAccountant]:
@@ -43,29 +48,27 @@ def _make_test_accountants(
   def test_make_test_accountants(self):
     self.assertNotEmpty(self._make_test_accountants())
 
-  def test_unsupported(self):
-
-    class UnknownDpEvent(dp_event.DpEvent):
-      pass
-
+  @parameterized.product(
+      unsupported_event=(dp_event.UnsupportedDpEvent(), UnknownDpEvent()),
+      nest_fn=(
+          lambda event_: event_,
+          lambda event_: dp_event.ComposedDpEvent([event_]),
+          lambda event_: dp_event.SelfComposedDpEvent(event_, 10),
+      ),
+  )
+  def test_unsupported(self, unsupported_event, nest_fn):
+    event_ = nest_fn(unsupported_event)
     for accountant in self._make_test_accountants():
-      for unsupported_event in [
-          dp_event.UnsupportedDpEvent(),
-          UnknownDpEvent()
-      ]:
-        for nested_unsupported_event in [
-            unsupported_event,
-            dp_event.SelfComposedDpEvent(unsupported_event, 10),
-            dp_event.ComposedDpEvent([unsupported_event])
-        ]:
-          composition_error = accountant._maybe_compose(
-              nested_unsupported_event, count=1, do_compose=False)
-          self.assertIsNotNone(composition_error)
-          self.assertEqual(composition_error.invalid_event, unsupported_event)
-          self.assertFalse(accountant.supports(nested_unsupported_event))
-          with self.assertRaisesRegex(privacy_accountant.UnsupportedEventError,
-                                      'caused by subevent'):
-            accountant.compose(nested_unsupported_event)
+      composition_error = accountant._maybe_compose(
+          event_, count=1, do_compose=False
+      )
+      self.assertIsNotNone(composition_error)
+      self.assertEqual(composition_error.invalid_event, unsupported_event)
+      self.assertFalse(accountant.supports(event_))
+      with self.assertRaisesRegex(
+          privacy_accountant.UnsupportedEventError, 'caused by subevent'
+      ):
+        accountant.compose(event_)
 
   def test_no_events(self):
     for accountant in self._make_test_accountants():
@@ -96,11 +99,18 @@ def test_no_op(self):
         # Implementing `get_delta` is optional.
         pass
 
-  def test_non_private(self):
+  @parameterized.parameters(
+      dp_event.NonPrivateDpEvent(),
+      dp_event.ComposedDpEvent([dp_event.NonPrivateDpEvent()]),
+      dp_event.ComposedDpEvent(
+          [dp_event.NoOpDpEvent(), dp_event.NonPrivateDpEvent()]
+      ),
+      dp_event.SelfComposedDpEvent(dp_event.NonPrivateDpEvent(), 10),
+  )
+  def test_non_private(self, non_private_event):
     for accountant in self._make_test_accountants():
-      event = dp_event.NonPrivateDpEvent()
-      self.assertTrue(accountant.supports(event))
-      accountant.compose(event)
+      self.assertTrue(accountant.supports(non_private_event))
+      accountant.compose(non_private_event)
       self.assertEqual(accountant.get_epsilon(0.99), float('inf'))
       self.assertEqual(accountant.get_epsilon(0), float('inf'))
       self.assertEqual(accountant.get_epsilon(1), float('inf'))