Wrapping query perturbation mechanisms in ggplot stats: Towards privacy-integrated plotting

Graphical displays are a major element in exploratory data analysis and important for communicating one's findings. But when working with confidential data, like personally identifying information, or data that is subject to business secrets, compliant plotting can be difficult, as graphed data points can reveal information on individual observations from the data. Analysts working in research data centers, for instance, are commonly subject to the following "rule of thumb" (taken from Bond et al., 2015, sect. 2.3):

"Graphs in themselves are not permitted as output. The underlying (aggregated) data may be submitted as output, which when cleared may be used to reconstruct graphs in the researcher’s own environment."

What this means is that computing aggregates and graphing them is split in separate procedures. This slows down data exploration and may discourage the use of graphs overall, potentially risking the loss of their analytical power.

Would it not be nice of one could wrap the computation of confidentialised statistics directly in the graphics software? So that one can use established EDA workflows with minimal effort added for ensuring output confidentiality? 

In this short proof of concept I show how the use of many of the common plot types from popular R package ggplot2 (Wickham et al.) can be made privacy-compliant by wrapping query perturbation mechanisms from the Differential Privacy (DP) and Cell Key Method (CKM) toolkits.

A general strategy: Treating collective geoms like tabular aggregates

While it seems sensible that individual geoms, which put an identifiable data point per unit, are considered a-priori non-confidential, the use of collective geoms, which plot an aggregate or statistical summary metric, are strongly analogus to the counts and sums shown in tables, for which established confidentiality mechanisms exist. 

This analogy extends our initial goal to two specific desiderata:

1. The workflow for producing the plot should be encumbered as little as possible. Ideally, we don't perform any additional data transformation by hand and just have to change the plot call in minor ways.
2. The mechanism used for ensuring confidentiality (typically a perturbation mechanism) replicates that used by the data owner for other outputs (e.g. tables) and produces results that are consistent with them.

[The second point rarely shows up in the data privacy literature, but is an important requirement in production settings for canonical publication, e.g. by state statistical offices.]

The general tactic is binning data points, but, as Avraam et al. (2021) point out, that is not enough. Geoms representing small counts or those whose values can be combined with other published aggregates to reveal information on small numbers of individuals, are still problematic.

Query perturbation then, in essence, means that instead of the real bin size we aim to plot the bin size plus a small-ish (pseudo-)random perturbation value, drawn from a designated unbiased noise distribution. For this PoC, two popular noise mechanisms are considered: 

• the Laplace mechanism, associated with DP (Dwork et al., 2006), and
•  the sticky-noise mechanism associated with CKM (Fraser & Wooton, 2005).

[Another alternative would be a random sample mechanism, as described in this previous post.]

A preliminary example: count_dp & count_ckm stats

Let's start with a simple bar chart. The following is a standard example, adapted from the ggplot2 reference:

ggbar <- ggplot(data = mpg, aes(class))
ggbar + geom_bar()

Adding some labels to make the class sizes easier to read:

ggbar + geom_bar() +
    geom_text(stat = "count", aes(label = after_stat(count)), vjust = 1.5, color = "white")

ggplot builds on stats, which are computations done on the data during the creation of the plot. Here, it simply counts the number of observations by class. We see reference to this "count" stat explicitly in the call to geom_text above (which allows us to plot the number of observations, as they have to be computed from the micro data).

But already the command

ggbar + geom_bar()

is equivalent to

ggbar + geom_bar(stat = "count")

which, if not specified explicitly, ggplot does under the hood.

The stats are our ticket to privacy-compliant plotting. Instead of using the default count stat, we can create a custom stat that first computes the count per class (as before), but then applies a perturbation to each result before displaying it in the chart.

[The details on how the new stats were created are omitted here. I am indebted to the chapter on 'Building a Stat' from the book by Peng et al. (2020), to the 'Extending ggplot2' section from the ggplot2 online documentation as well as to the open code on the ggplot2 GitHub.]

For a differentially private option, we may consider wrapping a call to the Laplace mechanism from the package DPpack in our new "count_dp" stat. This requires us to pass a parameter $\epsilon$ (also called the privacy budget for this query) to DPpack's LaplaceMechanism function. For convenience, we allow passing $\epsilon$ directly in the ggplot call. For reproducibility, we add the option to pass a seed for the RNG as well (important if we want the numbers shown to represent the actual bar heights, rather than a different, independent perturbation of the original count).

ggbar + geom_bar(stat = "count_dp", eps = 1, seed = 42) +
    geom_text(stat = "count_dp", aes(label = after_stat(count)), eps = 1, seed = 42, 
    		  vjust = 1.5, color = "white")

What is happening under the hood? 

In the normal bar chart, all observations $i = 1, \ldots, N$ are grouped in the assigned classes and the count $N_C = |\{i \in C\}|$ for class $C$ determines the height of the bar. With a query perturbation in the loop, however, we instead display $N_C + \Delta$, where for the Laplace mechanism $\Delta \sim \text{Laplace}(0, s / \epsilon)$. $s$ is the query sensitivity, set here to the default in DPpack's tableDP function for count queries. $\epsilon$ determines the signal-to-noise ratio in the differentially private output, with small $\epsilon$ resulting in noisier answers.

The CKM noise mechanism, on the other hand, puts a stronger focus on consistency of outputs across different queries and releases. The idea is to add a random key $r_i \in [0, 1]$ permanently to each record in the data set (called a record key). When calculating $N_C$, our new "count_ckm" stat must also derive a new compound key (cell key) as $r_C = (\sum_{i \in C} r_i) \text{ mod } 1$ in the background. The noise distribution is then sampled precisely at point $r_C$ to derive $\Delta$, which is added to $N_C$ as before. This leads to queries of the same contributing records automatically receiving the same noise. To use this mechanism we begin by affixing record keys to the input data frame.

mpg$rk <- runif(nrow(mpg), 0 , 1)

Our count_ckm stat should also allow us to pass a custom noise distribution. These can be created using the ptable package.

library(ptable)
ptab <- create_cnt_ptable(D = 5, V = 2.01, js = 2)@pTable

[The parameters D, V and js control the shape of the noise distribution, including the amount of noise added. For their interpretation, check out the ptable vignette.] 

The benefit of all this preparation is that now, instead of having to set seeds again and again, we simply pass the name of the record key variable to our plotting function.

ggbar <- ggplot(data = mpg, aes(class, rkey = rk))
ggbar + geom_bar(stat = "count_ckm", pTable = ptab) +
  geom_text(stat = "count_ckm", aes(label = after_stat(count)), pTable = ptab,
            vjust = 1.5, color = "white")

Naturally, we can extend the binning to include subject-level groupings. In this case, each group will be perturbed independently.

# stacked bar chart - original
ggplot(data = mpg, aes(class, fill = drv)) +
  geom_bar() +
  geom_text(stat = "count", aes(label = after_stat(count)), 
            vjust = 1.5, color = "white", position = "stack") +
  ggtitle("unperturbed")

# stacked bar chart - CKM
ggplot(data = mpg, aes(class, fill = drv, rkey = rk)) +
  geom_bar(stat = "count_ckm", pTable = ptab) +
  geom_text(stat = "count_ckm", aes(label = after_stat(count)), pTable = ptab, 
            vjust = 1.5, color = "white", position = "stack") +
  ggtitle("perturbed by CKM")

 

General histograms: bin_db & bin_ckm stats

We are not confined to fixed classes though. The ideas work just as well for binning with heterogeneous bin widths, as used in histograms. After writing "bin_dp" and "bin_ckm" stats to replace the original ggplot bin stat, the modifications needed in the plot call are the same as before:

data("diamonds")
diamonds$rk <- runif(nrow(diamonds)) # adding record keys
price_bins <- quantile(diamonds$price, probs = seq(0, 1, length = 11))

# histogram - original
ggplot(diamonds, aes(price, after_stat(count / width))) +
  geom_histogram(breaks = price_bins, fill = "white") +
  geom_text(stat = "bin", breaks = price_bins, aes(label = after_stat(count))) +
  ggtitle("unperturbed")

# histogram - DP
ggplot(diamonds, aes(price, after_stat(count / width))) +
  geom_histogram(stat = "bin_dp", breaks = price_bins, fill = "white",
                 eps = 1, seed = 23) +
  geom_text(stat = "bin_dp", breaks = price_bins, aes(label = after_stat(count)),
            eps = 1, seed = 23) +
  ggtitle("perturbed by Laplace")  
  
# histogram - CKM
ggplot(diamonds, aes(price, after_stat(count / width), rkey = rk)) +
  geom_histogram(stat = "bin_ckm", breaks = price_bins, fill = "white", pTable = ptab) +
  geom_text(stat = "bin_ckm", breaks = price_bins, aes(label = after_stat(count)),
            pTable = ptab) +
  ggtitle("perturbed by CKM")

The same stat can, of course, be used for frequency polygons:  

# frequency polygon - original
ggplot(diamonds, aes(price)) +
  geom_freqpoly(bins = 10, lty = "dotted", color = "blue") +
  geom_text(stat = "bin", bins = 10, aes(label = after_stat(count))) +
  ggtitle("unperturbed")

# frequency polygon - DP
ggplot(diamonds, aes(price)) +
  geom_freqpoly(stat = "bin_dp", bins = 10, lty = "dotted", color = "blue",
                eps = .5, seed = 23) +
  geom_text(stat = "bin_dp", bins = 10, aes(label = after_stat(count)),
            eps = .5, seed = 23) +
  ggtitle("perturbed by Laplace")

# frequency polygon - CKM
ggplot(diamonds, aes(price, rkey = rk)) +
  geom_freqpoly(stat = "bin_ckm", bins = 10, lty = "dotted", color = "blue",
                pTable = ptab) +
  geom_text(stat = "bin_ckm", bins = 10, aes(label = after_stat(count)),
            pTable = ptab) +
  ggtitle("perturbed by CKM")

What to do about scatter plots? bin2d_dp, bin2d_ckm, binhex_dp, binhex_ckm, sum_dp & sum_ckm stats

Scatter plots are great. Unfortunately, they are also patently privacy-incompatible, since each data point shows the exact measurements associated with a single (potentially identifiable) unit of observation. The binning + perturbing strategy needs to apply here too. In fact, binning to what is, in essence, a 2D histogram, is already an established approach, though not necessarily in connection with ensuring confidentiality, but rather to handle overplotting. Carr et al. (1987) are often credited with pioneering the idea.

In ggplot2, we have to draft new stats (called bin2d_dp and bin2d_ckm) from the bin2d stat. 

# 2D histogram - original
ggplot(mpg, aes(cty, hwy)) +
  geom_bin_2d() +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("unperturbed")

# 2D histogram - DP
ggplot(mpg, aes(cty, hwy)) +
  geom_bin_2d(stat = "bin2d_dp", eps = .5, seed = 1234) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("perturbed by Laplace")

# 2D histogram - CKM
ggplot(mpg, aes(cty, hwy, rkey = rk)) +
  geom_bin_2d(stat = "bin2d_ckm", pTable = ptab) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("perturbed by CKM")

  

A maybe more elegant version of binning in two dimensions is the bin_hex stat, which creates honeycomb-like bin shapes. Thompson et al. (2013) suggest this plot type as a surrogate for scatter plots in sensitive data settings. My new stats are called binhex_dp and binhex_ckm:

# hex histogram - original
ggplot(mpg, aes(cty, hwy)) +
  geom_hex() +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("unperturbed")

# hex histogram - DP
ggplot(mpg, aes(cty, hwy)) +
  geom_hex(stat = "binhex_dp", eps = .5, seed = 1234) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("perturbed by Laplace")

# hex histogram - CKM
ggplot(mpg, aes(cty, hwy, rkey = rk)) +
  geom_hex(stat = "binhex_ckm", pTable = ptab) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("perturbed by CKM")

If the measurement accuracy is limited, so that several points already share positions, no extra binning may be needed and the count of data points at each position could be perturbed directly. To do this, the new stats sum_dp and sum_ckm are derived from the sum stat, used with geom_count:

# point counts - original
ggplot(mpg, aes(cty, hwy)) +
  geom_count() +
  scale_size_area(breaks = c(3, 6, 12)) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("unperturbed")

# point counts - DP
ggplot(mpg, aes(cty, hwy)) +
  geom_count(stat = "sum_dp", seed = 1000, eps = 0.5) +
  scale_size_area(breaks = c(3, 6, 12)) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("perturbed by Laplace")

# point counts - CKM
ggplot(mpg, aes(cty, hwy, rkey = rk)) +
  geom_count(stat = "sum_ckm", pTable = ptab) +
  scale_size_area(breaks = c(3, 6, 12)) +
  xlim(c(8, 38)) + ylim(c(10, 45)) +
  ggtitle("perturbed by CKM")

Outlook

If confidentiality-preserving data work is to catch on further, it has to play within the confines of established data science practices like exploratory data analysis (Snoke et al., 2024). Ideally, making plots respect the data owner's confidentiality rules should become seamless - or at least as close to seamless as possible. Allowing well-established plotting packages to create provably privacy-respecting charts with minimal effort would be an important step towards that goal.

Implementation remarks

Stats were adapted from the original ggplot2 source code available on GitHub, by adding the relevant perturbation steps to internal compute_group or compute_panel functions. It is planned to provide the new stats in a designated R package for privacy-preserving plotting in the future.

Literature

D. Avraam, R. Wilson, O. Butters, T. Burton, C. Nicolaides, E. Jones, A. Boyd, P. Burton, "Privacy preserving data visualizations," EPJ Data Science, vol.10, no.2, 2021.

S. Bond, M. Brandt, P.-P de Wolf, "Guidelines for the checking of output based on microdata research," Data without Boundaries (DwB) Project Deliverable, WP11: Improved Methodologies for Managing Risks of Access to Detailed OS Data, 2015.

D.B. Carr, R.J. Littlefield, W.L. Nicholson, J.S. Littlefield, "Scatterplot matrix techniques for large N," Journal of the American Statistical Association, vol.82, no.398, pp.424-436, 1987.

C. Dwork, F. McSherry, K. Nissim, A. Smith, "Calibrating Noise to Sensitivity in Private Data Analysis," Third Theory of Cryptography Conference TCC '06, New York, USA, 2006.

B. Fraser & J. Wooton, "A proposed method for confidentialising tabular output to protect against differencing," Joint UNECE / Eurostat work session on statistical data confidentiality, Geneva, Switzerland, 2005.

R.D. Peng, S. Kross, B. Anderson, Mastering Software Development in R, 2020.

J. Snoke, C. McKay Bowen, A.R. Williams, A.F. Barrientos, "Incompatibilities between current practices in statistical data analysis and Differential Privacy," Journal of Privacy and Confidentiality, vol.14, no.3, 2024.

G. Thompson, S. Broadfoot, D. Elazar, "Methodology for the automatic confidentialisation of statistical outputs from remote servers at the Australian Bureau of Statistics," Joint UNECE/Eurostat Work Session on Statistical Data Confidentiality '13, Ottawa, Canada, 2013.

H. Wickham, D. Navarro, T. Lin Pedersen, ggplot2: Elegant Graphics for Data Analysis (3rd Ed.), in progress.