Differential privacy for aggregates - add dp_laplace and dp_gaussian built-in functions with budget accounting

benchmark/dp/README.md

@@ -0,0 +1,67 @@
+<!--
+{% comment %}
+Licensed to the Apache Software Foundation (ASF) under one or more
+contributor license agreements.  See the NOTICE file distributed with
+this work for additional information regarding copyright ownership.
+The ASF licenses this file to you under the Apache License, Version 2.0
+(the "License"); you may not use this file except in compliance with
+the License.  You may obtain a copy of the License at
+
+http://www.apache.org/licenses/LICENSE-2.0
+
+Unless required by applicable law or agreed to in writing, software
+distributed under the License is distributed on an "AS IS" BASIS,
+WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+See the License for the specific language governing permissions and
+limitations under the License.
+{% end comment %}
+-->
+
+# DP-FedAvg Benchmark
+
+Benchmarks the `dp_gaussian` built-in by sweeping the privacy budget ε for
+federated logistic regression (DP-FedAvg) on the UCI Adult dataset, and
+plotting the accuracy/privacy trade-off.
+
+## Setup
+
+Run these commands from the repository root. They create a virtual
+environment named `python_venv` and install the Python dependencies listed
+in [scripts/requirements.txt](scripts/requirements.txt):
+
+```bash
+python3 -m venv benchmark/dp/scripts/python_venv
+source benchmark/dp/scripts/python_venv/bin/activate
+pip install -r benchmark/dp/scripts/requirements.txt
+```
+
+Build SystemDS before running the benchmark, if you haven't already:
+
+```bash
+mvn clean package -DskipTests
+```
+
+## Running
+
+With `python_venv` activated, run the benchmark from the repository root:
+
+```bash
+bash benchmark/dp/scripts/run_benchmark.sh
+```
+
+This prepares the dataset, starts the federated workers, sweeps
+epsilon over {0.5, 1, 4, 8} plus a non-private baseline, stops the workers, and
+generates the plots.
+
+## Outputs
+
+- Prepared dataset and per-worker federated shards: [benchmark/dp/data/](data/)
+- Trained models, accuracy logs, and results table: [benchmark/dp/results/](results/)
+- Accuracy vs. epsilon plot: [benchmark/dp/results/accuracy_vs_epsilon.png](results/accuracy_vs_epsilon.png)
+- Utility cost of privacy plot: [benchmark/dp/results/privacy_cost.png](results/privacy_cost.png)
+
+Deactivate the virtual environment when done:
+
+```bash
+deactivate
+```

benchmark/dp/scripts/benchmark_utilities.py

@@ -0,0 +1,30 @@
+#-------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+#-------------------------------------------------------------
+
+
+"""
+Common constants
+"""
+import pathlib
+
+RESULTS_DIR = pathlib.Path("benchmark", "dp", "results")
+DATA_DIR = pathlib.Path("benchmark", "dp", "data")
+

benchmark/dp/scripts/collect_results.py

@@ -0,0 +1,62 @@
+#-------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+#-------------------------------------------------------------
+
+"""
+Parse per-run accuracy files into a single results.csv.
+
+Output columns: label, epsilon, private, accuracy
+"""
+import pathlib, csv
+
+from benchmark_utilities import RESULTS_DIR
+
+rows = []
+
+def parse_acc(path: pathlib.Path) -> float:
+    txt = path.read_text().strip()
+    # SystemDS writes a bare float.
+    return float(txt)
+
+# Non-private baseline.
+baseline_path = RESULTS_DIR / "acc_baseline.txt"
+if baseline_path.exists():
+    rows.append(dict(label="baseline", epsilon="inf",
+                     private=0, accuracy=parse_acc(baseline_path)))
+
+# DP runs.
+for eps in [0.5, 1, 4, 8]:
+    p = RESULTS_DIR / f"acc_eps_{eps}.txt"
+    if p.exists():
+        rows.append(dict(label=f"epsilon={eps}", epsilon=eps,
+                         private=1, accuracy=parse_acc(p)))
+    else:
+        print(f"Warning: {p} not found - skipping")
+
+out = RESULTS_DIR / "results.csv"
+with open(out, "w", newline="") as f:
+    w = csv.DictWriter(f, fieldnames=["label","epsilon","private","accuracy"])
+    w.writeheader()
+    w.writerows(rows)
+
+print(f"Wrote {out}")
+for r in rows:
+    print(f"  {r['label']:12s}  acc={r['accuracy']:.4f}")
+

benchmark/dp/scripts/eval.dml

@@ -0,0 +1,43 @@
+#-------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+#-------------------------------------------------------------
+
+# eval.dml — compute binary classification accuracy on held-out test set.
+# Arguments: data_dir, model_path, out_acc
+data_dir   = $data_dir;
+model_path = $model_path;
+out_acc    = $out_acc;
+
+X_test = read(data_dir + "/X_test.csv",
+              data_type="matrix", value_type="double", format="csv");
+y_test = read(data_dir + "/y_test.csv",
+              data_type="matrix", value_type="double", format="csv");
+w      = read(model_path,
+              data_type="matrix", value_type="double", format="csv");
+
+scores    = X_test %*% w;
+preds     = (scores > 0.0);          # threshold at 0 (log-odds)
+correct   = sum(preds == y_test);
+n_test    = nrow(y_test);
+accuracy  = correct / n_test;
+
+print("Accuracy: " + accuracy);
+write(accuracy, out_acc, format="csv");
+

benchmark/dp/scripts/fedavg_dp.dml

@@ -0,0 +1,151 @@
+#-------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+#-------------------------------------------------------------
+
+# ── fedavg_dp.dml ────────────────────────────────────────────────────────────
+# Arguments (passed via -nvargs):
+#   data_dir   : path to benchmark/data/
+#   epsilon    : DP privacy budget ε  (use 9999 for non-private baseline)
+#   delta      : DP delta             (1e-5 can be used)
+#   clip_norm  : per-example gradient L2-norm clip bound (default 4.0) —
+#                sensitivity = clip_norm / n follows from this, since
+#                clipping each record's gradient contribution to clip_norm
+#                is what makes that bound actually hold.
+#   n_rounds   : number of FedAvg rounds (default 50)
+#   lr         : learning rate          (default 0.1)
+#   w1_rows    : rows in worker 1 shard
+#   w2_rows    : rows in worker 2 shard
+#   w3_rows    : rows in worker 3 shard
+#   w4_rows    : rows in worker 4 shard
+#   n_features : number of features
+#   out        : path to write final weights
+#   private    : 1 = apply DP noise (default), 0 = non-private baseline
+
+data_dir   = $data_dir;
+epsilon    = $epsilon;
+delta      = $delta;
+clip_norm  = ifdef($clip_norm,  4.0);
+n_rounds   = ifdef($n_rounds,   50);
+lr         = ifdef($lr,         0.1);
+private    = ifdef($private,    1);
+out_path   = $out;
+
+# Per-round epsilon: dp_gaussian is called once per round, and
+# DPBudgetAccountant composes the cost of every release against the single
+# total budget set by dp_set_budget(). Spending the full epsilon on every
+# round would exhaust the budget almost immediately, so split it evenly
+# across rounds instead.
+round_epsilon = epsilon / n_rounds;
+
+w1r = $w1_rows;
+w2r = $w2_rows;
+w3r = $w3_rows;
+w4r = $w4_rows;
+d   = $n_features;
+n   = w1r + w2r + w3r + w4r;
+
+# Sensitivity of the released (mean) gradient to a single record changing:
+# with each record's contribution clipped to L2-norm <= clip_norm below,
+# the sum can move by at most clip_norm, so the average moves by clip_norm/n.
+sensitivity = clip_norm / n;
+
+# ── Build federated matrix from 4 local workers ───────────────────────────
+# Row ranges are 0-based [begin, end) for each partition.
+r1s=0;        r1e=w1r;
+r2s=w1r;      r2e=w1r+w2r;
+r3s=w1r+w2r;  r3e=w1r+w2r+w3r;
+r4s=w1r+w2r+w3r; r4e=n;
+
+X = federated(
+    addresses=list(
+        "localhost:8301/" + data_dir + "/worker1/X_train.csv",
+        "localhost:8302/" + data_dir + "/worker2/X_train.csv",
+        "localhost:8303/" + data_dir + "/worker3/X_train.csv",
+        "localhost:8304/" + data_dir + "/worker4/X_train.csv"),
+    ranges=list(
+        list(r1s, 0), list(r1e, d),
+        list(r2s, 0), list(r2e, d),
+        list(r3s, 0), list(r3e, d),
+        list(r4s, 0), list(r4e, d)));
+
+y = federated(
+    addresses=list(
+        "localhost:8301/" + data_dir + "/worker1/y_train.csv",
+        "localhost:8302/" + data_dir + "/worker2/y_train.csv",
+        "localhost:8303/" + data_dir + "/worker3/y_train.csv",
+        "localhost:8304/" + data_dir + "/worker4/y_train.csv"),
+    ranges=list(
+        list(r1s, 0), list(r1e, 1),
+        list(r2s, 0), list(r2e, 1),
+        list(r3s, 0), list(r3e, 1),
+        list(r4s, 0), list(r4e, 1)));
+
+# ── Initialise weights ────────────────────────────────────────────────────
+w = matrix(0, rows=d, cols=1);
+
+if (private == 1) {
+        eps = dp_set_budget($epsilon, $delta);
+}
+
+# ── Training rounds ───────────────────────────────────────────────────────
+for (round in 1:n_rounds) {
+
+    # Forward pass — executes on federated workers.
+    scores = X %*% w;                           # (n × 1), federated
+    probs  = 1.0 / (1.0 + exp(-scores));        # (n × 1), federated
+
+    residuals = probs - y;                       # (n × 1), federated
+
+    # DP noise injection — only when private=1.
+    if (private == 1) {
+        # Per-example L2-norm clipping: each record's raw gradient
+        # contribution is X[i,:] * residual_i, with norm
+        # ||X[i,:]||_2 * |residual_i|. Scaling it down to clip_norm
+        # whenever it exceeds that bound is what makes
+        # sensitivity = clip_norm / n (set above) an actual, provable
+        # bound instead of an arbitrary constant.
+        row_norms         = sqrt(rowSums(X^2));               # (n × 1) ||X[i,:]||_2
+        contrib_norms     = abs(residuals) * row_norms;        # (n × 1) ||X[i,:]*residual_i||_2
+        clip_scale        = clip_norm / pmax(contrib_norms, clip_norm); # (n × 1), in (0,1]
+        clipped_residuals = residuals * clip_scale;            # (n × 1)
+
+        # Gradient aggregation — t(X) %*% clipped_residuals is a (d × 1)
+        # local sum that the coordinator collects in one federated
+        # aggregate instruction.
+        grad = t(X) %*% clipped_residuals / n;                # (d × 1), LOCAL after agg, clipped
+
+        noisy_grad = dp_gaussian(grad,
+                                 "identity",
+                                 sensitivity=sensitivity,
+                                 epsilon=round_epsilon,
+                                 delta=delta);
+        w = w - lr * noisy_grad;
+    } else {
+        # Gradient aggregation — t(X) %*% residuals is an (d × 1) local sum
+        # that the coordinator collects in one federated aggregate instruction.
+        grad = t(X) %*% residuals / n;                        # (d × 1), LOCAL after agg
+        w = w - lr * grad;
+    }
+}
+
+# ── Write model weights ───────────────────────────────────────────────────
+write(w, out_path, format="csv");
+print("FedAvg done. epsilon=" + epsilon + "  rounds=" + n_rounds);
+