Differential Privacy Over Riemannian Manifolds
Year: 2022

Abstract:
In this work we consider the problem of releasing a differentially private statistical summary that resides on a Riemannian manifold. We present an extension of the Laplace or K-norm mechanism that utilizes intrinsic distances and volumes on the manifold. We also consider in detail the specific case where the summary is the Fréchet mean of data residing on a manifold. We demonstrate that our mechanism is rate optimal and depends only on the dimension of the manifold, not on the dimension of any ambient space, while also showing how ignoring the manifold structure can decrease the utility of the sanitized summary. We illustrate our framework in two examples of particular interest in statistics: the space of symmetric positive definite matrices, which is used for covariance matrices, and the sphere, which can be used as a space for modeling discrete distributions.

heading: 1 Introduction
text: Over the last decade we have seen a tremendous push for the development and application of methods in data privacy. This surge has been fueled by the production of large sophisticated datasets alongside increasingly complex data gathering technologies. One theme that has emerged with the proliferation of highly structured and dynamic data is the importance of exploiting underlying structures in the data or models to maximize utility while controlling disclosure risks. In this paper we consider the problem of achieving pure Differential Privacy, DP, when the statistical summary to be released takes values on a complete Riemannian manifold.
Riemannian manifolds are used extensively in the analysis of data or parameters that are inherently nonlinear, meaning, either addition or scalar multiplication may cause the summary to leave the manifold, or such operations are not even well defined. Classic examples of such objects include spatio-temporal processes, covariance matrices, projections, rotations, compositional data, densities, and shapes. Traditional privacy approaches for handling such objects typically consist of utilizing an ambient or embedding space that is linear so that standard DP tools can be employed. For example, Karwa and Slavković [2016] considered the problem of releasing private degree sequences of a graph which required them to project back onto a particular convex hull as a post-processing step. Such an approach is a natural starting point and reasonable so long as the space doesn’t exhibit too much curvature. However, there are several interrelated motivations for working with the manifolds directly. First, if one employs an ambient space (known as taking an extrinsic approach), then calculations such as the sensitivity may depend on the dimension of the ambient space, which will in turn impact
35th Conference on Neural Information Processing Systems (NeurIPS 2021).
the utility of the private statistical summary. For example, minimax rates in DP typically scale polynomially in the dimension [e.g. Hardt and Talwar, 2010, Bun et al., 2018, Kamath et al., 2019]. Second, the Whitney embedding theorem states that, in the worst case, to embed a manifold in Euclidean space requires a space that is twice the dimension of the manifold. Third, if the manifold exhibits substantial curvature, then even small distances in the ambient space may result in very large distances on the manifold. Lastly, the choice of the ambient space may be arbitrary and one would ideally prefer if this choice did not play a role in the resulting statistical analysis.
Related Literature: To the best of our knowledge, general manifolds have not been considered before in the DP literature. The closest works come from the literature on private covariance matrix estimation and principal components [Blum et al., 2005, Awan et al., 2019, Amin et al., 2019, Kamath et al., 2019, Biswas et al., 2020, Wang and Xu, 2020]. While not always explicitly described in some works [Wang et al., 2013, Wei et al., 2016], these objects lie in nonlinear manifolds, namely, the space of symmetric positive definite matrices (SPDM) and the space of projections matrices respectively, called the Stiefel manifold. For example, in Chaudhuri et al. [2013] they consider the problem of generating a synthetic PCA projection by using the matrix Bingham distribution, a distribution over the Stiefel manifold [Khatri and Mardia, 1977, Hoff, 2009]. In contrast, producing private covariance matrix estimates usually involves adding noise in a way that preserves symmetry, but does not use any deeper underlying manifold structure. A related problem comes from the literature on private manifold learning [Choromanska et al., 2016, Vepakomma et al., 2021], though this is entirely distinct from the present work, which assumes the underlying manifold is known, usually because of some physical constraints on the data or statistical summaries.
Contributions: In this paper we utilize tools from Differential Geometry that allow us to extend the Laplace mechanism for -Differential Privacy to general Riemannian manifolds. Under this framework, we consider the problem of privately estimating the Fréchet mean of data lying on a d-dimensional manifold. We are able to bound the global sensitivity of the mean and provide bounds on the magnitude of the privacy noise, as measured using the distance on the manifold, that match the optimal rates derived in Euclidean spaces. However, we demonstrate the influence of curvature of the space in understanding the sensitivity of the mean, and how the situation becomes especially challenging on positively curved spaces. We conclude by providing two specific numerical examples that elucidate this phenomenon: the first considers data coming from the space of positive definite matrices equipped with a geometry that results in negative curvature, while the second example considers data lying on the sphere, which has constant positive curvature.

heading: 2 Notation and Background
text: In this section we provide the basic notation, terminology, and mathematical concepts needed from differential geometry. For a more detailed treatment of differential geometry and Shape Analysis there are many excellent classic texts, e.g., Gallot et al. [1990], Lang [2002], Dryden [2014], Srivastava and Klassen [2016], Lee [2018], while an overview of DP can be found in Dwork and Roth [2014].
Throughout the paper we let M denote a d-dimensional complete Riemannian manifold. For m ∈ M, denote the corresponding tangent space as TmM. We assume M is equipped with a Riemannian metric {〈·, ·〉m : m ∈ M}, which is a collection of inner products over the tangent spaces {TmM : m ∈M} that vary smoothly in m. Two quantities that will be used extensively in this work are that of the distance and volume induced from the Riemannian metric. Consider two points m1,m2 ∈M and a smooth path γ : [0, 1]→M such that γ(0) = m1 and γ(1) = m2. The derivative γ̇(t) represents the velocity of γ as it passes through the point γ(t) and can thus be identified as an element of the tangent space Tγ(t)M. We define the length of the curve as
L(γ) := ∫ 1 0 〈γ̇(t), γ̇(t)〉1/2γ(t) dt.
The distance between m1 and m2 is the infimum over all possible paths connecting the two points ρ(m1,m2) := inf
γ:γ(0)=m1 γ(1)=m2
L(γ).
If this distance is achieved by a particular path, γ, then we say that γ is a geodesic. Geodesics generalize the concept of a straight line to nonlinear spaces, and are thus a natural tool to consider
when generalizing DP perturbation mechanisms. The distance ρ(·, ·) defines a valid distance metric overM; we say thatM is complete if it is complete as a metric space. By the Hopf-Rinow theorem [Lee, 2018, Theorem 6.19], this is equivalent to saying that to every pair of points there exists a minimizing geodesic, though if the points are far enough apart it need not be unique.
In the next section we will use the exponential map, expm : TmM → M, which is a means of moving between the manifold and tangent spaces. If there exists a unique geodesic, γ, between between points γ(0) = m1 and γ(1) = m2, then the exponential map is defined as the mapping expm1(γ̇(0)) = m2. In other words, expm maps initial velocities to points onM through the use of geodesic curves. If the manifold is complete, then the exponential map is locally diffeomorphic, meaning that for anym ∈M there exists an open neighborhood Ur(m) ofM that is diffeomorphic to an open ball Br(0) centred at the origin of TmM. This ensures that the inverse exp−1m : Ur → Br(0) known as the inverse-exponential (also known as logarithm) map exists locally. The injectivity radius of m is the supremum of r 7→ Ur over all such radii r. The injectivity radius ofM, denoted by injM, is defined to be the infimum of the injectivity radii of all points m ∈M. The Riemannian metric can be used to define a notion of volume, called the Riemannian volume measure denoted as µ, which acts analogously to Lesbesgue measure in Euclidean space. To define the measure, it helps to employ a chart, (U, φ), although the final definition will not depend on which chart we choose. Since φ : U → φ(U) ⊂ Rd is a homeomorphism, at each m ∈ U the inverse φ−1 induces a basis, ∂x1, . . . , ∂xd on TmM and a corresponding dual basis dx1, . . . , dxd on TmM∗, the dual space of TmM. Then the Riemannian volume form is defined as √ |g|dx1 ∧ · · · ∧ dxd, where gij = 〈∂xi, ∂xj〉m and | · | is the absolute value of the determinant, which can be shown to be invariant to the choice of chart. This induces a volume over the set U and upon employing a partition of unity, one can define µ over the entire manifoldM, equipped with the Borel σ-algebra.

heading: 3 Manifold Perturbations and Differential Privacy
text: Denote the dataset as D = {x1, . . . , xn} with the data coming from points xi collected from an arbitrary set X . We aim to release a statistical summary f(D) which takes values onM. Defining differential privacy over a Riemannian manifold presents no major challenge since it is a well defined concept over any measurable space [Wasserman and Zhou, 2010, Awan et al., 2019], which includes Riemannian manifolds equipped with the Borel σ-algebra. Denote the (random) sanitized version of f(D) as f̃(D). We can then define what it means for f̃(D) to satisfy -DP.
Definition 1. A family of randomized summaries, {f̃(D) ∈ M : D ∈ Xn}, is said to be - differentially private with > 0, if for any adjacent database D′, denoted as D ∼ D′, differing in only one record we have P (f̃(D) ∈ A) ≤ e P (f̃(D′) ∈ A), for any measurable set A.
In a similar fashion, we can extend the notion of sensitivity to Riemannian manifolds using the distance function. However, this definition is by no means canonical and intimately connected to the type of noise one intends to use [Mirshani et al., 2019]. Definition 2. A summary f is said to have a global sensitivity of ∆ <∞, with respect to ρ(·, ·), if for any two adjacent databases D and D′ we have
ρ(f(D), f(D′)) ≤ ∆.
With a bounded sensitivity, it still isn’t obvious how to produce a differentially private summary. In particular, we no longer have linear perturbations or classic noise mechanisms such as Laplace or Gaussian. However, the Riemannian structure allows us to use the volume measure as a base measure, similar to the Lebesgue measure on Euclidean spaces. We can then define a new distribution that can be viewed as a generalization of the Laplace distribution over manifolds [e.g. Hajri et al., 2016, for SPDM]. It is worth noting that this distribution is not equivalent to the multivariate Laplace in Euclidean settings, instead it is an instantiation of the K-norm distribution [Hardt and Talwar, 2010]. Definition 3. A probability measure P overM is called a Laplace distribution with footpoint η ∈M and rate σ > 0 if for any measurable set A we have
P (A) = ∫ A C−1η,σe −ρ(η,m)/σ dµ(m),
where 0 < Cη,σ <∞ is the normalizing constant and µ is the Riemannian volume measure overM.
Unlike with Euclidean distance, the normalizing constant might only be finite for certain values of σ. The foot-point represents the center of the distribution but, in general, need not be equal to the mean of the distribution (which also might not exist). For data privacy, one often has to restrict the data or parameter space anyway, in which case the Laplace distribution can be restricted so that the normalizing constant, mean, etc are all finite and well defined. The advantage in using such a mechanism is that sensitivity can be readily transferred into differential privacy. Theorem 1. Let f : Xn →M be a summary with global sensitivity ∆. Then the Laplace mechanism with footpoint f(D) and rate σ = 2∆/ satisfies -differential privacy. If the normalizing constant, Cη,σ does not depend on the footpoint, η, then one can take σ = ∆/ .
Proof. (Sketch) The proof follows from a direct verification via the triangle inequality.
In principle, it is possible to sample from the Laplace distribution using Markov Chain Monte Carlo (MCMC). One can start the chain at η and then make small proposed steps by randomly selecting a direction and radius on the tangent space. The resulting tangent vector can be pushed to an element ofM via the exponential map. Alternatively, if enough structure is known, as we illustrate in Section 5, one may be able to sample from the distribution directly or if the space is bounded then one can use rejection sampling.
An interesting alternative to our approach that is still inherently intrinsic is to instead generate the sanitised summary on a particular tangent space and then map it to the manifold using the exponential map. On the surface, this seems like a reasonable idea, however there are some subtle technicalities that would have to be overcome. In particular, one has to choose which tangent space to work with. Ideally one would work with the tangent space at the summary of interest, but that isn’t private. If another plane is used, then there is the chance for more serious distortions from the noise. Likely these issues could be overcome, but would require additional work.

heading: 4 Differentially Private Fréchet Means
text: As before, suppose the data consists of D = {x1, . . . , xn}, but now with xi ∈ M. The sample Fréchet mean, x̄, is defined to be the global minimizer of the energy functional
M3 x 7→ F2(x) := 1
2n n∑ i=1 ρ2(x, xi) ,
which is a natural generalization of the Euclidean mean to manifolds. Conditions that ensure existence and uniqueness of x̄ have been extensively studied since its inception in the 1970s [Karcher, 1977, Kendall, 1990]. Even when the mean is unique, the following example shows that the sensitivity need not decrease with the sample size, which produces sanitized estimates with low utility. Example 1. Let x1, . . . , xn be points on the unit circle M = S1 = {x ∈ R2 : ‖x‖ = 1} with arc-length distance ρ, represented as angles such that xi = 2πin−1 , i = 1, . . . , n− 1 and xn = xi′ for some i′ ∈ [n − 1]. Then minimizing F2 occurs when we take x̄ = xi. So, we can make the mean any of the xi by shifting a single point. If n is even, the furthest any two points can be is π and the resulting sensitivity is π, which clearly does not decrease with n.
The above example illustrates that one must have some additional structure to ensure that the sample Fréchet mean as a statistic is stable and that the sensitivity is properly decreasing with the sample size. The first requirement is that the sample Fréchet mean is unique; this imposes strong constraints on the spread of the data onM given by its curvature. Denote by Br(m) the open geodesic ball at m of radius r inM. For a given dataset D with points inM we make the following assumption. Assumption 1. The data D ⊆ Br(m0) for some m0, where r < r∗ := 12 min{injM, π 2κ −1/2} and κ > 0 is an upper bound on the sectional curvatures ofM.
For flat and negatively curved manifolds, Assumption 1 only says that the data lies in some bounded ball. In that case κ ≤ 0 and we can interpret κ−1/2 to be +∞. Furthermore, the injM can be arbitrarily large. Thus the radius of this ball only impacts the sensitivity, which is very common in
data privacy. However, it is, in general, difficult to relax Assumption 1 for positively curved manifolds even when privacy isn’t a concern. For example, it suffices to use the slightly weaker assumption r < 12 min{injM, πκ
−1/2} to ensure that (i) the closure B̄r(m0) is geodesically convex; (ii) x̄ exists and is unique; and (iii) x̄ belongs to the closure of the convex hull of points in D [Afsari, 2011]. For the unit sphere Sd−1 in Example 1 we have injM = π and κ = 1, we require r < π/2 so that a dataset D lying within a hemisphere on Sd−1 will have a unique sample Fréchet mean only if D contains no point lying on the equator. However, we need the stronger Assumption 1 to ensure that (x, y) 7→ ρ2(x, y) is convex along geodesics for every x and y inside Br(m0) [Le, 2001], which is required to determine the sensitivity of x̄.
In Theorem 2 we provide a bound on the global sensitivity of the Fréchet mean. The bound depends on the sample size n, the radius r of the ball that contains the data, and a function h(r, κ) which depends only on r and on the upper bound κ of the sectional curvatures ofM. For flat or negatively curved manifolds, we will see that h(r, κ) = 1, which matches classical results for the Euclidean space, owing to the classical Hadamard-Cartan theorem that states that a simply connectedM with non-negative sectional curvatures is diffeomorphic to Rd. However, the situation is more subtle for positively curved manifolds where h can no longer be ignored. Theorem 2. Under Assumption 1 consider two datasets D = {x1, . . . , xn−1, xn} and D′ = {x1, . . . , xn−1, x′n} differing by only one element. If x̄ and x̄′ are the two sample Fréchet means of D and D′ respectively, then
ρ(x̄, x̄′) ≤ 2r(2− h(r, κ)) nh(r, κ) , h(r, κ) =
{ 2r √ κ cot( √ κ2r) κ > 0;
1 κ ≤ 0 .
Proof. Consider the energy functionals F2 and F̃2 for the datasets D and D′ with unique sample means x̄ and x̄′, respectively. SinceM is complete the exponential map expx : TxM → M is surjective and under Assumption 1 the log map or inverse exponential mapM3 y 7→ exp−1x (y) ∈ TxM is well-defined for every x ∈ Br(m0). Here x 7→ ρ2(x, y) is twice continuously differentiable, and under Assumption 1 the function x 7→ ρ(x, y) is strictly convex for all x, y ∈ Br(m0) [Karcher, 1977, Afsari, 2011]. Consider an arc length parameterized, unit speed minimizing geodesic γ between x̄ and x̄′ such that γ(0) = x̄, γ(b) = x̄′ with b = ρ(x̄, x̄′). The composition, G2 := F2 ◦ γ : [0, b] → R is now a twice continuously differential real-valued function with derivatives Ġ2 and G̈2, and thus
Ġ2(b) = Ġ2(0) + bG̈2(t0) = ρ(x̄, x̄ ′)G̈2(t0),
for some 0 ≤ t0 ≤ b since Ġ2(0) = 0.
To determine G̈2(t0), we need to calculate the second derivative of ρ(γ(t0 + ), q)2 evaluated at = 0 and for an arbitrary q ∈ Br(m0), which equals
2 ( d d ρ(γ(t0 + ), q) ∣∣ =0 )2 + 2ρ(γ(t0), q) d2 d 2 ρ(γ(t0 + ), q) ∣∣ =0 .
Let βq be the angle between γ̇(t0) and α̇q(t0) formed in TzM, where αq is a minimizing geodesic from q to γ(t0); this implies that dd ρ(γ(t0 + ), q) ∣∣ =0
= 〈∇ρ(γ(t0), q), γ̇(t0)〉z = cosβq, with∇ as the Riemannian gradient, since ddtρ(γ(t), q)|t=0 = γ̇(ρ(x̄, q)) and γ is a unit-speed geodesic. As a consequence, with minimizing geodesics αxi from xi to z = γ(t0) and corresponding angles βxi ,
G̈2(t0) = d2
d 2 F2(γ(t0 + ))
∣∣ =0 = 1
n n∑ i=1 [ cos2 βxi + ρ(z, xi) d2 d 2 ρ(γ(t0 + ), xi) ∣∣ =0 ] .
Note that if z = xi for any i, the angle βxi is not well-defined, but regardless of the chosen path αxi the contribution to the sum from the particular xi will be one. The Hessian d 2 d 2 ρ(γ(t0 + ), xi) ∣∣ =0 of the distance function can be lower bounded using the Hessian comparison theorem [e.g. Lee, 2018, Theorem 11.7] to obtain
G̈2(t0) ≥ 1
n n∑ i=1 [cos2 βxi + a(ρ(z, xi), κ) sin 2 βxi ],
where
a(s, κ) =  s √ κ cot( √ κs) κ > 0; s−1 κ = 0;
s √ |κ| coth( √ |κ|s) κ < 0 .
If κ > 0, then (s, κ) 7→ a(s, κ) ≤ 1 and decreasing; on the other hand if κ < 0, (s, κ) 7→ a(s, κ) ≥ 1 and increasing. We hence have that
G̈2(t0) ≥ h(r, κ) :=
{ 2r √ κ cot( √ κ2r) κ > 0;
1 κ ≤ 0, (1)
since for κ = 0, a(s, κ) ≥ (2r)−1 and we can choose r = 1/2 since it is effectively unconstrained in this setting. The lower bound on G̈2(t0) thus depends on whetherM is positively or non-negatively curved depending on the sign of κ. This results in
ρ(x̄, x̄′) ≤ Ġ2(b) h(r, κ) = 1 h(r, κ) [Ġ2(b)− ˙̃ G2(b)],
where ˙̃G2(b) = 0 since∇F̃2(x̄′) = 0. For any x ∈M , in normal coordinates, the gradients
∇F̃2(x) = − 1
n [ n−1∑ i=1 exp−1x (xi) + exp −1 x (x ′ n) ] , ∇F2(x) = − 1 n [ n−1∑ i=1 exp−1x (xi) + exp −1 x (xn) ] ,
belong to TxM. This leads to the desired result since
ρ(x̄, x̄′) ≤ 1 nh(r, κ) ∥∥exp−1x̄′ (xn)− exp−1x̄′ (x′n)∥∥x̄′ ≤ 2r(2− h(r, κ))nh(r, κ) , using Lemma 1 in the Supplemental based on Jacobi field estimates [Karcher, 1977].
In our next Theorem we provide a guarantee on the utility of our mechanism. We demonstrate that, in general, the magnitude of the privacy noise added is O(dr/n ). Classic results on -DP [e.g. Hardt and Talwar, 2010] in Rd typically do not calculate sensitivity based on a Euclidean ball, instead focusing on privatizing each coordinate separately, in which case the optimal privacy noise is O(d3/2r/n ). To reconcile the two, the classic rate can equivalently be thought of as calculating sensitivity using an `∞ ball. It is easy to verify that it requires an `2 ball of radius r √ d to cover an `∞ ball of radius r, in which case the two rates agree, meaning that our mechanism is rate optimal. Theorem 3. Let the Assumptions of Theorem 2 hold. Let x̃ denote a draw from the Laplace mechanism conditioned on being in Br(m0). Assume that n and are such that σ → 0. Then x̃ is -DP and furthermore
E ρ(x̃, x̄)2 = O
( d2r2
n2 2
) .
Proof. First, x̃ conditioned on falling in Br(m0) guarantees the existence of the Laplace distribution (over Br(m0)) since it is now clearly integrable. That x̃ is DP follows from the same arguments as Theorem 1.
Turning to our utility guarantee, first notice that r was chosen such that there exists a set Ar ⊂ Tmo such that the expm0 : Ar → Br(m0) is a diffeomorphism. Identifying Tm0M with R
d, we also have that Ar is a ball centered at 0 with radius r (as measured using the inner product 〈·, ·〉m0 ). A change-of-variables, expm0(ṽ) = x̃, implies that ṽ has density equal to
f(v) = c−1f,σe −|v|/σ|Jv|
with support on Ar, where |Jv| is the determinant of the Jacobian of expm0 . This is not the K-norm distribution over Rd unless the Jacobian is constant in v. Since the set Ar is compact, the determinant
of the Jacobian is bounded from above and from below (away from 0), we can find constants c1 and c2, independent of σ, satisfying
c1e −|v|/σ ≤ e−|v|/σ|Jv| ≤ c2e−|v|/σ.
We can use this to bound the desired expected value as
E ρ(x̃, x̄)2 = c−1f,σ ∫ Ar |v|2e−|v|/σ|Jv| dv ≤ [ c1 ∫ Ar e−|v|/σ dv ]−1 c2 ∫ Ar |v|2e−|v|/σ dv.
Using a change of variables with u = σ−1v we have that
E ρ(x̃, x̄)2 ≤ c2σ 2
c1 [∫ Ar/σ e−|u| du ]−1 ∫ Ar/σ |u|2e−|u| du.
This, however, is simply σ22c2/c1 multiplied by the expected squared norm of a Euclidean Laplace that is conditioned on falling within Ar/σ. If we remove the condition that the Laplace falls within Ar/σ the value necessarily increases, thus we have that
E ρ(x̃, x̄)2 ≤ c2σ 2
c1 [∫ Rd e−|u| du ]−1 ∫ Rd |u|2e−|u| du.
Using a change of variables in both integrals to spherical coordinates and noting σ = O(r/n ) as in Theorems 1 and 2, we get that
E ρ(x̃, x̄)2 ≤ c2σ 2
c1 [∫ ∞ 0 yd−1e−y dy ]−1 ∫ ∞ 0 yd+1e−y dy = c2σ 2d(d− 1) c1 = O ( d2r2 n2 2 ) .
Our final Theorem focuses on the case of linear manifolds to highlight mathematically the benefit of constructing the privacy mechanism directly on the manifold as opposed to an ambient space and then projecting back onto the manifold. Practically, the variance of the mechanism is inflated by a factor of D/d where d and D are the dimensions of the manifold and ambient space respectively with d ≤ D. Intuitively, if the ambient space is of a higher dimension than the manifold, then one has to expend additional privacy budget to privatise the additional dimensions, which our approach avoids. Since the mechanism concentrates around a single point as the sample size grows (and thus one can use a single tangent space to parameterise the problem), one should be able to extend this to more general manifolds, though for ease of exposition we focus on the linear case. Theorem 4. LetM⊂ RD be a d-dimensional linear subspace of RD equipped with the Euclidean metric. Assume the assumptions of Theorem 3 hold. Let x̃D denote the private summary generated from the Laplace over RD with scale σ, in the sense of Definition 3 (equivalently, this is the K-norm mechanism with the `2 norm). Then
E ‖PMx̃D − x̄‖2 = O ( dDr2
n2 2
) ,
where PM is the projection operator ontoM.
Proof. Choose {v1, . . . , vD} as an orthonormal basis of RD such that the matrix PM = V V T is an orthogonal projector ontoM, where V = [v1, . . . , vd] ∈ RD×d and V TV = Id. Let 〈·, ·〉 be the usual inner product on RD with norm ‖ · ‖. Note that x̃D = x̄ + σRU, where R distributed as Gamma(D, 1), U is uniform on SD−1 (i.e. 〈U,U〉 = 1) independent of R, and x̄ ∈ M (see Supplemental 1.4 for details). Then
E ‖PMx̃D − x̄‖2 = E ‖σRPMU‖2 = σ2 E[R2] d∑ i=1 E[〈U, vi〉2].
Since U is uniform on the sphere, it follows that the vector (〈U, v1〉2, . . . , 〈U, vD〉2) follows a Dirichlet distribution with concentration parameters all equal to 1. Thus, for each i, 〈U, vi〉2 is distributed as Beta(1, D − 1). This completes the proof as
E ‖PMx̃D − x̄‖2 = dσ2(D +D2) 1
D = σ2d(D + 1).

heading: 5 Examples
text: In this section we numerically explore two examples that are common in statistics. In the first example, we consider the space of symmetric positive definite matrices (SPDM) equipped with the Rao-Fisher affine invariant metric, under which the space is a negatively curved manifold. In the second example we consider data lying on the sphere, which can be used to model discrete distributions or compositional data, as an example of a positively curved manifold. In both cases we demonstrate substantial gains in utility when privatizing the Fréchet mean using our proposed mechanism against using a more standard Euclidean approach that utilizes an ambient space. Other examples of Riemannian manifolds that commonly arise in statistics include the Steifel manifold for PCA projections, quotient spaces for modeling shapes, and hyperbolic spaces for modelling phylogenetic structures. Simulations are done in Matlab on a desktop computer with an Intel Xeon processor at 3.60GHz with 31.9 GB of RAM running Windows 10. Additional details on each example are also provided in the supplemental.

heading: 5.1 SPDM
text: Let P(k) denote the space of k× k symmetric positive definite matrices. In addition to being used for modeling covariance matrices, this space is widely used in engineering of brain-computer interfaces [Congedo et al., 2017], computer vision [Zheng et al., 2014], and radar signal processing [Arnaudon et al., 2013].
We consider P(k) equipped with the Riemannian metric 〈v, u〉p = Tr(p−1up−1v), known as the Rao-Fisher or Affine Invariant metric where u, v ∈ TpP(k) are symmetric matrices. The metric makes P(k) into a manifold with negative sectional curvature [e.g., Helgason, 2001]. Under this metric expp(v) = p 1/2Exp ( p−1/2vp−1/2 ) p1/2 and exp−1q (p) = q 1/2Log ( q−1/2pq−1/2 ) q1/2, where Exp and Log are matrix exponential and logarithm respectively. The (squared) distance between q, p ∈ P(k) is given in closed form by ρ2(q, p) = Tr[Log(q−1/2pq−1/2)2]; these expressions are widely available [e.g., Hajri et al., 2016, Said et al., 2017]. To calculate the Fréchet mean of a sample we use a gradient descent algorithm as proposed in Le [2001] (Supplemental material 1.1).
We generate samplesD = {x1, · · · , xn−1, xn} from P(k) using the Wishart distribution as discussed in the Supplemental material 1.2.1. In the first and second panels of Figure 1 we show simulation results which illustrate Theorems 2 and 3 and compare the utility of the Euclidean counterpart. In the first panel we illustrate the sensitivity by plotting ρ(x̄, x̄′), for neighboring databases, as blue dots as well as the theoretical bound The blue line is the average distance at each sample size.
We compute the Fréchet mean x̄ and then privatize the mean using two approaches: (i) we generate the privatized mean x̃ by sampling using the Laplace distribution from Definition 3 defined directly on P(k) with footpoint x̄ using the algorithm in Hajri et al. [2016]; (ii) we use the embedding of P(k) into the set of k × k symmetric matrices, isomorphic to Rk(k+1)/2, to represent x̄ as a vectorized matrix vech(x̄) in R3 (i.e., k = 2) and obtain a privatized mean vech(x̃E) by adding to vech(x̄) a vector drawn from the standard Euclidean Laplace. Then x̃E is obtained by reverting to the matrix representation, which is not guaranteed to stay in P(k) but is symmetric by construction. In the second panel we plot the average, across repetitions, of the distances ‖vech(x̄) − vech(x̃)‖ (blue) and ‖vech(x̄)− vech(x̃E)‖ (red). Since the Euclidean summary need not belong toM, using the Euclidean distance enables a common comparison between the two methods. The shaded regions around the lines correspond to ±2SE, where SE is the standard error of the average distances. Examining panel 2 in Figure 1, we see that our approach has better utility even when calculated using the Euclidean distance. Here, the ambient space approach does not increase the dimension of the statistic since P(k) and the space of symmetric matrices are of the same dimension, thus the gain in utility appears to be primarily due to respecting the geometry of the problem. Furthermore, as expected for smaller sample sizes, approach (ii) can produce summaries that are not positive definite, with about 25% not being in P(k) at sample sizes 20-40.

heading: 5.2 Spheres
text: Let Sdκ denote a d-dimensional sphere of radius κ−1/2 parameterized such that the sectional curvature is constant κ > 0. Identifying the sphere as a subset of Rd+1, the tangent space at p is TpSdκ =
{v ∈ Rd+1 : 〈v, p〉 = 0}. The exponential map, defined on all of TpSdκ, is given by expp(v) = cos(‖v‖)p+ κ−1/2 sin(‖v‖)v/‖v‖ and expp(0) = p. The inverse exponential map exp−1p : Sdκ → TpSdκ is defined only within the ball of radius strictly smaller than π/2 around p and is given by exp−1p (q) = θ sin(θ) (q − cos(θ)p) with q 6= {p,−p}. The corresponding distance function is ρ(p, q) = θ where θ = cos−1(〈p, q〉) and p, q ∈ Sdκ. In similar fashion to Section 5.1, in the third and fourth panels of Figure 1 we show simulation results which illustrate Theorems 2 and 3, and compare utility to its Euclidean counterpart. For our simulations we fix d = 2 and κ = 1. Consistent with Assumption 1 on support of the data, we choose a ball Br(m0) of radius r = π/8 < r∗ = π/2 and take m0 as the north pole. We generate random samples as shown in the Supplemental material 1.3.1. The red line again corresponds to the theoretical bound from Theorem 2. The (unique) Fréchet mean x̄ is computed using a gradient descent algorithm.
In the last panel we compare utility of the privatized means, again obtained using two approaches: (i) exactly as in approach (i) with SPDM; (ii) using the embedding of x̄ into R3 (Cartesian coordinates) to represent x̄ and obtaining a private x̃E by adding to x̄ a draw from the Euclidean Laplace. We display the average Euclidean distance between the mean and private mean ±2SE, where SE is the standard error of the distances at each sample size; we use 1000 replicates at each sample size. The blue band is obtained using approach (i) and the red band using approach (ii). While the contrast between (i) and (ii) is not as stark as with SPDM, our approach still produces about 15% less noise, with an average of 16.8% reduction in the smaller sample sizes and 12% reduction in the larger sample sizes. Furthermore, unlike in the SPDM case, the Euclidean private summary is never on the manifold sinceM as a subset of Rd+1 has measure zero.

heading: 6 Conclusions and Future Work
text: In this paper we have demonstrated how to achieve pure differential privacy over Riemannian manifolds by relying on the intrinsic structure of the manifold instead of the structure induced by a higher-dimensional ambient or embedding space. Practically, this ensures that the private summary preserves the same geometric properties as the non-private one. Theorem 3 shows that our mechanism matches the known optimal rates for linear spaces for the Fréchet mean, while Theorem 4 highlights the benefit of the intrinsic approach in contrast to projection-based ones using an ambient space.
The benefits of directly working on the manifold come at the expense of a more complicated mathematical and computational framework, challenges with which vary between different manifolds depending on availability of closed-form expressions for geometric quantities. Conversely, working in a linear ambient space is usually computationally simpler, although other issues abound: projecting onto the manifold may not be possible (e.g. SPDM under negative curvature), or the projection may lead to poor utility in high-curvature places on the manifold.
As is well appreciated in the mathematics/statistics literature, there are challenges that are unique to positively curved spaces, which we also encounter here. The central issue is that the squared distance function need not be convex over large enough areas, and strong restrictions on the spread of the data are required to ensure that underlying summaries are unique and well defined. This phenomenon manifests in the form of a correction term h in our results, which impacts sensitivity of the statistics. Numerical illustrations in Section 5 show that this is not just a technical oddity or gap in our proofs since our empirical sensitivity can be seen to be quite close to our theoretical bound.
As this is the first paper we are aware of in DP over general manifolds, there are clearly many research opportunities. A deeper exploration over positively curved spaces would be useful given how common they are in practice (e.g., landmark shape spaces). A class of spaces unexplored in this paper, but well worth investigating, are Hadamard or (complete) CAT(0) spaces. They are non-positively curved metric spaces that need not be manifolds, on which geodesics and Fréchet means can be computed, and represent a natural geometric setting for graph and tree-structured data. One could also extend any number of privacy tools to manifolds including the Gaussian mechanism, exponential mechanism, objective perturbation, K-norm gradient mechanism, approximate DP, concentrated DP, and many others.
Acknowledgments and Disclosure of Funding
This work was funded in part by NSF SES-1853209, the Simons Institute at Berkeley and their 2019 program on Data Privacy to MR; and, NSF DMS-2015374, NIH R37-CA214955 and EPSRC EP/V048104/1 to KB. We thank Huiling Le for helpful discussions.

Review summaries:
This paper is concerned with extending the Laplace mechanism (and its guarantees) beyond Euclidean spaces, in order to achieve pure differential privacy over manifolds without resorting to first embed the space into an Euclidean space. This seems like an interesting idea, which can lead to some savings in the parameters (by constant factors, see below) at the price of some computational issues and (as far as I can tell) a more restrictive assumption on the data (Assumption 1).  The authors specifically consider the question of bounding the sensitivity of the Frechet mean (Theorem 2), in order to compute privately (pure privacy) the Frechet mean of a dataset (Theorem 3) matching the corresponding guarantee for Euclidean spaces.

This paper studies the problem of releases a differentially private statistic when the data and the statistic lie on a Riemannian manifold. They give a formulation of the Laplacian mechanism in this setting. They then consider the specific example of.a Frechet mean to display the benefits of using the manifold structure rather than adding noise in the ambient (higher dimensional) space, then projecting back to the manifold.   The authors show through specific examples and experiments the benefit of using the manifold structure rather than adding noise in the ambient space. As a first paper on the topic this seems like a nice proof of concept. I would have liked it to go deeper in discussing the complexity of actually sampling from the Laplace mechanism on a manifold. In particular, assumption 1 felt like it simplified the problem considerably.

There are two main contributions in the paper: - First, the language of differential privacy is adapted to the framework of Riemannian manifolds. In particular, an adaption of the Laplace mechanism on manifolds is proposed. - Second, the Laplace mechanism is applied to produce differentially private Fréchet means on a manifold. The key result to do so is to provide a stability result for Fréchet means. The curvature of the manifold is shown to play a central role in this stability property.