ABSTRACT Matrix mechanisms are often used to provide unbiased differentially private query answers when publishing statistics or creating synthetic data. Recent work has developed matrix mechanisms, such as ResidualPlanner and Weighted Fourier Factorizations, that scale to high dimensional datasets while providing optimality guarantees for workloads such as marginals and circular product queries. They operate by adding noise to a linearly independent set of queries that can compactly represent the desired workloads. In this paper, we present QuerySmasher, an alternative scalable approach based on a divide-and-conquer strategy. Given a workload that can be answered from various data marginals, QuerySmasher splits each query into subqueries and re-assembles the pieces into mutually orthogonal subworkloads. These subworkloads represent small, low-dimensional problems that can be independently and optimally answered by existing low-dimensional matrix mechanisms. QuerySmasher then stitches these solutions together to answer queries in the original workload. We show that QuerySmasher subsumes prior work, like ResidualPlanner (RP), ResidualPlanner+ (RP+), and Weighted Fourier Factorizations (WFF). We prove that it can dominate those approaches, under sum squared error, for all workloads. We also experimentally demonstrate the scalability and accuracy of QuerySmasher.

INTRODUCTION

Confidential datasets about individuals and businesses provide a wealth of information that is useful for scientific, economic, and policy purposes. Clearly, they cannot be released without modification and so a variety of perturbed summary statistics are released instead. The fundamental theorem of database reconstruction [5] demonstrates limits on how many approximate query answers can be released before privacy breaks down under any reasonable notion, and differential privacy [4, 6, 8, 23] has become a de facto standard for navigating the privacy/utility tradeoff. Although the datasets can be high-dimensional, the queries of interest [27–29] often involve a small set of variables at a time. Hence, the queries of interest can be answered from a collection of data marginals. For example, age range queries for each race can be answered by a 2-way marginal on age and race, 2-d range queries

on age and income for each gender can be answered by a 3-way marginal on age, gender, and income, etc. These are examples of linear queries, and they can be answered using a class of differentially private mechanisms known as matrix mechanisms [9, 15, 17, 20, 32, 33, 35, 36]. Given a workload of linear queries, these methods create a set of strategy queries that can be answered with low noise. The original workload queries can then be answered by linear postprocessing of the strategy query answers. This is a powerful framework that has been used in the 2020 Census disclosure avoidance system [1] and is also a key element of more complex algorithms [10, 19, 22, 24, 38, 39]. One of the key challenges of the matrix mechanism line of work is scalability to high-dimensional datasets. If a dataset has attributes 𝐴1, , . . . , 𝐴𝑘 with domains sizes 𝑑 1, . . . , 𝑑 𝑘 , then optimizing a matrix mechanism to minimize the sum squared error of an arbitrary query workload can as a semi-definite programming be formulated 2 Î𝑘 problem with 𝑂 variables [35]. Solving such an op𝑖=1 𝑑𝑖 timization problem is only feasible for low-dimensional datasets that have a small number of attributes. However, most applications require answers to structured query workloads rather than arbitrary workloads. Specifically, each linear query accesses a small number of attributes (different queries can access different attributes) and so the entire workload can be computed using a collection of lower dimensional marginals [27–29]. McKenna et al. [20] first showed how to scale matrix mechanism construction by exploiting Î such 𝑘 structure. However, their algorithm still required 𝑂 𝑖=1 𝑑𝑖 space, which becomes a bottleneck for scalability. A recent line of work called ResidualPlanner (RP) [33], ResidualPlanner+ (RP+) [32], and Weighted Fourier Factorizations (WFF) [15] overcame this bottleneck by proposing matrix mechanisms whose time and space complexity scales with the output size (i.e., size of the marginals needed to answer the workload queries) rather than the exponentially sized data domain. This is enabled by special types of linearly independent queries, such as the residual basis [33] or multidimensional Fourier basis [15], that can compactly represent any collection of marginals and therefore they can compactly represent any linear workload that is computable from collections of

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

marginals (such as range queries, prefix-sum queries, etc.). Surprisingly, a proper choice of basis also guarantees optimality (among matrix mechanisms) for expressive classes of queries like marginals [33] or, more generally, circular product queries [15]. In this paper, we present QuerySmasher, a novel approach to matrix mechanism scalability that is provably more accurate than this prior work [15, 32, 33] under the weighted sum squared error metric. It breaks apart each query into subqueries and groups those pieces into subworkloads. These subworkloads are mutually orthogonal, so that matrix mechanisms for each subworkload can be constructed independently. The subworkloads are also low dimensional, so that optimal low-dimensional techniques [35] can be used to generate the matrix mechanisms for each subworkload. Afterwards, these matrix mechanisms are combined together to create a mechanism for the original workload. This divide-andconquer approach ensures scalability. QuerySmasher can also be viewed from an alternate lens of automatically searching for a good workload-dependent compact linear basis that can represent marginals, instead of always using the residual basis [33] or Fourier basis [15]. This ensures high accuracy. To summarize, our contributions are the following: • We propose a new matrix mechanism, called QuerySmasher, for computing privacy-preserving answers to any workload of linear queries that are computable from sets of marginals. QuerySmasher uses a novel query decomposition approach to achieve both scalability and high accuracy. • We mathematically prove that QuerySmasher can dominate the scalable compact-linear-basis approaches RP [33], RP+ [32], and WFF [15], in terms of weighted sum squared error, for every workload. Specifically, if the optimal low-dimensional matrix mechanism construction [37] for each subworkload is replaced by approximate methods, the result becomes mathematically equivalent to RP, RP+, and WFF. • We empirically evaluate QuerySmasher against competing scalable matrix mechanisms [15, 20, 32, 33] to show that its theoretical accuracy and scalability advantages translate into practice. This paper is organized as follows. In Section 2, we introduce our notation. In Section 3 we explain background on differential privacy and matrix mechanisms. In Section 4, we discuss related work. We present an intuitive overview of QuerySmasher in Section 5 and the technical details in Section 6. In Section 7, we prove that QuerySmasher dominates RP [33], RP+ [32], and WFF [15], in terms of sum squared error, for every workload. We present experiments in Section 8 and conclusions in Section 9. Due to space constraints, proofs can be found in the appendix.

NOTATION

We use the following linear algebra notation. Let ⊗ be the Kronecker product; e.g., if the matrix V =

𝑣0,0 𝑣0,1 𝑣0,2 𝑣1,0 𝑣1,1 𝑣1,2 𝑣2,0 𝑣2,1 𝑣2,2

and W is another matrix 𝑣 W 𝑣 W 𝑣 W 0,0 0,1 0,2 then V⊗W can be written as the block matrix 𝑣1,0 W 𝑣1,1 W 𝑣1,2 W . 𝑣2,0 W 𝑣2,1 W 𝑣2,2 W

Let I𝑑 be the 𝑑 × 𝑑 identity matrix. Let 1𝑑 be the d-dimensional column vector whose entries are all 1. Data. Let 𝐴1, . . . , 𝐴𝑘 be attributes having discrete domains. For 𝑖 = 1, . . . , 𝑘, let 𝑑 𝑖 denote the domain size of attribute 𝐴𝑖 and let its

domain elements be 𝑎 0(𝑖 ) , 𝑎 1(𝑖 ) , . . . , 𝑎𝑑(𝑖 )−1 . A value 𝑎 (𝑖𝑗 ) for attribute 𝑖 𝐴𝑖 can also be represented in one-hot notation as a 𝑑 𝑖 -dimensional column vector e (𝑖𝑗 ) , where the 𝑗 th component is 1 and all others are 0. We index vectors using square brackets (e.g., e (𝑖𝑗 ) [ 𝑗] = 1). A record 𝑟 on attributes 𝐴1, . . . , 𝐴𝑘 is a tuple that specifies a value for each attribute. It can also be represented as a vector of diÎ mension 𝑘𝑖=1 𝑑 𝑖 by taking the Kronecker product ⊗ of the one-hot encodings of its attribute values. For example, if the domain of 𝐴1 is {low, medium, high}, the domain of 𝐴2 is {yes, no}, and the domain of 𝐴3 is {0, 1, 2, 3, 4}, then the record 𝑟 = (medium, yes, 3) can be represented as e1(1) ⊗ e0(2) ⊗ e3(3) = [0, 1, 0] ⊤ ⊗[1, 0] ⊤ ⊗[0, 0, 0, 1, 0] ⊤ . A dataset D = {𝑟 1, . . . , 𝑟 𝑛 } is a set of records, where each record 𝑟 𝑖 corresponds to one individual. In addition to this set view, a dataset D can also be viewed as a tensor/multidimensional array T, in which the entry T[ 𝑗1, 𝑗2, . . . , 𝑗𝑘 ] is the number of records for (2) (𝑘 ) which 𝐴1 = 𝑎 (1) 𝑗1 , 𝐴2 = 𝑎 𝑗2 , . . . , 𝐴𝑘 = 𝑎 𝑗𝑘 . A dataset D can also be represented as a column vector x® which is the sum of the vector representations of the records in the dataset. There is a simple relationship between the tensor view T and vector view x® : if T is stored in row-major order (the default in numpy) then x® = T.flatten() in Python. Marginals. Given a subset A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ } of attributes, a marginal on A can also be represented in 3 different ways: as a set, as a tensor, and as a vector. The set view D A is obtained from D by projecting on the attributes in A (i.e., removing attributes not in A from each record). If A = {} is the empty set, then D { } = {⊥, ⊥, . . . , ⊥} and therefore the only information it retains is the total number of records. The tensor view T A is obtained from D A in the obvious way. That is, T A is created from D A in the same way as T is created from D. Equivalently, if A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ }, then T A is created from T by summing up over all dimensions that are not in {𝑖 1, 𝑖 2, . . . , 𝑖 ℓ }. In the special case of A = {} then T { } = 𝑛 (the number of records). Finally, the vector view x® A is obtained by flattening T A . Once again, when A = {} then x® { } = 𝑛. Linear Queries. A linear query over a marginal can be represented as a multidimensional array q if we want to use it with the tensor view of a marginal. It can also be represented as a row vector q® if we want to use it with the vector view of a marginal. If T A is a marginal with the same shape as q then the query answer is calculated by taking the element-wise product between q and T A and then computing the sum (in the notation of numpy: (q ∗ T A ).sum()). In the vector view q® , we have q® = q.flatten() and the query answer can be computed as the dot product q® · x® A . If all entries of q (or q® ) are in {0, 1}, then we call it a counting query. A set of queries is called a workload and is denoted by WKLoad. For convenience, our notation is summarized in Table 1.

BACKGROUND

Differential Privacy. Differential privacy refers to a family of privacy definitions whose goal is to provide plausible deniability about the presence and contents of any record in a dataset. The most popular versions that support privacy accounting for Gaussian noise are: (𝜖, 𝛿)-DP [7] (also known as approximate DP), Gaussian DP [6], zCDP [4], and

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

Table 1: Table of Notation D: 𝑛: 𝑘: 𝐴𝑖 : 𝑎 (𝑖𝑗 ) : 𝑑𝑖 : T: x® : A: DA: TA : x® A : q: q® : WKLoad: M:

Dataset of records {𝑟 1, . . . , 𝑟 𝑛 }. number of records in D The number of attributes in the data schema. Attribute 𝑖 of the data schema. the 𝑗 th domain element of attribute 𝐴𝑖 . Size of domain of Attribute 𝐴𝑖 . Tensor view of the data. Vector view of the data (column vector) A subset of attributes. Projection of D onto attributes in A. Tensor view of D A . Vector view of D A . Linear query represented as a tensor. Linear query represented as a row vector. A workload of queries. A privacy mechanism.

Renyi DP [23]. QuerySmasher supports all of these definitions, but for simplicity, we focus on approximate DP and Gaussian DP. The goal of these privacy definitions is to add noise into a mechanism’s computation so that any two neighboring datasets D and D ′ (i.e., datasets differing on the presence/absence of 1 record) become nearly indistinguishable. That is, if a mechanism M is run on either D or D ′ , an attacker that observes the output of M would have difficulty in determining which of the two databases was the input. Approximate and Gaussian DP only differ in the mathematical conditions they use to enforce indistinguishability: Definition 1 (Approximate DP [7]). Given privacy parameters 𝜖 ≥ 0 and 𝛿 ∈ [0, 1], a mechanism M satisfies (𝜖, 𝛿)-DP if for all pairs of neighbors D, D ′ and all 𝑆 ⊆ 𝑟𝑎𝑛𝑔𝑒 (M), 𝑃 (M (D) ∈ 𝑆) ≤ 𝑒 𝜖 𝑃 (M (D ′ ) ∈ 𝑆) + 𝛿 Definition 2 (Gaussian DP [6]). Given privacy parameter 𝜇 ≥ 0, a mechanism M satisfies 𝜇-Gaussian DP if for all pairs of neighbors D, D ′ and all 𝑆 ⊆ 𝑟𝑎𝑛𝑔𝑒 (M), Φ−1 (𝑃 (M (D) ∈ 𝑆) ≤ Φ−1 (𝑃 (M (D ′ ) ∈ 𝑆)) + 𝜇 where Φ is the CDF of the standard Gaussian. These privacy definitions can be satisfied by the Correlated Gaussian Linear Mechanism, defined as follows. Definition 3 (Correlated Gaussian Linear Mechanism, Privacy Cost [34]). Given a dataset represented as a vector x® of counts, a Correlated Gaussian Linear Mechanism has the form M (®x) = B®x + 𝑁 (0, Σ), where B is a matrix referred to as the strategy matrix and Σ controls the covariance of the Gaussian noise added to the vector of query answers B®x. The privacy cost of M is defined as max𝑖 (B⊤ Σ−1 B) [𝑖, 𝑖] (i.e., largest diagonal entry of B⊤ Σ−1 B). Computing the approximate DP or Gaussian DP privacy parameters for a collection of correlated Gaussian linear mechanisms is done using the following result: Theorem 1 ([34]). Let M 1, . . . , M 𝑘 be a collection of correlated Gaussian linear mechanisms with respective strategy

matrices B1, . . . , B𝑘 and noise covariance matrices Σ1, . . . , Σ𝑘 . Define the overall privacy cost as 𝑝𝑐𝑜𝑠𝑡 (M 1, . . . , M 𝑘 ) = Í𝑘 ⊤ −1 max𝑖 𝑗=1 B 𝑗 Σ 𝑗 B 𝑗 [𝑖, 𝑖]. If 𝑝𝑐𝑜𝑠𝑡 (M 1 , . . . , M 𝑘 ) ≤ 𝛽 then: • Releasing the outputs of M satisfies (𝜖, 𝛿)√︁ √︁1, . . . , M 𝑘 collectively √︁ √︁ DP with 𝛿 = Φ( 𝛽/2 − 𝜖/ 𝛽) − 𝑒 𝜖 Φ(− 𝛽/2 − 𝜖/ 𝛽). • Releasing the outputs√︁of M 1, . . . , M 𝑘 collectively satisfies 𝜇Gaussian DP with 𝜇 = 𝛽. An important feature of differential privacy and its variations is postprocessing invariance. If M satisfies (𝜖, 𝛿)-DP (resp., 𝜇Gaussian DP) and G is a possibly randomized algorithm whose domain contains the range of M, then the algorithm G ◦ M (which outputs G(M (D))) satisfies (𝜖, 𝛿)-DP (resp., 𝜇-Gaussian DP) with the same values of the privacy parameters. Matrix Mechanisms. The matrix mechanism framework [9, 17] uses the correlated Gaussian linear mechanism1 to answer a workload WKLoad of linear queries. The workload is typically represented as a matrix W, whose rows are the vector representations of queries in WKLoad. The weight of each query q® , denoted as weight(®q), indicates its importance. The selection step chooses matrices L, B, and covariance matrix Σ with the condition that W = LB. The measure step uses the correlated Gaussian linear mechanism determined by the strategy matrix B and covariance matrix Σ to get noisy answers 𝑧® to the strategy queries: 𝑧® = B®x + 𝑁 (0, Σ). The reconstruction phase computes L® 𝑧 , which are noisy answers to the original workload queries W. These noisy answers have the following properties: • 𝐸 [𝐿® 𝑧 ] = W®x, so that the reconstructed answers are unbiased. • Covar[𝐿® 𝑧 ] = LΣL⊤ . Hence the covariance between the 𝑖th query th and 𝑗 query is the (𝑖, 𝑗) entry in this matrix and the (𝑖, 𝑖) entry is the variance of the 𝑖 th query. • The privacy cost is: max𝑖 (B⊤ Σ−1 B) [𝑖, 𝑖]. Hence to minimize the weighted sum squared error of the reconstructed query answers, the selection phase should set L = WB† (the pseudoinverse of B) [17] and then choose B and Σ to minimize the weighted sum of the diagonal entries of LΣL⊤ while (1) keeping the privacy ≤ some user-defined threshold 𝛽 and (2) ensuring the rowspace of B contains the rows of W. When W⊤ W fits in memory, these requirements can be satisfied by solving a semidefinite programming problem using the algorithm of Yuan et al. [35] or the accelerated variation used within HDMM [20].

RELATED WORK

There is a rich and ever-growing literature on differentially private query-answering. Some of the most successful approaches [21] combine ideas from the MWEM [12] and matrix mechanism [17] line of work. MWEM-style algorithms [2, 10–12, 18, 22, 24, 38] maintain a privacy-preserving data synopsis and iteratively identify which queries are poorly answered by it. They obtain noisy answers to those queries (or to alternative strategy queries [10, 24]) and then update the data synopsis. These query answers typically have bias and data-dependent noise, but tend to be very accurate. The matrix mechanism line of work [9, 13, 15, 17, 20, 25, 26, 30– 33, 35, 36] answers workload queries by first finding surrogate 1 There also exist matrix mechanisms that use Laplace noise [17, 36].

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

queries, known as strategy queries, that can be answered with small amounts of noise. Workload query answers are then reconstructed from noisy answers to the strategy queries. They are unbiased and have data-independent variances. Matrix mechanisms are often used inside of other differentially private query-answering systems [1, 10, 19, 22, 24, 26, 39]. Although some matrix mechanisms support pure 𝜖-differential privacy and use Laplace noise [17, 36], better privacy utility trade-offs are obtained with relaxed variants of differential privacy [4, 6, 7, 23] and Gaussian noise. Most work on matrix mechanisms tries to optimize for the weighted sum squared error of workload queries, although optimizing a weighted max [9, 25, 31] is also important in practice. Yuan et al. [35] provide an optimal algorithm for optimizing the weighted sum-squared error in low dimensional settings, while Xiao et al. [31] provide a corresponding algorithm for the weighted max objective. RP [33] and WFF [15] can support a much wider range of convex objectives. Workloads that can be answered from a collection of data marginals, such as marginal and range queries, are important in practice [1] and HDMM was the first scalable matrix mechanism that can optimize for them. However, its runtime and memory requirements were limited by the size of the data domain. RP [33] overcame this limitation by adding noise to carefully chosen strategy queries called the residual basis. The scale of the noise added to the residual basis queries can be tuned to optimize workload accuracy. RP is optimal for marginal queries and scales with the size of the workload rather than the exponentially sized data domain. RP+ [32] is an extension of RP that allows a user to customize which basis should be used, but it does not provide a way of choosing a good basis. WFF [15] proposed an alternative basis that is optimal for marginals and circular product queries (e.g., circular range queries). It computes the multidimensional Fourier transform of data marginals and optimizes for workload accuracy by tuning the variance of the noise added to the Fourier coefficients. Our proposed mechanism QuerySmasher maintains the scalability of [15, 32, 33], while provably dominating them under sum-squared-error on any workload.

OVERVIEW OF QuerySmasher

We first provide an overview of, and intuition behind, QuerySmasher, and then present the full algorithm in Section 6. Afterwards, we show how prior work like RP, RP+, and WFF can be recast as special cases of QuerySmasher, and consequently QuerySmasher dominates them on all workloads under sum-squared error. Let x® be the vector representation of the data, let A ⊆ {𝐴1, . . . , 𝐴𝑘 } be a subset of attributes, and let x® A be the vector representation of the marginal on A (as in Section 2). The vector x® A can (𝑚𝑎𝑟𝑔) be computed from x® by the matrix multiplication: 𝑄 A x® = x® A , where (𝑚𝑎𝑟𝑔)

= V1 ⊗V2 ⊗ · · · ⊗V𝑘 ( I𝑑 if 𝐴𝑖 ∈ A where V𝑖 = 𝑇 𝑖 1𝑑 𝑖 otherwise 𝑄A

and so a linear query q® over x® A (whose true answer is q® · x® A ) can (𝑚𝑎𝑟𝑔) be turned into a linear query q® ′ over x® by setting q® ′ = q®𝑄 A .

(𝑚𝑎𝑟𝑔)

Given two subsets of attributes A and A ′ , the matrices 𝑄 A

(𝑚𝑎𝑟𝑔) and 𝑄 A ′ have linear dependencies between them. So, Xiao et al. proposed the concept of a residual space Resid A on A [33]. (𝑚𝑎𝑟𝑔) Resid A is a linear subspace of the row space of 𝑄 A and any two residual spaces Resid A and Resid A ′ (for A ≠ A ′ ) are mutually

orthogonal. Intuitively, Resid A corresponds to information that is unique to the marginal on A [33]. Formally, Resid A is defined as: Resid A = rowspan (W1 ⊗W2 ⊗ · · · ⊗W𝑘 ) ( if 𝐴𝑖 ∈ A I𝑑 − 1 1𝑑 1𝑇 where W𝑖 = 𝑇 𝑖 𝑑 𝑖 𝑖 𝑑 𝑖 otherwise 1𝑑 𝑖

(1)

RP [33] and WFF [15] proposed specific linearly independent bases for each Resid A and used those bases as the rows in the strategy matrix B of their Gaussian linear mechanisms (Definition 3). In the case of RP+ [32], the user was asked to supply the bases. QuerySmasher takes a different approach: given a workload WKLoad of linear queries that can be answered from a collection of marginals, QuerySmasher breaks apart each q® ∈ WKLoad by projecting it onto the different residual spaces. Subqueries over the same residual space are grouped together into a subworkload. This results in a collection of low-dimensional subworkloads that are small enough to be independently optimized by an exact non-scalable algorithm [37]. The matrix mechanisms created for each subworkload are stitched together to create a matrix mechanism for the original workload. Specifically, the steps are: (1) First, queries in WKLoad are organized by the marginals that can answer them. For example, a query that uses only attributes 𝐴2 and 𝐴8 is answerable from the marginal x® A on A = {𝐴2, 𝐴8 }. (2) Each q® ∈ WKLoad that can be answered by x® A is projected onto the residual spaces Resid A ′ for A ′ ⊆ A. For example, if query q® can be answered from marginal x® {𝐴2 ,𝐴8 } , it will be decomposed into 4 pieces that we denote as q® ⇒{ } , q® ⇒{𝐴2 } , q® ⇒{𝐴8 } , and q® ⇒{𝐴2 ,𝐴8 } . Subquery q® ⇒{ } can be answered from x® { } , subquery q® ⇒{𝐴2 } can be answered from x® {𝐴2 } , etc. The answer to the original query q® can computed by adding up the answers to the subqueries. The details are presented in Section 6. (3) After this decomposition, QuerySmasher creates subworkloads by grouping together subqueries that belong to the same residual space. That is, subworkload WKLoad⇒A consists of all subqueries that result from projections onto the residual space Resid A . For each subworkload WKLoad⇒A , QuerySmasher can use any matrix mechanism construction technique, including the optimal low-dimensional algorithm from [35], to create a matrix mechanism M ⇒A that minimizes the sum-squared error over the subworkload at a privacy cost of 1. (4) Given the mechanisms M ⇒A for different sets A, QuerySmasher rescales the noise that they use, in order to optimize for the global weighted sum squared error of queries in WKLoad subject to a given bound 𝛽 on the privacy cost. Then QuerySmasher runs those mechanisms to obtain noisy answers for the subqueries in each subworkload. (5) To reconstruct an answer to a query q® ∈ WKLoad, 𝑄𝑢𝑒𝑟𝑦𝑆𝑚𝑎𝑠ℎ𝑒𝑟 simply adds up the noisy answers to its subqueries.

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

THE QuerySmasher ALGORITHM

q : yes no

We present the pseudocode for QuerySmasher in the next subsection and then explain the key details in the following subsections.

6.1

The Overall Algorithm

Algorithm 1: Pseudocode for QuerySmasher Input: D, WKLoad, weight for each q® ∈ WKLoad, target privacy cost 𝛽 1 for q ® ∈ WKLoad do 2 A ← attributes used by q® . 3 for A ′ ⊆ A do // Decompose q® 4 Create q® ⇒A ′ using the method in Section 6.2. 5 weight(®q⇒A ′ ) ← weight(®q). 6 Add q® ⇒A ′ to the subworkload WKLoad⇒A ′ . // Optimize subworkloads independently. for each non-empty subworkload WKLoad⇒A ′ do 8 Optimize WKLoad⇒A ′ using a low-dimensional solver. 9 M A ′ is the Gaussian Linear Mechanism provided by the solver. I.e., M A ′ (®x A ′ ) ≡ B A ′ x® A ′ + 𝑁 (0, Σ A ′ ). 10 Recon A ′ is a linear reconstruction function provided by the solver. It can take a query and the output of M A ′ and return an unbiased estimate of the query answer. 11 Var A ′ is the variance function provided by the solver. It takes a query and returns that variance of the reconstructed query answer given by Recon A ′ . 12 L A ′ is the weighted sum square variance of queries in WKLoad⇒A ′ . // Assemble & run overall mechanism for WKLoad. Í √ LA′ A′ 13 𝛾 ← 𝛽 14 for each non-empty subworkload WKLoad⇒A ′ do √︁ 2 = 𝛾/ L ′ . 15 𝜎A ′ A 2 Σ ′) 16 Define M ∗A ′ (®x A ′ ) ≡ B A ′ x® A ′ + 𝑁 (0, 𝜎 A ′ A ∗ 17 𝜔 A ′ ← M A ′ (®x A ′ )// run the mechanism 7

// Reconstruction phase (query answering) for q® ∈ WKLoad do 19 answer(®q) ← 0 and variance(®q) ← 0 20 A ← attributes used by q® . 21 for A ′ ⊆ A do 22 answer(®q) += Recon A ′ (®q⇒A ′ , 𝜔 A ′ ) 23 variance(®q) += Var A ′ (®q⇒A ′ )

a 0 0

b 1 0

c 1 1

q® : [0, 1, 1, 0, 0, 1]

Figure 1: A query answerable by the marginal on A = {𝐴1, 𝐴2 } with domain(𝐴1 )={yes, no} and domain(𝐴2 )={a,b,c}. The tensor form q (left) and vector form q® (right) are shown. The query asks for the number of records that satisfy (𝐴1 = yes, 𝐴2 = b) ∨ (𝐴1 = yes, 𝐴2 = c) ∨ (𝐴1 = no, 𝐴2 = c). The first for-loop, covering Lines 1 through 6, decomposes each query q® into its subqueries by projecting it onto the appropriate residual spaces. If a query accesses a subset A of attributes (i.e., q® is answerable from x® A ), then there is one subquery for each A ′ ⊆ A. Thus, for example, a 2-d range query on income and education level for a specific gender accesses 3 attributes (income, education level, gender) and would be decomposed into 8 subqueries. Then each subquery is assigned to the appropriate subworkload. The weight of each subquery is equal to the weight of the original query. The formula for the decomposition is in Section 6.2. The next for-loop, covering Lines 7 through 12 iterates through the subworkloads. For each subworkload WKLoad⇒A ′ , it uses a low-dimensional algorithm (e.g., [37]) to create a matrix mechanism for the queries in WKLoad⇒A ′ that minimizes their weighted sumsquared error at a privacy cost of 1. Each subworkload can be processed independently and in parallel. Note: data are not accessed in this step. The details are presented in Section 6.3. Also note that algorithms equivalent to RP, RP+, and WFF can be recovered as special cases by using weaker matrix mechanisms than [37] in this step (we explain this in Section 7). A matrix mechanism created in this step consists of the following components (1) the gaussian linear mechanism M A ′ (®x A ′ ) = B A ′ x® A ′ + 𝑁 (0, Σ A ′ ), (2) a linear reconstruction function Recon A ′ that can answer a query q⇒A ′ ∈ WKLoad⇒A ′ if one has the noisy output of M A ′ (®x A ′ ), and (3) a function Var A ′ such that Var A ′ (q⇒A ′ ) is the variance of the reconstructed answer to q⇒A ′ . Neither Recon A ′ nor Var A ′ require access to data, so that one can calculate the variance of any query before seeing the data. Next, the amount of noise used by each M A ′ needs to be rescaled so that the total privacy cost is ≤ 𝛽. Lines 13 and 15 compute the optimal rescaling (Section 6.4) and Line 16 uses it to adjust each mechanism’s noise. Line 17 then runs the mechanisms (this is the only line that accesses data). Lines 18 through 23 reconstruct the answer to every q® ∈ WKLoad and compute its variance by reconstructing the answers to the subqueries of q® (along with their variances) and adding them together.

return reconstructed answer to each query and its variance

6.2 The pseudocode for QuerySmasher is shown in Algorithm 1. Every q® ∈ WKLoad has a weight, denoted as weight(®q). The goal is to provide unbiased noisy answers to these queries while minimizing the weighted sum of their variances: ∑︁ L= weight(®q) ∗ Variance(noisy answer to q® ) (2) q® ∈WKLoad

while keeping the privacy cost bounded by 𝛽.

Query Decomposition Details

We next explain how QuerySmasher decomposes each q® . Suppose a query q® uses the attribute set A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ }. For each A ′ ⊆ A, define the subquery q® ⇒A ′ as: (3) q® ⇒A ′ = q® W𝑖 1 ⊗W𝑖 2 ⊗ · · · ⊗W𝑖 ℓ  1  if 𝑑 𝑖 𝑗 ∉ A ′  𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗  where W𝑖 𝑗 = 1 𝑇  otherwise  I𝑑 𝑖 𝑗 − 𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗 

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

q⇒{ } : q® ⇒{ } : q⇒{𝐴1 } : q® ⇒{𝐴1 } : q⇒{𝐴2 } : q® ⇒{𝐴2 } : q⇒{𝐴1 ,𝐴2 } : q® ⇒{𝐴1 ,𝐴2 } :

. . . . . . . . . . . . . . . . . . . . . . . . . . . . 1/2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1/2 yes 1/6 .................... no -1/6 . . . . . . . . . . . . . . . . . . . . . [1/6, −1/6] a b c ................ -1/2 0 1/2 . . . . . . . . . . . . . . . . . . [−1/2, 0, 1/2] a b c . . . . . yes −1/6 1/3 −1/6 no 1/6 −1/3 1/6 [−1/6, 1/3, −1/6, 1/6, −1/3, 1/6]

Figure 2: The decomposition of q from Figure 1 into the subqueries q⇒{ } , q⇒{𝐴1 } , q⇒{𝐴2 } , and q⇒{𝐴1 ,𝐴2 }

The last two properties in the list above are the reason that reconstruction (Lines 18 through 23 in Algorithm 1) is correct, since subqueries are placed in different subworkloads, and the subworkloads are optimized separately. Theorem 2. Let A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ } be a subset of attributes. Let q® be a query that is answerable from the marginal x® A . For any A ′ ⊆ A, define Q A ′ using Equation 3 so that q® ⇒A ′ = q® Q A ′ . Then q® · x® A = Í q® ⇒A ′ · x® A ′ . Furthermore, if A ′ ≠ A ′′ then Q⊤A ′ Q A ′′ = 0. A′ ⊆ A

Proof. We prove orthogonality first. Represent the matrices in their kronecker product forms: Q A ′ = W1′ ⊗W2′ ⊗ · · · ⊗Wℓ′ and Q A ′′ = W1′′ ⊗W2′′ ⊗ · · · ⊗Wℓ′′ . Then since A ′ ≠ A ′′ , there exists an index 𝑗 such that W′𝑗 ≠ W′′𝑗 . Specifically, one of those matrices must be equal to 𝑑1𝑗 1𝑑 𝑖 𝑗 and the other must be equal to I𝑑 𝑖 𝑗 − 𝑑1𝑗 1𝑑 𝑖 𝑗 1𝑇𝑑𝑖 . 𝑗

In this case, clearly (W′𝑗 )𝑇 W′′𝑗 = 0 and (W′′𝑗 )𝑇 W′𝑗 = 0. Thus, For example, consider the query whose vector representation is shown in Figure 1 (right). The dataset can have many attributes, but the query can be answered from x® {𝐴1 ,𝐴2 } , the vector representation of the marginal on A = {𝐴1, 𝐴2 }. Attribute 𝐴1 has domain size 𝑑 1 = 2 and 𝐴2 has domain size 𝑑 2 = 3. This query is decomposed into 4 subqueries: q® ⇒{ } , q® ⇒{𝐴1 } , q® ⇒{𝐴2 } , and q® ⇒{𝐴1 ,𝐴2 } . To calculate q® ⇒{ } using Equation 3, we multiply q® on the right by ( 21 12 )⊗( 13 13 ), which is a 6-dimensional vector whose entries are all 1/6. Thus q® ⇒{ } is a scalar, the average of the entries in q® . To calculate q® ⇒(𝐴2 ) , we multiply q® on the right by ( 12 12 )⊗(I 3 − h i 2/3 −1/3 −1/3 1 𝑇 ), i.e., 1/2 ⊗ −1/3 2/3 −1/3 , and so on. The four queries 1 1 1/2 3 3 3 −1/3 −1/3

2/3

that q® is decomposed into are shown in Figure 2. Decomposed queries have the following important properties:

Q𝑇A ′ Q A ′′ = ((W1′ )𝑇 W1′′ )⊗ · · · ⊗((W′𝑗 )𝑇 W′′𝑗 )⊗ · · · ⊗((Wℓ′ )𝑇 Wℓ′′ ) = 0 Next, we note that for any A ′ ⊆ A, then x® A ′ can be obtained from x® A by left multiplication by a matrix. Specifically, x® A ′ = V1′ ⊗V2′ ⊗ · · · ⊗Vℓ′ x® A   if 𝐴𝑖 ∉ A ′  1𝑇𝑑𝑖  𝑗 where V′𝑗 =   I𝑑 𝑖 𝑗 otherwise  Let Q A ′ = W1′ ⊗ · · · ⊗Wℓ′ be the matrix that turns q® into q® ⇒A ′ by right multiplication as in Equation 3 and let V A ′ = V1′ ⊗ · · · ⊗Vℓ′ be the matrix that converts x® A into x® A ′ by left multiplication. It is clear that W′𝑗 = 𝑑1𝑖 1𝑑 𝑖 𝑗 if and only if V′𝑗 = 1𝑇𝑑𝑖 . Therefore 𝑗

(1) Each subquery q® ⇒A ′ is answerable from the marginal x® A ′ by taking dot products: q® ⇒A ′ · x® A ′ . (2) If A ′ ≠ A ′′ then the decompositions into q® ⇒A ′ and q® ⇒A ′′ are orthogonal. That is, if q® ⇒A ′ = q® Q1 and q® ⇒A ′′ = q® Q2 then the columns of Q1 are orthogonal to the columns of Q2 . (3) The answer to q® is equal to the sum of the answers of the Í subqueries. That is q® · x® A = q® ⇒A ′ · x® A ′ . A′ ⊆ A

(4) If the subqueries are answered by matrix mechanisms that are independent of each other, the variance of q® equals the sum of the variances of its subqueries. The first and fourth properties are obvious and the other two are proved in the following Theorem 2. The second property means that the information contained in one subquery q® ⇒A ′ is independent of information contained in a different subquery q® ⇒A ′′ . That is, if both are rewritten as query vectors over the full data vector x® , then those query vectors are orthogonal. Since q® ⇒A ′ is later placed into subworkload WKLoad⇒A ′ and q® ⇒A ′′ is placed into WKLoad⇒A ′′ , this means that subworkloads are orthogonal and linearly independent (i.e., if their subqueries are rewritten as query vectors over the full data vector x® , the subqueries in WKLoad⇒A ′ are orthogonal to subqueries in WKLoad⇒A ′′ ). Hence, noisy answers to subqueries in WKLoad⇒A ′ carry no information about noisy answers to subqueries in WKLoad⇒A ′′ and hence subworkloads can be optimized independently of each other.

𝑗

if 𝐴𝑖 𝑗 ∉ A ′

 1 𝑇   𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗  W′𝑗 V′𝑗 =  I𝑑 𝑖 − 1 1𝑑 𝑖 1𝑇   𝑗 𝑑𝑖 𝑗 𝑗 𝑑𝑖 𝑗 Therefore ∑︁

otherwise

q® ⇒A ′ · x® A ′

A′ ⊆ A

! =

∑︁

q® Q A ′ V A ′ x® A = q®

A′ ⊆ A

= q® = q®

∑︁

Q A ′ V A ′ x® A

A′ ⊆ A

∑︁ Ìℓ A′ ⊆ A Ìℓ

𝑗=1

1 𝑇 if 𝐴𝑖 𝑗 ∉A ′ 𝑑 𝑖 1𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗 𝑗 1 otherwise I𝑑 𝑖 − 𝑑 1𝑑 𝑖 1𝑇 𝑖𝑗 𝑗 𝑑𝑖 𝑗 𝑗

I𝑑 𝑖 𝑗 − 𝑑1𝑖 1𝑑 𝑖 𝑗 1𝑇𝑑𝑖 𝑗=1 𝑗 𝑗 𝑇 1 + 𝑑 𝑖 1𝑑 𝑖 𝑗 1𝑑 𝑖 x® A 𝑗

!! x® A

𝑗

= q® x® A where the second-to-last equality is obtained by transforming the Í expression A ′ ⊆ A 𝑓0 (A ′ ) (here 𝑓0 is that Kronecker product on Í the 3rd equation line) into A ′ ⊆ A\{𝐴𝑖 } 𝑓1 (A ′ ), where 𝑓1 (A ′ ) = 1 𝑓 (A ′ ) + 𝑓0 (A ′ ∪ {𝐴𝑖 1 }), then converting that expression into Í0 ′ ′ ′ ′ A ′ ⊆ A\{𝐴𝑖 1 ,𝐴𝑖 2 } 𝑓2 (A ), where 𝑓2 (A ) = 𝑓1 (A ) + 𝑓1 (A ∪ {𝐴𝑖 2 }), etc. □

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

Tensor View of Query Decomposition As some readers find tensor representations of data and queries more intuitive, we explain query decomposition from this viewpoint. Let A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ } be a subset of attributes and let T A be a marginal on those attributes, expressed as a tensor. Let q be a query (in tensor form) that is answerable from T A (e.g., left part of Fig 1). Given an A ′ ⊆ A, we construct q⇒A ′ from q as follows: (1) Add up the tensor representation of q over the dimensions not in A ′ , to get an intermediate tensor X1 . The size of the query Î Î representation is 𝑑 𝑖 𝑗 while the size of X1 is 𝑑𝑖 𝑗 . 𝐴𝑖 ∈ A

𝐴 ∈ A′

𝑖𝑗 𝑗 Î Î (2) Multiply X1 by the scalar ( 𝐴𝑖 𝑗 ∈ A ′ 𝑑 𝑖 𝑗 )/( 𝐴𝑖 𝑗 ∈ A 𝑑 𝑖 𝑗 ), to get the intermediate tensor X2 . Thus X2 is equivalent to averaging q over the dimensions not included in A ′ . (3) q⇒A ′ is obtained from X2 by centering over each dimension.

As an example, consider the query in Figure 1. To obtain q⇒{ } , we sum up q along each dimension to get 3. Since the size (number of cells) of q is 6 and the size of q⇒{ } is 1, we multiply the 3 by 1/6 to get q⇒{ } = 1/2, which is the average of the entries in q. To obtain q⇒{𝐴2 } , we sum up along the dimension for 𝐴1 to get 0 1 2 . This result has size 3 while q has size 6, so we multiply by 3/6 to get 0 1/2 1 . To center this row (so it has mean 0), we need to subtract 1/2 from each entry to get -1/2 0 1/2 , which is equal to the q⇒{𝐴2 } shown in Figure 2. Finally, to obtain q⇒{𝐴1 ,𝐴2 } (i.e., A ′ = A), no summation needs to be performed in Step 1. Step 2 multiplies every entry by 1. Step 3 requires us to iteratively center each row and then each column (one can alternatively center columns first and then rows – this doesn’t change the outcome). We center the first row by subtracting 2/3 from each entry in row 1, and center the second row by subtracting -2/3 1/3 1/3 1/3. This gives . Now each row adds up to 0. To -1/3 -1/3 2/3 center the columns, we add 1/2 to the first column, 0 to the second, -1/6 1/3 -1/6 and subtract 1/2 from the third. This yields , 1/6 -1/3 1/6 which is q⇒{𝐴1 ,𝐴2 } . Note the rows and columns add up to 0. Thus, we see that in the tensor view, we average along the dimensions we eliminate, and then center along each dimension we keep so that the average along those dimensions is 0.

6.3

Optimizing a Workload Partition

Once QuerySmasher has performed an orthogonal decomposition of a workload WKLoad into independent subworkloads WKLoad⇒A ′ , it can optimize each subworkload independently and in parallel. That is, for each WKLoad⇒A ′ , it needs to create a matrix mechanism that can answer the queries in WKLoad⇒A ′ at a privacy cost of 1 while minimizing the weighted sum squared error over WKLoad⇒A ′ . QuerySmasher is free to use any technique in the literature to perform this task, even the optimal non-scalable low-dimensional algorithms [37] for matrix mechanism construction. The reason is that if all queries q® ∈ WKLoad only access a few attributes each, then the decomposed queries q® ⇒A ′ will also access only a few attributes each. The size of the vector representation of

Î q® ⇒A ′ is 𝐴𝑖 ∈ A ′ |𝐴𝑖 | and this is manageable when A ′ is small (otherwise, one can add approximation algorithms to QuerySmasher as discussed below). Given a subworkload WKLoad⇒A ′ that can be answered from a marginal x® A ′ , a matrix mechanism constructor [17] produces: • A gaussian linear mechanism M (®x A ′ ) ≡ B®x A ′ + 𝑁 (0, Σ). This mechanism is run once to get a noisy output 𝑧®. • Recon: a linear function such that Recon(®q, 𝑧®) is a noisy answer to q® ∈ WKLoad A ′ . • Var: a function such that Var(®q) is the variance of the reconstructed answer (i.e. Variance(Recon(®q, 𝑧®)). This variance is dataindependent hence Var does not need any data access. Optimal Matrix Mechanism Construction. Let W be the matrix whose rows are the query vectors in WKLoad⇒A ′ . Let D be a diagonal matrix, where entry D[𝑖, 𝑖] is the weight of the query that is in the 𝑖 th row of W. When the goal is to minimize the weighted sum squared error of the queries, Yuan et al. [37] formulated matrix mechanism construction in terms of the following semidefinite optimization problem (where † denotes pseudoinverse): arg min trace(W⊤ DWV† )

(4)

s.t. V is symmetric positive semidefinite the row space of V contains the row space of W and V[𝑖, 𝑖] ≤ 1 for all 𝑖 which can be solved using the algorithm of Yuan et al. [37] or a simplified version by McKenna et al. [20] or even off-the-shelf semidefinite program solvers. Once a solution V is obtained, the rest of the steps are: (1) Compute the Cholesky decomposition V = B𝑇 B and returns B. Optionally, project the rows of B onto the rowspace of W. (2) Given a desired privacy cost 𝛽, set 𝜎 2 = max𝑖 (B⊤ B) [𝑖, 𝑖]/𝛽. Set M (®x A ′ ) ≡ B®x A ′ + 𝑁 (0, 𝜎 2 I). (3) Define Recon(®q, 𝑧®) ≡ q® B†𝑧®, where B† is the pseudoinverse. (4) Define Var(®q) ≡ 𝜎 2 q® B† (B† ) ⊤ q® ⊤ (recall q® is a row vector). Approximate Matrix Mechanism Construction. Sometimes the semidefinite program in Equation 4 is difficult to solve because the main quantity W⊤ DW is very large. This can happen when A ′ has many variables (i.e., some queries do access many variables) or when some attributes 𝐴𝑖 ∈ A ′ already have large domains. The flexibility of QuerySmasher means that in those situations, approximate algorithms for creating matrix mechanisms can be added to QuerySmasher to handle those cases. Common approaches include (1) not running the solver for Equation 4 to completion, or (2) imposing additional structure on the optimization variables V. For example, one common approach used in practice is to fix a basis U whose rowspan includes the rows of W, and require that V = U𝑇 HU for some diagonal matrix H [16, 17, 30]. This results in fewer optimization parameters (the diagonal entries instead of a full matrix) and can often be solved in close form. This flexibility and modularity is responsible for the generality of QuerySmasher. In fact, we show in Section 7 that using approximate (instead of optimal) techniques can make QuerySmasher equivalent to RP, RP+, and WFF.

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

6.4

Assembling the Overall Solution.

After QuerySmasher is done constructing a matrix mechanism for each nonempty subworkload WKLoad⇒A ′ , it ends up with a linear gaussian mechanism for each A ′ , having the form M A ′ (®x A ′ ) ≡ B A ′ x® A ′ + 𝑁 (0, Σ A ′ ). Each such mechanism has privacy cost 1. Let 𝐿 A ′ be the weighted sum squared error of the matrix mechanism over the subworkload WKLoad⇒A ′ . If QuerySmasher rescales the noise used by M A ′ by some value 2 (i.e., redefining it as M ′ (® 2 Σ ′ )) 𝜎A ® A ′ + 𝑁 (0, 𝜎 A ′ ′ A A xA′ ) ≡ BA′ x 2 then the privacy cost becomes 1/𝜎 A ′ and the weighted sum-squared 2 𝐿 ′. error over the subworkload becomes 𝜎 A ′ A 2 , for each A ′ , so Thus, QuerySmasher searches for values of 𝜎 A ′ that the overall privacy cost is 𝛽, while making the total weighted sum squared error as small as possible. That is, it needs to solve the following optimization problem: ∑︁ ∑︁ 1 2 =𝛽 min 𝜎A s.t., ′ 𝐿 A′ 2 𝜎2 𝜎A ′ A′ A′ A′ The form of this optimization problem was studied in [32] (Lemma 2 are: 3) and the optimal choices of the 𝜎 A ′ Í √ 1 A ′′ 𝐿 A ′′ 2 𝜎A √ ′ = 𝛽 𝐿 A′ These scale factors are computed by Algorithm 1 in Line 15. This discussion proves the privacy properties of QuerySmasher: Theorem 3. The privacy cost of Algorithm 1 is ≤ the input parameter 𝛽.

Definition 4. A gaussian linear mechanism M respects residual spaces if there exists an A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ } such that: • M (®x A ) ≡ B®x A +𝑁 (0, Σ) for some B and Σ (i.e., M adds Gaussian noise to linear queries computable from the marginal on A). • Each row of B belongs to the row space of 1 1 I𝑑 𝑖 1 − 1𝑑 𝑖1 1𝑑⊤𝑖 ⊗ · · · ⊗ I𝑑 𝑖ℓ − 1𝑑 𝑖ℓ 1𝑑⊤𝑖 . ℓ 1 𝑑 𝑖1 𝑑 𝑖ℓ Intuitively, this means that rows of B belong to the same linear space as queries in the chosen subworkload WKLoad⇒A and hence M could be used to provide answers to them. Example 4. ResidualPlanner (RP) [33] uses gaussian linear mechanisms of the form: M (®x A ) ≡ B®x A + 𝑁 (0, Σ), where the matrix B = Sub𝑑 𝑖1 ⊗ · · · ⊗Sub𝑑 𝑖ℓ , where Sub𝑑 is a (𝑑 − 1) × 𝑑 matrix that can be written in block matrix form as [1𝑑 −I𝑑 ] (the first column is full of ones, entries of the form (𝑖, 𝑖 + 1) contain -1, and every 0 . It is easy to verify thing else is 0). For example Sub3 = 11 −1 0 −1 that the row space of each Sub𝑑 𝑖 𝑗 is a subspace of the rowspace of I𝑑 𝑖 𝑗 − 𝑑1𝑖 1𝑑 𝑖 𝑗 1𝑑⊤𝑖 and so B satisfies the condition of Definition 𝑗

Example 5. ResidualPlanner+ (RP+) [32] uses gaussian linear mechanisms of the form: M (®x A ) ≡ B®x A + 𝑁 (0, Σ), where the matrix B = S1 ⊗ · · · ⊗Sℓ , where each S𝑖 depends on a user’s input and is modified so that its row space belongs to the row space of I𝑑 𝑖 𝑗 − 𝑑1𝑖 1𝑑 𝑖 𝑗 1𝑑⊤𝑖 . Hence RP+ also respects residual spaces. 𝑗

After these calculations have been performed, the measure phase (which is the only step that accesses raw data, and hence the only step that incurs the privacy cost) can begin. QuerySmasher runs all of the linear gaussian mechanisms M A ′ (®x A ′ ) = B A ′ x® A ′ + 2 Σ ′ ) with this appropriately scaled noise. It collects their 𝑁 (0, 𝜎 A ′ A respective noisy answers 𝑧®A ′ and does not access the data anymore. To answer a query q® ∈ WKLoad, QuerySmasher re-uses the noisy answers it has collected. If q® accesses attributes in the set A, then QuerySmasher decomposes it into subqueries q® ⇒A ′ for all A ′ ⊆ A using Equation 3 from Section 6.2. It reconstructs an answer to each q® ⇒A ′ using the reconstruction function Recon A ′ provided by the matrix mechanism for the subworkload WKLoad⇒A ′ (i.e., QuerySmasher computes Recon A ′ (®q⇒A ′ , 𝑧®A ′ )). Since the answer to q® is the sum of the answers of its subqueries (Theorem 2), the answer to q® is calculated as: ∑︁ answer(®q) = Recon A ′ (®q⇒A ′ , 𝑧®A ′ ) A′ ⊆ A

ACCURACY OF QuerySmasher

Next we show that algorithms equivalent to RP, RP+, and WFF can be obtained as special cases when QuerySmasher uses approximate rather than optimal matrix mechanism construction techniques, which implies that QuerySmasher would dominate those methods on every workload, if it uses optimal solvers for the subworkloads. We do this by abstracting some common linear algebraic properties that they share. The first property is respecting residual spaces, which we define as follows:

𝑗

4. Therefore every gaussian linear mechanism used by RP respects residual spaces.

𝑗

Example 6. Weighted Fourier Factorization (WFF) [15] applies the multidimensional Fourier transform to different data marginals and adds independent noise to their coefficients. This is equivalent to creating linear mechanisms of the form M (®x A ) ≡ B®x A +𝑁 (0, Σ) where Σ is a diagonal matrix and each row of B is created by taking the real and imaginary components of vectors of the form 𝑣®1 ⊗® 𝑣 2 ⊗ · · · ⊗® 𝑣 ℓ , where each 𝑣®𝑗 has size 𝑑 𝑖 𝑗 and has the form: h i𝑡 −2𝜋 i𝑡 𝑣®𝑗 = exp −2𝜋 𝑑 𝑖 𝑗 · 0 , exp 𝑑 𝑖 𝑗 · 1 , i i𝑡 . . . , exp −2𝜋 𝑑 𝑖 𝑗 · (𝑑 𝑖 𝑗 − 1) √ where i = −1 and 𝑡 is an integer between 1 and |𝑑 𝑖 𝑗 | − 1. Varying the values of 𝑡 inside each of the 𝑣®𝑗 gives the different rows of B. Although this is a bit complex, one can see that each 𝑣®𝑗 is orthogonal to 1𝑑 𝑖 𝑗 and therefore B satisfies the conditions of Definition 4. The next condition is that if a matrix mechanism is given an overall target privacy cost of 𝛽, then the privacy costs of the gaussian linear mechanisms M 1, . . . , M𝑚 it uses add up to 𝛽. It is not a vacuous condition because if M 𝑖 uses strategy matrix B𝑖 and covariance matrix Σ𝑖 then its privacy cost is max 𝑗 B𝑖⊤ Σ𝑖−1 B𝑖 ) [ 𝑗, 𝑗] (Definition 3). However, the overall privacy cost (Theorem 1) is: ! 𝑚 ∑︁ ⊤ −1 𝛽 = max B𝑖 Σ𝑖 B𝑖 [ 𝑗, 𝑗] 𝑗 𝑖=1

≤

𝑚 ∑︁

max B𝑖⊤ Σ𝑖−1 B𝑖 [ 𝑗, 𝑗] 𝑗

𝑖=1

(5)

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

𝑚 ∑︁

𝑝𝑐𝑜𝑠𝑡 (M 𝑖 )

𝑖=1

Thus the overall privacy cost, in general, can be smaller than the sum of the individual privacy costs of the mechanisms. Definition 5. A matrix mechanism is splittable if, when given a target privacy cost bound 𝛽, it creates a collection of gaussian linear mechanisms and the sum of their individual privacy costs is 𝛽. The RP+ algorithm is explicitly designed to be splittable. The linear mechanisms used by RP all have constant diagonals for B𝑖⊤ Σ−1 B𝑖 and so clearly Equation 5 holds with equality, so RP is splittable. When WFF adds noise to Fourier coefficients, the complex part of the coefficient gets the same variance as the real part of the coefficient. As a result, the privacy cost matrix of each gaussian linear mechanism it uses has a constant diagonal and therefore, like RP, it is splittable. Our final theoretical result is that if a competing matrix mechanism is splittable and every gaussian linear mechanism it uses respects residual spaces then the weighted sum squared error over WKLoad that it achieves cannot be smaller than the weighted sum squared error of QuerySmasher (provided that QuerySmasher uses an optimal solver for Equation 4). In particular each gaussian linear mechanisms it uses is a feasible (and potentially suboptimal) solution to the matrix mechanism construction semidefinite program (Equation 4) for some subworkload WKLoad⇒A . Theorem 7. Given a workload WKLoad of linear queries that can be answered from data marginals, and a privacy cost bound 𝛽, let Compete denote a matrix mechanism that can answer the queries in WKLoad, is splittable, and which uses gaussian linear mechanisms that all respect residual spaces. Then: • The weighted sum-squared error of Compete over WKLoad is equal to the sum of the weighted sum-squared errors of Compete on the subworkloads created by QuerySmasher. • For any given A, let M 1, . . . , M𝑚 be the subset of Gaussian linear mechanisms of Compete that are answerable from the marginal x® A . Then M 1, . . . , M𝑚 together correspond to a feasible (but potentially suboptimal) solution to the optimization problem in Equation 4 for subworkload WKLoad⇒A (i.e., the W in that equation is formed from the query vectors in WKLoad⇒A ). • When QuerySmasher uses solvers that solve Equation 4 to optimality, the weighted sum squared error of QuerySmasher over WKLoad is ≤ the weighted sum-squared error of Compete over WKLoad. Proof sketch. (see the appendix for the full proof). The main steps are to show that for any non-empty subworkload generated by QuerySmasher, a subset of the strategy queries used by Compete form a linear basis for the space spanned by that subworkload. Furthermore, the strategy queries from Compete that are assignable to one subworkload are linearly independent from those assignable to a different subworkload. Then, like with Theorem 2, answering a query q® using Compete is equivalent to using Compete to answer the subqueries that q® is decomposed into, and adding up the answers. Since QuerySmasher is optimal for the subqueries, it then has lower error than Compete. □

EXPERIMENTS

In this section, we present our experimental results. Baselines. We compare QuerySmasher (QS) to ResidualPlanner (RP) [33], ResidualPlanner+ (RP+) [32], HDMM [20], and Weighted Fourier Factorization (WFF) [15]. Public implementations of prior work do not support all of the query workloads we consider, so we take the following approach. RP only supports marginals but is completely subsumed by RP+ (hence comparison to RP can be dropped). Furthermore, we modifed the RP+ implementation to support a wider set of queries. WFF does not have a public implementation so we re-implemented it (a practical application of Theorem 7). We did not modify HDMM code and so query workloads that are not supported by the implementation are marked as NI (not implemented). Implementation details. We implemented QuerySmasher in Python 3.11 using NumPy, SciPy, and PyTorch. All computations are CPU-only; our experiments are run on a dual-socket AMD EPYC 7763 server (128 cores, 1 TB RAM) running Linux and the residual spaces are optimized sequentially.

8.1

Workload queries.

We evaluate the matrix mechanisms on a wider variety of queries than prior work. Specifically, we consider the following: • Marginal queries. One-way marginals have the form “how many records have 𝐴𝑖 = 𝑐” (for each choice of c), two-way marginals have the form “how many records have 𝐴𝑖 = 𝑐 1 and 𝐴 𝑗 = 𝑐 2 ”, etc. Note that RP, RP+, WFF, and QuerySmasher are all optimal for pure marginal workloads, so this serves as a sanity check. • Prefix-sum queries: these are queries of the form “how many records have 𝐴𝑖 ≤ 𝑐” for different values of 𝑐. Two-dimensional prefix-sum queries have the form “how many records have 𝐴𝑖 ≤ 𝑐 1 and 𝐴 𝑗 ≤ 𝑐 2 ,” etc. • Range queries: these are axis-aligned queries of the form “how many records have 𝐴𝑖 ∈ [𝑐 1, 𝑐 2 ]” and are extended to multiple dimensions in the obvious way (e.g., “how many records have 𝐴𝑖 ∈ [𝑐 1, 𝑐 2 ] and 𝐴 𝑗 ∈ [𝑐 3, 𝑐 4 ]”). • Circular range queries: these are range queries that wrap around like on a clock (e.g., “how many checkins are between “10pm and 1am”) and are extended to multiple dimensions in the obvious way. WFF [15] is provably optimal for these queries and therefore so is QuerySmasher. • Affine queries: these are two-dimensional queries that are not axis-aligned. They are answerable from two-way marginals and have the form “how many records have 𝐴𝑖 + 𝐴 𝑗 ≤ 𝑐” for different choices of 𝑐. These queries allow for cross-attribute comparisons, which prior work did not consider. • Abs queries: these are also two-dimensional queries that are not axis-aligned. They allow cross-attribute comparisons and have the form “how many records have |𝐴𝑖 − 𝐴 𝑗 | ≤ 𝑐” for different choices of 𝑐. • Random: finally we consider random counting queries that are formed by taking a query vector q® , initially filled with zeros, and then independently flipping each bit with probability 𝑝. This query type simulates queries without any particular structure. Noting that there is no universally accepted choice of 𝑝, we chose 𝑝 = 0.3 to get somewhat sparse query vectors.

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

8.2

Metrics

In our experiments, each matrix mechanism is required to optimize its workload at a privacy cost 𝛽 = 1. We measure the root-meansquared error:

v u t 𝑅𝑀𝑆𝐸 =

8.3

Var(®q)

q® ∈WKLoad

|WKLoad|

Executive Summary

For the accuracy experiments, we evaluate all methods on real data in Section 8.4 and provide more extensive experiments on synthetic data in Section 8.5. Scaling experiments are provided in Section 8.6. The high-level summary of our observations from those sections are: • Among the baselines (RP+, WFF, HDMM), each has its own specialty: HDMM performs very well on marginal queries, RP+ performs well on marginals and prefix sums. WFF is provably optimal for marginals and circular range queries [15]; it also performs well for range queries and random queries (but struggles with prefix queries). • QuerySmasher consistently has the lowest RMSE for all workloads. It matches WFF in places where it is optimal (marginals and circular range queries) and is strictly better than all methods on the rest of the queries. • None of the competitors (RP+, WFF, HDMM) perform well on affine and abs queries, and QuerySmasher shows significant improvements for those query types. It also outperforms HDMM and WFF on prefix queries by a large margin and outperforms RP+ on prefix queries by a smaller margin. • In terms of scaling, the dominant cost arises in optimizing a subworkload, as deploying the optimal solver requires solving a semi-definite optimization program. The resource requirements only grow with the size of the output (e.g., size of the workload query answers) rather than with the domain size (which is exponentially large). To put it into context, consider the following two results from Section 8.6. – Consider a complex workload consisting of all 1-d range queries, all 2-d affine queries, and all 3-d prefix queries. Consider a synthetic dataset with 10 attributes, with each attribute having a domain size of 40. Thus, the overall domain size of a tuple is 1040 . For the workload on such a dataset, there are 10 prefix queries (i.e., 120 ways of choos3 = 120 groups of 3-d ing 3 attributes), 10 = 45 groups of 2-d affine queries, and 2 10 = 10 groups of 1-d range queries. Hence, answering this 1 workload, even in the non-private setting, requires computing 120+45+10 = 175 views of the dataset. For QuerySmasher the solver time was 60 seconds and required 3.4 GB of memory. – If the number of dimensions is increased to 40, the overall domain size of a tuple becomes 4040 ≈ 1.2 × 1064 . Now there are 40 prefix queries (i.e., 9880 ways of 3 = 9880 groups of 3-d choosing 3 attributes), 40 = 780 groups of 2-d affine queries, 2 and 40 1 = 40 groups of 1-d range queries. Hence answering this workload requires computing 9880 + 780 + 40 = 10700

views of the data. This is now a 10700 175 ≈ 61 times larger problem. The solve time is 3032 seconds (which is only 3032 60 ≈ 51 times larger) and the memory requirement is 227 GB (which is only 227 3.4 ≈ 67 times larger). We note that in addition to parallelizing the solve phase of QuerySmasher, one can also replace the optimal low-dimensional subworkload optimizer (Section 6.3) with an approximate optimizer. Hence, any future research improvements in lowdimensional matrix mechanisms will directly benefit QuerySmasher as well. We next consider the experiments in detail.

8.4

Utility Experiments for Real Data

Following prior work [20, 32], we consider the following real dataset schemas. Adult (9 categorical attributes with domain sizes 42, 16, 15, 9, 7, 6, 5, 2, 2 and 5 numeric attributes with domain sizes 100, 100, 100, 99, 85) [3], CPS (3 categorical attributes with domain sizes 7, 4, 2 and 2 numeric attributes with sizes 50, 100) [20], and Loans (8 categorical attributes with domain sizes 51, 36, 15, 8, 6, 5, 4, 3 and 4 numeric attributes with sizes 101, 101, 101, 101) [14]. Following prior work, our first experiment uses the hybrid marginal/prefix workload. For example, a 1-way hybrid query on attribute 𝐴𝑖 has the form “how many records have 𝐴𝑖 op1 𝑐 1 ”, where op1 is “=” if 𝐴𝑖 is categorical and “<” if 𝐴𝑖 is numeric. Then 2-way hybrid queries have the form “how many records have 𝐴𝑖 op1 𝑐 1 and 𝐴 𝑗 op2 𝑐 1 ” and so on. The RMSE of the different methods is shown in Table 2. Table 2: RMSE for hybrid marginal/prefix workloads at privacy cost 1. Dataset Adult

CPS

Loans

Workload 1-way 2-way 3-way ≤ 3-way 1-way 2-way 3-way ≤ 3-way 1-way 2-way 3-way ≤ 3-way

#Queries 588 148,137 20,894,536 21,043,261 163 7,000 72,556 79,719 532 118,974 14,539,522 14,659,028

RP+ 5.11 17.63 47.19 48.90 3.18 6.36 8.12 8.39 4.73 14.91 36.11 36.65

HDMM 5.648 21.111 60.025 61.613 3.347 6.625 8.634 8.623 5.127 17.149 42.938 44.338

WFF 5.804 22.493 66.934 66.794 3.589 7.841 10.267 10.594 5.330 18.736 49.823 49.721

QS 5.047 17.632 47.055 47.853 3.135 6.194 7.903 8.140 4.670 14.822 36.095 36.410

This is a relatively easy workload, with QuerySmasher consistently having the lowest RMSE, RP+ being slightly worse, while HDMM and WFF struggle with 3-way queries on the larger dataset schemas (Adult and Loans). To accommodate additional query classes, we consider an experiment where all the attributes are treated as numeric and we evaluate the RMSE on all 2-way affine and 2-way Abs queries. The results are shown in Table 3. We see from Table 3 that QuerySmasher consistently outperforms the competing methods while reducing the noise between 15 − 21% compared to the next best method. As these are relatively small datasets, we next consider more systematic evaluations on synthetic data.

8.5

Utility Experiments for Synthetic Data

We present consider synthetic datasets whose schemas have the form [𝑛]𝑑 . That is, there are 𝑑 attributes and each attribute has

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

Table 3: RMSE for 2-way Affine and abs-diff workloads. QuerySmasher improves RMSE by 15 − 21%. Dataset CPS Adult Loans

Workload affine abs_diff affine abs_diff affine abs_diff

#Queries 805 731 99,830 99,200 66,318 65,682

RP+ 10.012 10.476 19.395 19.441 18.082 18.148

HDMM NI NI NI NI NI NI

WFF 7.147 7.387 21.559 21.587 18.671 18.714

QS 5.935 5.900 16.435 16.254 14.305 14.330

a domain size of 𝑛. Thus the overall domain of a tuple is 𝑛𝑑 . We evaluate all of the workload types as 𝑛 and 𝑑 vary. Table 4: RMSE of 1-way + 2-way workloads at privacy cost= 1, for different query types, on a synthetic dataset with 40 attributes (𝑑 = 40) as 𝑛 (domain size of an attribute) varies. NI indicates that competing methods do not have an implementation for the query type. For HDMM, we added a * to indicate that it can estimate its RMSE but will run out of memory when providing workload query answers. 𝑛 10

Workload marginal prefix range circular affine abs random marginal prefix range circular affine abs random marginal prefix range circular affine abs random marginal prefix range circular affine abs random marginal prefix range circular affine abs random

#Queries 78,400 78,400 2,361,700 7,804,000 15,220 8,200 235,200 312,800 312,800 34,406,400 124,816,000 31,220 16,400 938,400 703,200 703,200 168,674,100 631,836,000 47,220 24,600 2,109,600 1,249,600 1,249,600 524,504,800 1,996,864,000 63,220 32,800 3,748,800 1,952,000 1,952,000 1,268,038,500 4,875,100,000 79,220 41,000 5,856,000

RP+ 23.48 36.23 48.95 NI 62.48 85.99 NI 25.70 52.63 75.44 NI 110.30 152.41 NI 26.46 63.10 93.22 NI 149.68 207.22 NI 26.84 70.95 106.85 NI 184.38 255.57 NI 27.07 77.29 118.00 NI 215.98 299.62 NI

HDMM 23.48∗ 52.35∗ 68.98∗ NI NI NI NI 25.70∗ 80.01∗ 109.56∗ NI NI NI NI 26.46∗ 98.35∗ 137.92∗ NI NI NI NI 26.84∗ 112.33∗ 160.16∗ NI NI NI NI 27.07∗ 123.73∗ 178.65∗ NI NI NI NI

WFF 23.48 39.70 41.36 39.77 45.23 64.11 107.81 25.70 62.95 63.58 63.01 70.31 102.54 234.89 26.46 79.09 79.02 79.14 87.62 129.31 362.86 26.84 91.67 91.11 91.72 101.10 150.21 490.66 27.07 102.09 101.15 102.13 112.23 167.51 618.61

QS 23.48 33.70 41.08 39.77 28.25 35.85 104.43 25.70 49.51 63.32 63.01 35.71 39.49 233.69 26.46 60.81 78.79 79.14 44.36 48.14 362.01 26.84 68.78 90.91 91.72 69.62 49.83 490.65 27.07 75.26 100.97 102.13 79.33 52.80 616.54

Table 4 considers the RMSE of 1-way and 2-way workloads, while Table 5 considers 1-way + 2-way + 3-way workloads of mixed

query types. We note that HDMM and RP+ not have implementation support for some of the workloads, and we mark such cases in a table cell as NI (not implemented). We also note that HDMM cannot fully run on datasets of this scale – it runs out of memory before it can reconstruct noisy workload query answers (since it scales as 𝑂 (𝑛𝑑 )). Nevertheless, it can calculate the RMSE it would have in the ideal case of unbounded computational resources. Hence we add a * to cell entries where the RMSE is predictable but the algorithm cannot run to completion. Table 4 shows that QuerySmasher has the lowest RMSE among all methods across all workloads. Among the competitors, WFF is provably optimal (among matrix mechanisms) for range queries and circular queries [15] and is nearly optimal for range queries. Meanwhile, RP+ has strong performance for prefix queries (as well as being optimal for marginals [32]). QuerySmasher matches the optimal performance on marginals and circular queries, and is strictly better on all of the other query types. The largest sources of improvements are for Abs and Affine. Table 5: RMSE on mixed-query workloads of 1-way range queries, 2-way affine queries, 3-way prefix queries at privacy cost 1 as 𝑑 (number of dataset attributes) and 𝑛 (domain size of each attribute) vary. 𝑛 10

𝑑 10 20 30 40 50 10 20 30 40 50 10 20 30 40 50

RP+ 22.49 57.17 99.66 148.46 202.72 36.82 97.22 172.22 259.00 355.93 47.08 126.20 224.99 339.60 467.86

HDMM NI NI NI NI NI NI NI NI NI NI NI NI NI NI NI

WFF 25.37 65.13 113.98 170.18 232.72 47.04 125.94 224.36 338.53 466.26 64.37 175.18 314.28 476.12 657.56

QS 20.41 51.63 93.50 138.38 187.24 34.60 95.63 167.16 249.29 340.55 44.46 126.19 221.80 331.86 454.37

Table 5 presents a more challenging workload that combines all 1-way range queries, all 2-way affine queries, and all 3-way prefix-sum queries. As the data dimensionality 𝑑 and attribute size 𝑛 vary, QuerySmasher shows consistent and large improvements in RMSE.

8.6

Scalability Experiments

QuerySmasher has a query decomposition phase, a subworkload optimization phase, and the assembly phase which computes privacypreserving queries over the data and then reconstructs the workload query answers. We next test the scalability of these phases with results shown in Tables 6 and 7. In these experiments, we consider datasets with 𝑑 attributes, where each attribute has a domain size 𝑛 (we vary 𝑑 and 𝑛). Hence, the tuple domain is 𝑛𝑑 . Table 6 considers the easier workload that consists of all 1-way and 2-way prefix sum queries. Table 7 considers a more challenging workload that consists of all 1-way range queries, all 2-way affine queries, and all 3-way prefix sum queries. The tables report the running time in seconds for the decomposition (Decomp.), subworkload optimization (Solve) and the assembly

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

(Assem.) phases. They also report memory usage. To help understand scaling behavior, we also included a #Views column which refers to the number of data views that would be necessary for answering the workload even in the non-private case. For example, when the number of attributes 𝑑 = 40, then there are 40 3 = 9880 types of three-way prefix-sum queries (i.e., 9880 ways of choosing 3 40 attributes), 40 2 = 780 types of two-way affine queries and 1 = 40 types of one-way range queries. Hence 9880 + 780 + 40 = 10700 dataset views (group-by count queries) on those different attributes are needed to answer the workload even in the non-private case. In Tables 6 and 7 we see that the resources (memory and running time) scale well with #Views. For example, consider the two rows in Table 7 corresponding to 𝑛 = 40, 𝑑 = 10 and 𝑛 = 40, 𝑑 = 40. The number of attributes increases by a factor of 4 but the number of views increases from 175 to 10,700 — a factor of ≈ 61×. Meanwhile, the solve time increases from 59.6 seconds to 3031.6 (a factor of ≈ 51×) and the memory requirements increase from 3.4 GB to 227.4 GB (a factor of ≈ 67×). Table 6 exhibits substantially lower resource usage because it only needs support solve 1-way and 2-way workloads. To get an idea of the necessary memory requirements, suppose ℓ queries in a workload can be answered using a marginal with 𝑐 total cells. A single linear query over that marginal is a linear combination of the values in the 𝑐 cells and hence is represented as a 𝑐-dimensional vector. Materializing the workload then requires ℓ𝑐 space (typically ℓ ≥ 𝑐). An optimal matrix mechanism for the part of the workload answerable by this marginal would generally ask for noisy answers to 𝑐 linear strategy queries, hence requiring 𝑂 (𝑐 2 ) space. This is a fundamental bottleneck that is independent of QuerySmasher. The only way around it would be to use sub-optimal matrix mechanisms that are restricted to strategy queries that have a sparse representation and/or to cache optimal solutions (so that marginals with the same shape can re-use the solution). The modular design of QuerySmasher would allow optimal solvers to be replaced with sub-optimal solvers to save resources.

Table 6: Timing experiments on 1-way + 2-way prefix queries as 𝑑 (number of attributes in the dataset) and 𝑛 (domain size of each attribute) vary. #Views refers to the number of data views that would be needed to answer the queries even in the non-private setting.

𝑛 10

Acknowledgments. This work was partially supported by NSF awards CNS-2349610, CNS-2317232, and CNS-2317233

REFERENCES [1] John M. Abowd, Robert Ashmead, Ryan Cumings-Menon, Simson Garfinkel, Micah Heineck, Christine Heiss, Robert Johns, Daniel Kifer, Philip Leclerc, Ashwin Machanavajjhala, Brett Moran, William Sexton, Matthew Spence, and Pavel Zhuravlev. forthcoming. Preprint https://www.census.gov/library/working-

#Views 55 210 465 820 1275 55 210 465 820 1275 55 210 465 820 1275 55 210 465 820 1275

QuerySmasher Phase Time Decomp Solve Assem <0.1s 2.8s <0.1s <0.1s 4.4s <0.1s <0.1s 2.6s <0.1s 0.1s 7.1s <0.1s 0.2s 9.5s 0.1s <0.1s 2.7s <0.1s <0.1s 4.2s <0.1s <0.1s 6.4s <0.1s 0.1s 8.1s <0.1s 0.1s 10.1s <0.1s <0.1s 3.6s <0.1s <0.1s 7.0s <0.1s 0.1s 9.9s <0.1s 0.1s 12.7s <0.1s 0.1s 16.4s 0.1s <0.1s 4.9s <0.1s <0.1s 9.3s <0.1s 0.1s 14.8s <0.1s 0.1s 20.5s <0.1s 0.2s 37.5s 0.1s

Total 2.8s 4.4s 2.6s 7.2s 9.7s 2.7s 4.3s 6.4s 8.2s 10.2s 3.6s 7.1s 10.0s 12.8s 16.6s 4.9s 9.3s 14.9s 20.6s 37.7s

Mem <1 GB <1 GB <1 GB <1 GB <1 GB <1 GB <1 GB <1 GB <1 GB 1 GB 1 GB 1 GB 1 GB 1 GB 1 GB 1 GB 1 GB 1 GB 1 GB 2 GB

Table 7: Timing experiments on mixed workloads containing all 1-way range queries, all 2-way affine queries, and all 3way prefix sum queries, as 𝑑 (number of attributes in the dataset) and 𝑛 (domain size of each attribute) vary. #Views refers to the number of data views that would be needed to answer the queries even in the non-private setting.

𝑛 10

CONCLUSIONS AND FUTURE WORK

In this paper, we introduced QuerySmasher, a matrix mechanism that uses a divide-and-conquer strategy to answer linear workloads that are computable from collections of marginals. This approach allows optimal matrix mechanism constructors to be re-used on small sub-problems. Future research includes improving optimal and approximate low-dimensional matrix mechanisms as the resource savings would automatically transfer to QuerySmasher. Other areas of future research include integration with MWEM style algorithms. Such an integration would combine their respective strengths: the ability to plan strategy queries with data-dependent exploration of which queries need more accuracy.

𝑑 10 20 30 40 50 10 20 30 40 50 10 20 30 40 50 10 20 30 40 50

𝑑 10 20 30 40 50 10 20 30 40 50 10 20 30 40 50 10 20 30 40

#Views 175 1350 4525 10700 20875 175 1350 4525 10700 20875 175 1350 4525 10700 20875 175 1350 4525 10700

QuerySmasher Phase Time Decomp Solve Assem <0.1s 4.5s <0.1s 0.5s 28.1s 0.8s 1.2s 25.6s 2.6s 2.7s 43.9s 4.7s 5.5s 88.7s 10.0s 0.1s 19.8s 1.1s 0.6s 46.6s 4.1s 1.6s 152.1s 13.3s 3.9s 366.2s 41.3s 7.2s 751.0s 79.3s 0.2s 22.7s 2.5s 1.3s 145.7s 22.3s 3.4s 499.2s 90.0s 7.2s 1174.0s 202.7s 12.6s 2531.9s 429.7s 0.5s 59.6s 8.4s 2.5s 410.6s 78.1s 7.5s 1326.4s 313.8s 12.9s 3031.6s 650.2s

Total 4.6s 29.4s 29.4s 51.3s 104.2s 20.9s 51.3s 167.1s 411.4s 837.5s 25.4s 169.3s 592.5s 1383.9s 2974.1s 68.6s 491.3s 1647.7s 3694.7s

Mem 0.7 GB 1.0 GB 2.1 GB 4.2 GB 7.7 GB 1.0 GB 3.8 GB 12.1 GB 28.4 GB 55.8 GB 1.8 GB 11.5 GB 39.1 GB 93.9 GB 186.0 GB 3.4 GB 26.2 GB 91.5 GB 227.4 GB

papers/2022/adrm/CED-WP-2022-002.html. The 2020 Census Disclosure Avoidance System TopDown Algorithm. Harvard Data Science Review (forthcoming. Preprint https://www.census.gov/library/working-papers/2022/adrm/CED-WP2022-002.html). [2] Sergul Aydore, William Brown, Michael Kearns, Krishnaram Kenthapadi, Luca Melis, Aaron Roth, and Ankit A Siva. 2021. Differentially private query release through adaptive projection. In International Conference on Machine Learning.

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

PMLR, 457–467. [3] Barry Becker and Ronny Kohavi. 1996. Adult. UCI Machine Learning Repository. DOI: https://doi.org/10.24432/C5XW20. [4] Mark Bun and Thomas Steinke. 2016. Concentrated Differential Privacy: Simplifications, Extensions, and Lower Bounds. In Proceedings, Part I, of the 14th International Conference on Theory of Cryptography - Volume 9985. [5] Irit Dinur and Kobbi Nissim. 2003. Revealing information while preserving privacy. In PODS. [6] Jinshuo Dong, Aaron Roth, and Weijie J. Su. 2022. Gaussian differential privacy. Journal of the Royal Statistical Society: Series B (Statistical Methodology) 84, 1 (2022), 3–37. https://doi.org/10.1111/rssb.12454 arXiv:https://rss.onlinelibrary.wiley.com/doi/pdf/10.1111/rssb.12454 [7] Cynthia Dwork, Krishnaram Kenthapadi, Frank McSherry, Ilya Mironov, and Moni Naor. 2006. Our Data, Ourselves: Privacy via Distributed Noise Generation. In EUROCRYPT. 486–503. [8] Cynthia Dwork, Frank McSherry, Kobbi Nissim, and Adam Smith. 2006. Calibrating Noise to Sensitivity in Private Data Analysis.. In TCC. [9] Alexander Edmonds, Aleksandar Nikolov, and Jonathan Ullman. 2020. The power of factorization mechanisms in local and central differential privacy. In Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing. 425–438. [10] Miguel Fuentes, Brett Mullins, Yingtai Xiao, Daniel Kifer, Cameron N Musco, and Daniel Sheldon. 2026. Fast Private Adaptive Query Answering for Large Data Domains. In The 29th International Conference on Artificial Intelligence and Statistics. https://openreview.net/forum?id=wkYTC6w2pK [11] Marco Gaboardi, Emilio Jesús Gallego Arias, Justin Hsu, Aaron Roth, and Zhiwei Steven Wu. 2014. Dual query: Practical private query release for high dimensional data. In International Conference on Machine Learning. PMLR, 1170– 1178. [12] Moritz Hardt, Katrina Ligett, and Frank McSherry. 2012. A simple and practical algorithm for differentially private data release. Advances in neural information processing systems 25 (2012). [13] Michael Hay, Vibhor Rastogi, Gerome Miklau, and Dan Suciu. 2010. Boosting the Accuracy of Differentially Private Histograms Through Consistency. Proceedings of the VLDB Endowment 3, 1 (2010). [14] Kaggle. 2021. Loan Prediction Problem Dataset. https://www.kaggle.com/ altruistdelhite04/loan-prediction-problem-dataset. Accessed: May 8th, 2023. [15] Christian Janos Lebeda, Aleksandar Nikolov, and Haohua Tang. 2025. Weighted Fourier Factorizations: Optimal Gaussian Noise for Differentially Private Marginal and Product Queries. arXiv:2512.21499 [cs.DS] https://arxiv.org/abs/2512. 21499 [16] Chao Li, Michael Hay, Gerome Miklau, and Yue Wang. 2014. A data- and workload-aware algorithm for range queries under differential privacy. Proc. VLDB Endow. (2014). [17] Chao Li, Gerome Miklau, Michael Hay, Andrew Mcgregor, and Vibhor Rastogi. 2015. The Matrix Mechanism: Optimizing Linear Counting Queries under Differential Privacy. The VLDB Journal 24, 6 (Dec. 2015), 757–781. https://doi.org/10.1007/s00778-015-0398-x [18] Terrance Liu, Giuseppe Vietri, Thomas Steinke, Jonathan Ullman, and Steven Wu. 2021. Leveraging public data for practical private query release. In International Conference on Machine Learning. PMLR, 6968–6977. [19] Ryan McKenna. 2025. Scaling up the Banded Matrix Factorization Mechanism for Large Scale Differentially Private ML. In The Thirteenth International Conference on Learning Representations. https://openreview.net/forum?id=69Fp4dcmJN [20] Ryan McKenna, Gerome Miklau, Michael Hay, and Ashwin Machanavajjhala. 2018. Optimizing error of high-dimensional statistical queries under differential privacy. Proceedings of the VLDB Endowment 11, 10 (2018). [21] Ryan McKenna, Gerome Miklau, and Daniel Sheldon. 2021. Winning the NIST Contest: A scalable and general approach to differentially private synthetic data. Journal of Privacy and Confidentiality 11, 3 (2021). [22] Ryan McKenna, Brett Mullins, Daniel Sheldon, and Gerome Miklau. 2022. AIM: an adaptive and iterative mechanism for differentially private synthetic data. In VLDB. [23] Ilya Mironov. 2017. Rényi Differential Privacy. In 30th IEEE Computer Security Foundations Symposium, CSF 2017, Santa Barbara, CA, USA, August 21-25, 2017. 263–275. [24] Brett Mullins, Miguel Fuentes, Yingtai Xiao, Daniel Kifer, Cameron Musco, and Daniel R Sheldon. 2024. Efficient and private marginal reconstruction with local non-negativity. Advances in Neural Information Processing Systems 37 (2024), 141356–141389. [25] Aleksandar Nikolov. 2014. New computational aspects of discrepancy theory. Ph.D. Dissertation. Rutgers University-Graduate School-New Brunswick. [26] Wahbeh Qardaji, Weining Yang, and Ninghui Li. 2013. Understanding Hierarchical Methods for Differentially Private Histograms. Proc. VLDB Endow. 6, 14 (2013). [27] U. S. Census Bureau. [n.d.]. American Community Survey Data. https://www. census.gov/programs-surveys/acs/data/summary-file.html.

[28] U. S. Census Bureau. 2020. Decennial Census Data Products. https://www.census. gov/programs-surveys/decennial-census/data-products.2020.html. [29] U.S. Bureau of Labor Statistics. (retrieved August 21, 2024). QCEW Public-Use Data Files. https://www.bls.gov/cew/downloadable-data-files.htm. [30] Xiaokui Xiao, Guozhang Wang, and Johannes Gehrke. 2011. Differential Privacy via Wavelet Transforms. IEEE Transactions on Knowledge and Data Engineering 23, 8 (2011), 1200–1214. [31] Yingtai Xiao, Zeyu Ding, Yuxin Wang, Danfeng Zhang, and Daniel Kifer. 2021. Optimizing fitness-for-use of differentially private linear queries. In VLDB. [32] Yingtai Xiao, Guanlin He, Levent Toksoz, Zeyu Ding, Danfeng Zhang, and Daniel Kifer. 2025. ResidualPlanner+: a scalable matrix mechanism for marginals and beyond. arXiv:2305.08175v3 [cs.DB] https://arxiv.org/abs/2305.08175v3 [33] Yingtai Xiao, Guanlin He, Danfeng Zhang, and Daniel Kifer. 2023. An optimal and scalable matrix mechanism for noisy marginals under convex loss functions. In Proceedings of the 37th International Conference on Neural Information Processing Systems (New Orleans, LA, USA) (NeurIPS). Curran Associates Inc., Red Hook, NY, USA, Article 901, 45 pages. [34] Yingtai Xiao, Guanhong Wang, Danfeng Zhang, and Daniel Kifer. 2023. Answering Private Linear Queries Adaptively Using the Common Mechanism. Proc. VLDB Endow. 16, 8 (April 2023), 1883–1896. https://doi.org/10.14778/3594512. 3594519 [35] Ganzhao Yuan, Yin Yang, Zhenjie Zhang, and Zhifeng Hao. 2016. Convex Optimization for Linear Query Processing under Approximate Differential Privacy. In Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining. [36] Ganzhao Yuan, Zhenjie Zhang, Marianne Winslett, Xiaokui Xiao, Yin Yang, and Zhifeng Hao. 2012. Low-Rank Mechanism: Optimizing Batch Queries under Differential Privacy. Proc. VLDB Endow. 5, 11 (July 2012), 1352–1363. https: //doi.org/10.14778/2350229.2350252 [37] Ganzhao Yuan, Zhenjie Zhang, Marianne Winslett, Xiaokui Xiao, Yin Yang, and Zhifeng Hao. 2015. Optimizing batch linear queries under exact and approximate differential privacy. ACM Transactions on Database Systems (TODS) 40, 2 (2015), 1–47. [38] Dan Zhang, Ryan McKenna, Ios Kotsogiannis, Michael Hay, Ashwin Machanavajjhala, and Gerome Miklau. 2018. Ektelo: A framework for defining differentiallyprivate computations. In Proceedings of the 2018 International Conference on Management of Data. 115–130. [39] Jun Zhang, Graham Cormode, Cecilia M. Procopiuc, Divesh Srivastava, and Xiaokui Xiao. 2017. PrivBayes: Private Data Release via Bayesian Networks. ACM Trans. Database Syst. 42, 4, Article 25 (oct 2017), 41 pages.

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

PROOFS FROM SECTION 6

Theorem 2. Let A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ } be a subset of attributes. Let q® be a query that is answerable from the marginal x® A . For any A ′ ⊆ A, define Q A ′ using Equation 3 so that q® ⇒A ′ = q® Q A ′ . Then q® · x® A = Í q® ⇒A ′ · x® A ′ . Furthermore, if A ′ ≠ A ′′ then Q⊤A ′ Q A ′′ = 0. A′ ⊆ A

Proof of Theorem 2. We prove orthogonality first. Represent the matrices in their kronecker product forms: Q A ′ = W1′ ⊗W2′ ⊗ · · · ⊗Wℓ′ and Q A ′′ = W1′′ ⊗W2′′ ⊗ · · · ⊗Wℓ′′ . Then since A ′ ≠ A ′′ , there exists an index 𝑗 such that W′𝑗 ≠ W′′𝑗 . Specifically, one of those matrices must be equal to 𝑑1𝑗 1𝑑 𝑖 𝑗 and the other must be equal to I𝑑 𝑖 𝑗 − 𝑑1𝑗 1𝑑 𝑖 𝑗 1𝑇𝑑𝑖 . In this case, clearly (W′𝑗 )𝑇 W′′𝑗 = 0 𝑗

and (W′′𝑗 )𝑇 W′𝑗 = 0. Thus, Q𝑇A ′ Q A ′′ = ((W1′ )𝑇 W1′′ )⊗ · · · ⊗((W′𝑗 )𝑇 W′′𝑗 )⊗ · · · ⊗((Wℓ′ )𝑇 Wℓ′′ ) = 0 Next, we note that for any A ′ ⊆ A, then x® A ′ can be obtained from x® A by left multiplication by a matrix. Specifically, x® A ′ = V1′ ⊗V2′ ⊗ · · · ⊗Vℓ′ x® A   if 𝐴𝑖 ∉ A ′  1𝑇𝑑𝑖  𝑗 where V′𝑗 =   I𝑑 𝑖 𝑗 otherwise  Let Q A ′ = W1′ ⊗ · · · ⊗Wℓ′ be the matrix that turns q® into q® ⇒A ′ by right multiplication as in Equation 3 and let V A ′ = V1′ ⊗ · · · ⊗Vℓ′ be the matrix that converts x® A into x® A ′ by left multiplication. It is clear that W′𝑗 = 𝑑1𝑖 1𝑑 𝑖 𝑗 if and only if V′𝑗 = 1𝑇𝑑𝑖 . Therefore 𝑗

𝑗

 1 𝑇   𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗  ′ ′ W𝑗 V𝑗 = 1 𝑇   I𝑑 𝑖 𝑗 − 𝑑 𝑖 1𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗

otherwise 𝑗

q® ⇒A ′ · x® A ′

A′ ⊆ A

! =

∑︁

q® Q

A′

x® A = q®

A′ ⊆ A

= q®

A′

x® A

A′ ⊆ A

∑︁ Ìℓ 𝑗=1

A′ ⊆ A

= q®

∑︁

Ìℓ

1 𝑇 if 𝐴𝑖 𝑗 ∉A ′ 𝑑 𝑖 1𝑑 𝑖 𝑗 1𝑑 𝑖 𝑗 𝑗 I𝑑 𝑖 − 𝑑 1 1𝑑 𝑖 1𝑇 otherwise 𝑖𝑗 𝑗 𝑑𝑖 𝑗 𝑗

I𝑑 𝑖 𝑗 − 𝑑1𝑖 1𝑑 𝑖 𝑗 1𝑇𝑑𝑖 𝑗 𝑗 𝑇 1 + 𝑑 𝑖 1𝑑 𝑖 𝑗 1𝑑 𝑖 x® A

!! x® A

𝑗=1

𝑗

• The weighted sum-squared error of Compete over WKLoad is equal to the sum of the weighted sum-squared errors of Compete on the subworkloads created by QuerySmasher. • For any given A, let M 1, . . . , M𝑚 be the subset of Gaussian linear mechanisms of Compete that are answerable from the marginal x® A . Then M 1, . . . , M𝑚 together correspond to a feasible (but potentially suboptimal) solution to the optimization problem in Equation 4 for subworkload WKLoad⇒A (i.e., the W in that equation is formed from the query vectors in WKLoad⇒A ). • When QuerySmasher uses solvers that solve Equation 4 to optimality, the weighted sum squared error of QuerySmasher over WKLoad is ≤ the weighted sum-squared error of Compete over WKLoad. Proof of Theorem 7. First, we note that if M 1, . . . , M𝑚 are answerable from marginal x® A then they can be turned into a single matrix mechanism M A that is answerable from x® A . Indeed, let B1, . . . , B𝑚 be their respective strategy matrices and Σ1, . . . , Σ𝑚 be their covariance matrices. Form the matrix B A by stacking the B𝑖 : B1 ©B ª 2® B A = . ®® .. ® «B𝑚 ¬ and form Σ A as the block diagonal matrix: Σ © 1 0 Σ A = ... «0

if 𝐴𝑖 𝑗 ∉ A ′



Therefore ∑︁

PROOFS FROM SECTION 7

Theorem 7. Given a workload WKLoad of linear queries that can be answered from data marginals, and a privacy cost bound 𝛽, let Compete denote a matrix mechanism that can answer the queries in WKLoad, is splittable, and which uses gaussian linear mechanisms that all respect residual spaces. Then:

𝑗

0 Σ2 .. . 0

··· .. . .. . ···

0 ª ® 0 ®® .. ®® . ® Σ𝑚 ¬

and define M A (®x A ) ≡ B A x® A + 𝑁 (0, Σ A ). Then it is easy to see that the distribution of the output of M A is the same as the joint distribution of the outputs of M 1, . . . , M𝑚 . Our next observation is that we can convert M A into a mechanism that has a diagonal covariance matrix. Define B∗A = Σ−1/2 A BA , 1/2 1/2 1/2 where Σ1/2 A is the matrix square root (Σ A Σ A = Σ A ). Σ A exists because Σ is symmetric positive definite as it is a covariance matrix for Gaussian noise. Define the gaussian linear mechanism M ∗A (®x A ) ≡ B∗A x® A + 𝑁 (0, I). Then it is easy to check that: ∗ • The output distributions of M A (®x A ) and Σ1/2 x A ) (i.e., A M A (® ∗ 1/2 run M A and multiply the result by Σ ) are the same. • The output distributions of Σ−1/2 x) and M ∗A (®x) are the A M A (® same.

Thus M A and M ∗A can be obtained from each other by postprocessing, and so by the postprocessing invariance of differential privacy, both mechanisms have the same privacy cost. Furthermore, in whatever situation one needs to use the output of M (e.g., in

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

the reconstruction function), one can instead use the output of M ∗ and left multiply it by Σ−1/2 A . Next, every M ∗A , whose input should be x® A can be converted 𝑓 𝑢𝑙𝑙 into a mechanism M A that operates on the full data x® . Specifically, define: (𝑚𝑎𝑟𝑔)

= V1 ⊗V2 ⊗ · · · ⊗V𝑘 ( I𝑑 if 𝐴𝑖 ∈ A where V𝑖 = 𝑇 𝑖 1𝑑 𝑖 otherwise 𝑄A

(𝑚𝑎𝑟𝑔)

M A (®x) = B∗A 𝑄 A 𝑓 𝑢𝑙𝑙

𝑓 𝑢𝑙𝑙

Since the M A have independent noise and independent strat𝑓 𝑢𝑙𝑙

egy matrices, then given their noisy outputs 𝑧®A = M A (®x) for A ∈ Good, the optimal linear reconstruction of the answer to a q® ∈ WKLoad that is answerable from some x® A ∗ is the query times the sum of the strategy pseudoinverse times the noisy answers: ! ∑︁ (𝑚𝑎𝑟𝑔) † 𝑓 𝑢𝑙𝑙 ∗ q® (B A 𝑄 A ) 𝑧®A A ∈Good

x® + 𝑁 (0, I)

Step 1: optimal reconstruction. We next consider what is the optimal linear reconstruction function. First, note that the row spaces of B A and B∗A are the same, and since residual subsapces are respected, those row spaces are are contained in the row space of 𝑃𝑟𝑜 𝑗 A ≡ I𝑑 𝑖1 − 𝑑1𝑖 1𝑑 𝑖1 1𝑑⊤𝑖 ⊗ · · · 1 1 ⊗ I𝑑 𝑖ℓ − 𝑑1𝑖 1𝑑 𝑖ℓ 1𝑑⊤𝑖 . ℓ

(𝑚𝑎𝑟𝑔)

Q𝑟𝑒𝑠𝑖𝑑 ≡ 𝑃𝑟𝑜 𝑗 A 𝑄 A = V1 ⊗ · · · ⊗V𝑘 A ( 1 ⊤ I𝑑 − 1𝑑 1 if 𝐴𝑖 ∈ A where V𝑖 = 𝑇 𝑖 𝑑 𝑖 𝑖 𝑑 𝑖 1𝑑 𝑖 otherwise ′ 𝑟𝑒𝑠𝑖𝑑 ⊤ and this means Q𝑟𝑒𝑠𝑖𝑑 A (Q A ′ ) = 0 when A ≠ A. ′ Combining these facts, for any A ≠ A, (𝑚𝑎𝑟𝑔) ⊤

)(B∗A ′ 𝑄 A ′

∑︁

(𝑚𝑎𝑟𝑔)

q®𝑄 A ∗

A ∈Good (𝑚𝑎𝑟𝑔) † ∗ † q®𝑄 A ∗ (Q𝑟𝑒𝑠𝑖𝑑 A ) (B A ) 𝑧®A

∑︁

† (B∗A Q𝑟𝑒𝑠𝑖𝑑 A ) 𝑧®A

A ∈Good

Next, we note that † (Q𝑟𝑒𝑠𝑖𝑑 A ) = V1 ⊗ · · · ⊗V𝑘 ( I𝑑 − 1 1𝑑 1⊤ where V𝑖 = 1 𝑖 𝑑 𝑖 𝑖 𝑑 𝑖 𝑑 𝑖 1𝑑 𝑖

if 𝐴𝑖 ∈ A otherwise

ℓ

where A = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ }. Since 𝑃𝑟𝑜 𝑗 A is an idempotent matrix (𝑃𝑟𝑜 𝑗 A 𝑃𝑟𝑜 𝑗 A = 𝑃𝑟𝑜 𝑗 A ), it means that B∗A 𝑃𝑟𝑜 𝑗 A = B∗A . Define

(𝑚𝑎𝑟𝑔)

(𝑚𝑎𝑟𝑔) † (B∗A 𝑄 A ) 𝑧®A

A ∈Good

and similarly any query q® over x® can be converted into a query (𝑚𝑎𝑟𝑔) q® 𝑓 𝑢𝑙𝑙 over D as follows: q® 𝑓 𝑢𝑙𝑙 = q®𝑄 A . (note q® · x® A = q® 𝑓 𝑢𝑙𝑙 · x® ).

(B∗A 𝑄 A

∑︁

(𝑚𝑎𝑟𝑔) q®𝑄 A ∗

)

(𝑚𝑎𝑟𝑔) (𝑚𝑎𝑟𝑔) ⊤ = (B∗A 𝑃𝑟𝑜 𝑗 A 𝑄 A )(B∗A ′ 𝑃𝑟𝑜 𝑗 A ′ 𝑄 A ′ ) ∗ 𝑟𝑒𝑠𝑖𝑑 ∗ 𝑟𝑒𝑠𝑖𝑑 ⊤ = (B A Q A )(B A ′ Q A ′ )

(𝑚𝑎𝑟𝑔)

† (Q𝑟𝑒𝑠𝑖𝑑 A ) = W1 ⊗ · · · ⊗W𝑘   1 if 𝐴𝑖 ∉ A ∗ ∪ A    0   if 𝐴𝑖 ∈ A \ A ∗ where W𝑖 = 1  if 𝐴𝑖 ∈ A ∗ \ A  𝑑 𝑖 1𝑑 𝑖    1 ⊤  I𝑑 𝑖 − 𝑑 1𝑑 𝑖 1𝑑 if 𝐴𝑖 ∈ A ∗ ∩ A 𝑖 𝑖 

𝑄 A∗

(𝑚𝑎𝑟𝑔)

∗ † and so 𝑄 A ∗ (Q𝑟𝑒𝑠𝑖𝑑 A ) equals the 0 matrix if A ⊈ A and otherwise is the matrix from Equation 3 used to project q® into q® ⇒A . Hence, the reconstruction of the answer to a q® ∈ WKLoad that is answerable from some x® A ∗ is ! ∑︁ (𝑚𝑎𝑟𝑔) † 𝑓 𝑢𝑙𝑙 ∗ q® (B A 𝑄 A ) 𝑧®A A ∈Good

=0 𝑓 𝑢𝑙𝑙 𝑓 𝑢𝑙𝑙 Thus M A and M A ′ use orthogonal (hence independent) strat-

egy matrices and hence access orthogonal information from x® . Next, if QuerySmasher generates a non-empty subworkload 𝑓 𝑢𝑙𝑙 WKLoad⇒A then Compete must have a corresponding M A (since a non-empty subworkload WKLoad⇒A implies q® 𝑓 𝑢𝑙𝑙 for some q® ∈ WKLoad must have a component in a direction orthogo𝑓 𝑢𝑙𝑙 nal to all of the strategy matrices of M A ′ for A ′ ≠ A, since the mechanisms of Compete respect residual spaces). Furthermore, if for some A there exists an empty subworkload WKLoad⇒A , 𝑓 𝑢𝑙𝑙 then there shouldn’t be any M A as its strategy matrix would be orthogonal to all queries (hence running such a mechanism would consume privacy budget for no improvement in accuracy). Hence, define Good to be the set of all A for which WKLoad⇒A is nonempty. The previous discussion means that Compete has a 𝑓 𝑢𝑙𝑙 mechanism M A if and only if A ∈ Good.

∑︁

(𝑚𝑎𝑟𝑔)

q®𝑄 A ∗

† ∗ † (Q𝑟𝑒𝑠𝑖𝑑 A ) (B A ) 𝑧®A

A ∈Good

∑︁

q® ⇒A (B∗A ) †𝑧®A

A ⊆ A∗

we also note that q® ⇒A (B∗A ) †𝑧®A is the optimal reconstruction (based on the pseudoinverse) of q® ⇒A based on the output 𝑧®A of M ∗A . Thus for Compete, just like for QuerySmasher, the reconstructed answer to q® is the sum of the reconstructed answers to its projected subqueries. Furthermore, since each projected subquery is answered using noise that is independent of the other subqueries (because it comes from a separate gaussian linear mechanism), the variance of q® is the sum of the variances of the subqueries. Since q® and its subqueries also have the same weights, then the weighted sum squared error of reconstructed answers to queries in WKLoad is equal to the sum of the weighted sum-squared errors of reconstructed answers to the subworkloads. This proves the first claim of the theorem.

Guanlin He, Yingtai Xiao, Jiamu Bai, Xin Gu, Zeyu Ding, Wenpeng Yin, and Daniel Kifer

Step 2: separate subworkload optimization. The previous results showed that the role of M ∗A is to serve as a gaussian linear mechanism from whose output the answers to queries in WKLoad⇒A are reconstructed. If 𝛽 A is the privacy cost of M ∗A √︁ then rescaling the entries of B∗A by 1/ 𝛽 A results in a gaussian linear mechanism whose privacy cost is 1. Hence (B∗A )𝑇 (B∗A )/𝛽 A is a feasible solution to the matrix V in Equation 4. This proves the second claim. Hence, if we obtain the optimal matrix mechanism for the subworkload WKLoad⇒A and rescale its noise variance by 1/𝛽 A (so that it has privacy cost 𝛽 A ), and then use that in place of M ∗A , the weighted sum-squared error over the subworkload cannot increase. If we then replace the 1/𝛽 A rescaling by the optimal rescaling used by QuerySmasher, the overall weighted sum-squared error cannot increase. This proves the third claim. □

INTEGRATION OF FOURIER BASIS INTO QuerySmasher

In this section, we explain how to integrate the Fourier basis [15] into QuerySmasher so that it can be used as a fast, but approximate, method to optimize a subworkload WKLoad⇒A ′ . If the Fourier basis is used for all subworkloads, the result becomes equivalent to WFF. Following the ideas of Lebeda et al. [15], the high level idea is to take the multidimensional Fourier transform of the tensor view of a marginal T A ′ , add different amounts of noise to each Fourier coefficient, and then take the inverse (multidimensional) Fourier transform of the noisy coefficients. This results in a privacyT A ′ that can be used to directly answer queries preserving marginal b in the subworkload. Thus, the goal is to determine how much noise different coefficients need. Expressing query variances in terms of tunable parameters. A marginal T A ′ , where A ′ = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ }, can be indexed using a tuple of ( 𝑗1, . . . , 𝑗ℓ ) of ℓ indexes (one for each dimension), where 𝑗1 range from 0 to 𝑑 𝑖 1 − 1 (inclusive), the index 𝑗2 ranges from 0 to 𝑑 𝑖 2 − 1, etc. The parameter set 𝜃 can be thought of as a similarly sized table, but only some entries correspond to actual parameters. Specifically, 𝜃 [ 𝑗1, . . . , 𝑗ℓ ] is a tunable parameter if 𝑗1 ∈ [1, . . . , 𝑑 𝑖 1 −1], . . . , 𝑗ℓ ∈ [1, . . . , 𝑑 𝑖 ℓ −1] and the tuple ( 𝑗1, . . . , 𝑗ℓ ) is ≤ than the tuple (𝑑 𝑖 1 − 𝑗1, 𝑑 𝑖 2 − 𝑗2, . . . , 𝑑 𝑖 ℓ − 𝑗ℓ ) in dictionary order. Given a query q ∈ WKLoad⇒A ′ expressed in tensor form, its data-independent variance in terms of the tunable parameters can be computed as follows: (1) Let 𝐹e be the inverse multidimensional discrete Fourier transform of q. (2) The variance is a linear combination of the parameters in 𝜃 , where the coefficient of the parameter 𝜃 [ 𝑗1, . . . , 𝑗ℓ ] is || 𝐹e[ 𝑗1, . . . , 𝑗ℓ ]|| 2 if ( 𝑗1, . . . , 𝑗ℓ ) = (𝑑 𝑖 1 −𝑗1, 𝑑 𝑖 2 −𝑗2, . . . , 𝑑 𝑖 ℓ −𝑗ℓ ) and 4|| 𝐹e[ 𝑗1, . . . , 𝑗ℓ ]|| 2 otherwise. Furthermore, if q′ is any circular-shift of q, then they have the same variance. Setting the tunable parameters. The weighted sum of variances of queries in WKLoad⇒A ′ can be computed by multiplying each

query’s variance by the query’s weight and adding them up. Thus, the weighted sum of variances is a linear combination of the parameters in 𝜃 . Let 𝑐 𝑗1 ,...,𝑗ℓ be the coefficient of parameter 𝜃 [ 𝑗1, . . . , 𝑗ℓ ]. Let 𝛾 be the sum of the square roots of the linear coefficients. Then √ the value of parameter is 𝜃 [ 𝑗1, . . . , 𝑗ℓ ] = 𝛾/ 𝑐 𝑗1 ,...,𝑗ℓ . How to run the mechanism. Setting the values for the parameter set 𝜃 defines a mechanism. This mechanism can be run as follows: (1) Create a tensor 𝑍 with entries initialized to 0. (2) Let 𝐹 be the multidimensional discrete Fourier transform of T A ′ . (3) For each valid parameter index ( 𝑗1, . . . , 𝑗ℓ ), let 𝜎 2 = 𝜃 [ 𝑗1, . . . , 𝑗ℓ ] and let 𝑎 + 𝑏 i = 𝐹 [ 𝑗1, . . . , 𝑗ℓ ]. (a) Add noise: 𝑎e = 𝑎 + 𝑁 (0, 𝜎 2 ) and 𝑏e = 𝑏 + 𝑁 (0, 𝜎 2 ). (b) Store: set 𝑍 [ 𝑗1, . . . , 𝑗ℓ ] = 𝑎e+𝑏ei and 𝑍 [𝑑 𝑖 1 − 𝑗1, . . . , 𝑑 𝑖 ℓ − 𝑗ℓ ] = 𝑎e − 𝑏ei. (4) Set e T A ′ to be the inverse multidimensional discrete Fourier transform of 𝑍 and return e TA′ . How to reconstruct a query answer. Any query q ∈ WKLoad⇒A ′ can be directly answered from e T A ′ in the usual way: perform elementwise multiplication of the tensor view of q with e T A ′ and then add up all the entries.

VARIATIONS OF THE MATRIX MECHANISM OPTIMIZATION PROBLEM

We present additional variations of the matrix mechanism construction optimization problem from Section 6.3, when the goal is to minimize weighted sum squared error. Let W be the matrix whose rows are the query vectors in WKLoad⇒A ′ . Let D be a diagonal matrix, where entry D[𝑖, 𝑖] is the weight of the query that is in the 𝑖 th row of W Suppose we have an independent linear basis U whose row span is equal to the row span of W and U† is the pseudoinverse of U. Then one can solve the following optimization problem: arg min trace((U† ) ⊤ W⊤ DWU† V −1 ) V

s.t. V is symmetric and invertible and (U⊤ VU) [𝑖, 𝑖] ≤ 1 for all 𝑖 This is a semidefinite program that can be solved using off-theshelf optimizers. Once a solution is obtained, B can be computed by taking the Cholesky decomposition V = L⊤ L and setting B = LU. The rest is the same as in Section 6.3. Some solvers for the optimization problem in Section 6.3 require W to be full-rank [37]. However, the W we get from a workload partition WKLoad⇒A ′ is typically not full rank. For example, in Figure 2 the query vectors all have entries that add up to 0. In the tensor view of each query, the sum along any dimension is clearly 0 (this is a direct consequence of the query decomposition construction). In fact, the rows of the query matrix W for the workload WKLoad⇒A ′ , where A ′ = {𝐴𝑖 1 , . . . , 𝐴𝑖 ℓ }, always belong to the subspace spanned by the rows of 1 1 I𝑑 𝑖 1 − 1𝑑 𝑖1 1𝑇𝑑𝑖 ⊗ · · · ⊗ I𝑑 𝑖ℓ − (6) 1𝑑 𝑖ℓ 1𝑇𝑑𝑖 ℓ 1 𝑑 𝑖1 𝑑 𝑖ℓ

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

Table 8: WRMSE at privacy cost 1 of weighted 1-way and 2-way marginal workload as 𝑑 (dataset dimensions) and 𝑛 (domain size of each attribute) vary. on [𝑛] 40 , varying domain size 𝑛. The RP+ and HDMM implementations do not support these queries. 𝑛 10

WFF QS 39.22 39.22 66.58 56.59 69.25 68.82 66.70 66.70 75.74 47.33 107.17 51.94 180.03 174.49 42.90 42.90 105.43 84.18 106.37 105.93 105.53 105.53 117.65 61.68 171.30 66.89 391.96 386.32 44.16 44.16 132.39 101.95 132.15 131.76 132.48 132.48 146.56 102.25 215.96 106.11 not finished 44.80 44.80 153.42 115.07 152.33 152.00 153.50 153.50 169.06 139.22 250.83 183.22 818.49 817.66

To sidestep this issue of a non-full-rank matrix, QuerySmasher can be modified in the following way: (1) If QuerySmasher is being used with a matrix mechanism constructor that needs a full rank query workload, one can change how WKLoad⇒A ′ is constructed. Given a query that is answerable from a marginal on A (with A ′ ⊆ A), the decomposition in Equation 3 would use the matrices I𝑑 𝑖 𝑗 instead of I𝑑 𝑖 𝑗 − 𝑑1𝑖 1𝑑 𝑖 𝑗 1𝑇𝑑𝑖 . 𝑗

𝑗

(2) After the convex matrix mechanism constructor has returned a Btmp , we project its rows onto the subspace spanned by the rows of Equation 6 by setting B′ to be Btmp times the matrix in Equation 6. (3) The rest is the same as in Section 6.3. (4) The reconstruction phase of Algorithm 1 uses the original query decomposition.

EVALUATION ON WEIGHTED QUERIES

Workload weights. We also experiment with weighted workloads to see if the results qualitatively change. For weighted workloads, we randomly sample an integer weight for each query from 1 to 5, to simulate the situation where some queries are more important than others. The weighted RMSE is defined as follows: v u t 𝑊 𝑅𝑀𝑆𝐸 =

weight(®q)Var(®q)

q® ∈WKLoad

|WKLoad| Table 8 shows the results. Qualitatively, little changes from Table 4, indicating that the latest scalable matrix mechanisms may be robust to weighted queries.

Accurate and Scalable Matrix Mechanisms via Divide and Conquer

Related documents