See without looking: joint visualization of sensitive multi-site datasets

1 Introduction

Large-scale datasets have proven to be unreasonably effective in facilitating solutions to many difficult machine learning problems (Halevy et al., 2009). High tolerance to mistakes and possible problems with individual data samples in applications relevant to internet businesses, 1 together with advances in machine learning methodologies (such as deep learning (Goodfellow et al., 2016)) are able to effectively average out problems with individual samples lead to improved performance in recognition tasks. The story is different in the domains working with biomedical data, such as neuroimaging, where a data sample is a magnetic resonance image (MRI) of the entire brain containing on the order of 100,000 volumetric pixels (voxels). There, the data collection process for each sample is expensive, considerations of data privacy often prevent pooling data collected at multiple places, thus the datasets are not as large. Yet they are large enough to be difficult to manually vet each sample. An incorrect data sample may still lead to wrong conclusions and quality control is an important part of every analysis. Methods for simultaneous embedding of multiple samples are welcomed, as it is very difficult to scan through each and are used in practice for quality control (Panta et al., 2016).

∗This work was supported by NIH R01DA040487; SP was in part supported by NSF IIS-1318759

1 It is expected that the mistakes average out and even if they do not the cost of displaying an image that the user did not request is low.

A common way of visualizing a dataset consisting of multiple high dimensional data points is embedding it to a 2 or 3-dimensional space. Such embedding can be an intuitive exploratory tool for quick detection of underlying structure and outliers. Although linear methods such as principal component analysis (PCA) provide the functionality they are not usually useful when there is a need to preserve and convey hidden nonlinear structure in the data. Many methods were developed for the task on nonlinear data embedding and visualization including Sammon mapping (Sammon Jr, 1969), curvilinear components analysis (Demartines and H´erault, 1997), Stochastic Neighbor Embedding (Hinton and Roweis, 2002), Isomap (Tenenbaum et al., 2000), Maximum Variance Unfolding (Weinberger and Saul, 2006), Locally Linear Embedding (Roweis and Saul, 2000), Laplacian Eigenmaps (Belkin and Niyogi, 2003). The problem with these approaches is in their inability to retain local and global structure in a single map. A method to handle this situation efficiently was alternatively proposed: t-distributed stochastic neighbor embedding (t-SNE) (Van der Maaten and Hinton, 2008). The embeddings resulting from t-SNE applications are usually intuitive and interpretable which makes it an attractive tool for domain scientists (Bushati et al., 2011; Panta et al., 2016).

All of these methods, however, are built on the assumption that the input dataset is locally available. t-SNE, for example, needs computation of the pairwise distances between all samples (points) in the dataset. Yet, in many situations it is impossible to pull the data together into a single location due to legal and ethical considerations. This is especially true for biomedical domain, where the risk of identifiability of anonymized data often prevents open sharing. Meanwhile, many systems allow virtually pooling datasets located at multiple research sites and analyzing them using algorithms that are able to operate on decentralized datasets (Carter et al., 2015; Gaye et al., 2014; Plis et al., 2016). The importance of operating on sensitive data without pooling it together and thus generating truly large-scale neuroimaging datasets is so high that researchers successfully engage into manually simulating a distributed system (Thompson et al., 2014). For all of the applications of the above systems quality control is essential and intuitive visualization of the complete virtual dataset physically spread across multiple locations is an important and much needed tool for filtering out participating sites with bad data, incorrect processing or simply mistakes in the input process.

In this paper we propose a way to embed into a 2D map a decentralized dataset that is spread across multiple locations such that the data at each location cannot be shared with others due to e.g. privacy concerns. Even the methods that are seemingly suitable to this setting (after a possible modification) do not seem to address the problem of inability to compute the distance between samples located at different sites. For example, Globerson et al. (2007) suggested a method of embedding multiple modalities into the same space. We could think of the modalities as our locations and modify their approach to our settings. This, however, is not straightforward as, again, the approach requires measuring co-occurrence, which transcends the borders of local sites.

We base our approach on availability of public anonymized datasets, which we use as a reference and build the overall embedding around it. This is most similar to the landmark points previously used for improving computational efficiency (De Silva and Tenenbaum, 2003, 2004). In fact, we start with a method that resembles the original landmark points approach. We show that it is not as flexible and does not produce good results. Then we introduce a dynamic modification that indeed is able to generate an embedding that reflects relations between points spread across multiple locations. Unlike the original landmark point approach, we use t-SNE as our base algorithm. We call the decentralized data algorithm dSNE. To evaluate the performance and to compare with the centralized version we use the MNIST dataset (LeCun et al., 1998) and taking advantage of the known classes of the samples introduce a metric of overlap and roundness to quantify the comparisons. We evaluate and compare our algorithms in a range of various settings, establishing 4 experiments with MNIST data. Furthermore, we apply our approach to a truly multi-site neuroimaging dataset: ABIDE (Di Martino et al., 2014). 2

2 Methods

In the general problem of data embedding we are given a dataset of N points X = [x1 . . . , xN], such that each point

xi ∈ Rn with the task of producing another dataset Y = [y1 . . . , yN], such that each point yi ∈ Rm, where m << n. Usually m = 2 for convenience of visualization. Of course, this is an incomplete definition of the problem, as for any interesting results Y must be constrained by X.

2.1 t-SNE

In t-SNE the distances between the points in Y must be as close to the distances between same points in X as possible given the weighted importance of preserving the relations with nearby points over those that are far. To achieve that, tSNE first converts the high dimensional Euclidean distances between datapoints into conditional probability that represent similarities (see Algorithm 1).

Algorithm 1 PairwiseAffinities

Input: p (site index), ρ (perplexity), X ∈ RN×C×Kp Output: P

1: Eq. (1) to compute pij with perplexity ρ

The similarity of datapoint xj to datapoint xi is the conditional probability, pj|i, that xi would pick xj as its neighbor if neighbors were picked in proportion to their probability density under a Gaussian centered at xi.

Algorithm 2 tSNE

Input:

Data: X = [x1, x2 . . . xN ], xi ∈ Rn

Objective parameters: ρ (perplexity)

Optimization parameters: T (number of iterations), η (learning rate), α (momentum)

1: {pij} = PairwiseAffinites 0, ρ,X

• 2: Y ∝ N(0, 10-4 I), I ∈ Rm×m initialize from Gaussian
• 3: for i = 1 to T do

4: Eq. (2) to compute low-dimensional affinities qij

• 5: Eq. (3) to compute δC / δyi

7: end for

At the same way we can compute qij from low dimensional output data. Algorithm 2 outlines the full procedure. In this algorithm, high dimensional and low dimensional pairwise affinities are formulated using equation (1) and (2) respectively. The gradient of the Kullback-Leibler divergence between P and the Student-t based joint probability distribution Q is expressed in (3).

Inspired by the overall quite satisfactory performance of t-SNE on a range of tasks, we base on it our algorithms for decentralized datasets.

2.2 Single-shot dSNE

In a scenario that we consider here, no data can leave a local site and thus it seems impossible to compute distances of samples across the sites. Without those distances (see equation (1)) we will not be able to obtain a common embedding. Fortunately, in neuroimaging not all data is private and unshareable. Public repositories of MRI data are popular and they provide diverse datasets for analyses (Castellanos et al., 2013; Hall et al., 2012; Ivory, 2015).

First we introduce some notation for sending and receiving messages. For a matrix X, X→p means that it is sent to site p and X←p means that it is received from site p. We assume that a shared dataset is accessible to all local sites and the sites have it downloaded.

Algorithm 3 singleshotDSNE

Input:

Objective parameters: ρ (perplexity)

Optimization parameters: T , η, α

Shared Data: Xs = [xs1 , xs2 . . . xsNs ], xsi ∈ Rn

Data at site p∀p: Xp = [xp1 , xp2 . . . xpNp ], xpi ∈ Rn

• 1: Ys ←tSNEXs , ρ, T, η, α
• 2: for p = 0 to P do
• 3: Y→p s
• 4: Run tSNE on [Xp,Xs]

▷ At the master node

• ▷ At local site p
• 5: At each iteration only update Yp ▷ At local site p
• 6: end for
• 7: Y← []
• ▷ At the master
• 8: for p = 0 to P do

▷ At the master

• 9: Y←p p
• 10: Y← [Y,Yp]
• 11: end for
• 12: Y← [Y,Ys]

For Single shot d-SNE (Algorithm 3) we at first pass the reference data from centralized site C to each local site.

Now each local site's data consists of two portions. One is its local dataset, for which we need to preserve privacy, and another one is the shared reference dataset both comprise the combined datasets. Each local site runs the t-SNE algorithm on this combine data and produces an embedding into a low dimensional space. However, while computing each iteration of tSNE a local site computes gradient based on combined data, but it only updates the embedding vectors y for local datasets. The embedding for the shared data has been precomputed at the master node and shared with each local site. Similarly to the landmark points approach of De Silva and Tenenbaum our method uses reference points to tie together data frommultiple sites. In practice the samples in the shared dataset are not controlled by the researchers using our method, and it is hard to assess the usefulness of each sample in the shared data in advance. In the end each local site obtains an embedding of its data together with the embedding of the shared dataset. Since the embedding points of the shared dataset did not change, all local embeddings are easily combined by aligning the points representing the shared data.

2.3 Multi-shot dSNE

Algorithm 4 GradStep

Input:

Data embeddings: Yp (local), Ys (shared), P Optimization parameters: η, α

Output: ˆYp (local), ˆYs (shared)

• 1: Eq. (2) to compute low-dimensional affinities qij
• 2: Eq. (3) to compute δC / δyi
• 3: ˆyi = η(δC/δyi) +α(yt-1 i -yt-2 i )
• 4: group ˆyi into ˆYp (local) and ˆYs (shared)

Singeshot dSNE, however, tends to produce results with significant overlap as the data from different sites does not affect each other and often end up in the same location despite belonging to different classes. To overcome this problem we developed an iterative algorithm: multishot d-SNE. Unlike the single shot version we allow the embedding of the shared dataset to change. However, we change it exactly the same way for each site. This is achieved by averaging the updates to shared data embedding based on the updates that each of the local sites requires to be made. This allows information about desirable embedding of the local data points to flow across sites and thus communicate preferred location through the influence on the shared data. Algorithm 6 details all of the steps.

Algorithm 5 UpdateStep

Input:

Data embeddings: Yp, Ys, ˆYp, ˆYs

Output: Yp, Ys

• 1: Yp =Yp + ˆYp
• 2: Ys =Ys + ˆYs

2.4 ComparisonMetrics

An objective comparison of inherently subjective visualization algorithms is difficult. We quickly found that the k-

Algorithm 6 multishotDSNE

Input:

Objective parameters: ρ (perplexity) Optimization parameters: T , η, α Shared Data: Xs = [xs1 , xs2 . . . xsNs ], xsi ∈ Rn Data at site p∀p: Xp = [xp1 , xp2 . . . xpNp ], xpi ∈ Rn Output: Y = {y1, y2, . . . , yN}, yi ∈ Rm,m << n,N = ∑p Np + Ns 1: Ys ∝ N(0, 10-4 I), I ∈ Rm×m initialize from Gaussian 2: for p = 0 to P do ▷ Initialize at sites 3: Y→p s 4: Pp ←PairwiseAffinities p, ρ, [Xp,Xs] 5: Yp ∝N(0, 10-4 I), I ∈ Rm×m 6: end for 7: for i = 0 to T do 8: for p = 0 to P do ▷ At local sites 9: ˆYp, ˆYs ←GradStep [Yp,Ys,Pp] 10: end for 11: ˆY ← 0 ▷ At the master 12: for p = 0 to P do ▷ At the master 13: ˆY←p s 14: ˆY← 1P ˆYs ▷ Average local ˆYs 15: end for 16: for p = 0 to P do ▷ At local sites 17: ˆY→p 18: Yp,Ys ←UpdateStep [Yp,Ys, ˆYp, ˆY] 19: end for 20: end for

means criterion is only weakly correlated with the usefulness of produced embeddings:

In an attempt to quantify perceptional quality of the resulting embeddings we have developed two additional metrics: overlap and roundness. To remove sensitivity to noise we first remove the outliers [see Figure 1] in each digit's cluster (Liu et al., 2012). Then compute the convex hull for each digit and use them to compute the measure of the overlap (the sum of all polytope areas minus the area of the union of all the polytopes normalized by this union's area) and the roundness (ratio of the area of each polytope to the area of the circumscribed circle).

3 Results

We base our experiments in this section on two datasets: MNIST (LeCun et al., 1998) for handwritten images of all digits in the 0 to 9 range, and Autism Brain Imaging Data Exchange (ABIDE) for fMRI data (Di Martino et al., 2014). MNIST data were taken from a Kaggle competition 3 which has 28, 000 gray-scale images of handwritten digits. Each image contains 28 × 28 = 784 pixels. Among these data, we randomly (but preserving class balance) pick 5,000 different samples from the data set for our needs. At first, we reduce dimension of the data from 784 to 50 using PCA. Then, dSNE and tSNE are used to generate (x, y) coordinates in two dimensional space. The second dataset is the ABIDE fMRI set, which contains data of 1153 subjects. The ABIDE data has been pre-processed down to multiple spatial and temporal quality control (QC) measures 4 . For ABIDE, because of the low dimension of QC measures, we do not use dimensionality reduction but directly run our tSNE and dSNE to produce the embeddings.

3 https://www.kaggle.com/c/digit-recognizer

3.1 MNIST Experiments

Experiment 1 (No diversity in the reference data): In this experiment, the reference dataset contains a single digit that is also present at all of the three (3) local sites. We run an experiment for each of the 10 digits. Each site contains 400 samples for each of its corresponding digits. Reference dataset contains 100 samples of its digit.

Experiment 2 (Effect of the sample size): Often centralized stores accumulate more data than any separate local site can contain. Our goal is to check the adaptability of our algorithm on this case. In this experiment, every local site contains only one digit and the reference dataset contains all digits (0-9). We consider 2 cases: when each site contains 400 samples and each digit in the reference consists of 100 samples; and the inverse case, when sites only have 100 samples, while the reference digits are represented by 400.

Experiment 3 (Missing digit in reference data): In this experiment we investigate the effect of the case when a digit is missing from shared data. This approximates the case of unique conditions at a local site. Each local site out of 10 contains a single digit. We run 10 experiments; in each, the reference dataset is missing a digit. For each of the experiments we have 2 conditions: in one the reference dataset is small (100 samples for each digit but the missing one) while the sites are large (400 samples per site); in another the reference data is large (400 samples per digit) and the sites are small (only 100 samples).

Experiment 4 (Effect of the number of sites): In this experiment, we investigate whether the overall size of the virtual dataset affects the result. Every local site, as well as the reference data, contains all digits (0-9) 20 samples per digit. We continuously increase the number of sites from 3 to 10. As a result, the total number of samples across all sites increases.

3.2 Single-shot

Figure 2 shows examples of the output generated by the single shot dSNE algorithm together with tSNE layout of the centralized data for comparison. We run single shot for the datasets of Experiment 2, Experiment 3, and Experiment 4. For all experiments, we are able to correctly group and embed same digits from different sites. However, the problem was that the digit clusters tended to heavily overlap. Even for experiment 4, where we found the best results, one can still observe heavy overlap of digit clusters.

4 The full list is available here: https://github. com/preprocessed-connectomes-project/ quality-assessment-protocol/tree/master/ normative_data

3.3 Multi-shot

The multi-shot algorithm was introduced to cope with the overlap problem. Figure 3 shows its performance in experiment 1. In this figure, the plots on the left show performance metrics in comparison to those of tSNE on centralized data. For each digit we rerun the experiment 10 times with a different seed value. On the right, the top figure presents the best, and the lower figure represents the worst performing run. The layouts are the same but colored by digits and by sites. From the analysis of the results, we find that each common digit from different sites is embedded in the plots perfectly. But, for many digits, the clusters are less separable.

Figure 4 represents the result of experiment 2. The comparison metrics show that when the shared portion contains large amount of data the metrics are better than in the case of smaller number of samples in the reference dataset. However, the cluster roundness degrades with the size of the sample in the shared data. dSNE clusters are less round compared to the centralized tSNE.

Figure 5 depicts the results of Experiment 3. Although there is some variability with respect to which digit is missing from the reference dataset, this is not the largest effect. Here, we observe as similar behavior as in experiment 2: larger size of the reference dataset leads to better results. Note, when the reference dataset contains a large number of samples, we always get better results for Kmeans and Intersection Ratio compared to the case of fewer samples in the reference. Nevertheless, a visual inspection of the best results in the case of smaller reference dataset size does not lead to conclusions that the results are unusable.

Figure 6 depicts the results of Experiment 4. The effect of sample size is the most pronounced for tSNE output, while the metrics remain fairly constant on the multishot dSNE embeddings regardless of the size of the virtual dataset. This is a surprising outcome. We present best and worst embedding based on the number of different local sites in our decentralized scenario.

3.4 Real Data

We investigate performance of multishot dSNE in comparison with the embedding produced by tSNE on the pulled data using the QC metrics of the ABIDE dataset. To simulate a consortium of multiple sites we randomly split these data into ten local and one reference datsets. Results show 10 different clusters for centralized data. For three random splits of our decentralized simulation we obtained 10 different clusters as well (see Figure 7). Notably, the split into the clusters in the embedding is stable regardless of the split into sites.

4 Conclusions

Our approach enables embedding and visualization of high dimensional private data spread around multiple sites. The data does not leave the sites and only minimal gradient information from the embedding space gets transferred across the sites. We consider this approach plausibly private as most of the information about individual samples was discarded. Extensive tests on MNIST demonstrate the usefulness of the approach and high quality of the obtained embeddings over a variety of settings. Although multi-shot dSNE is quite robust to various conditions, such as changes in the number of sites, rare or missing data etc., the best performance is achieved when the reference dataset is dense. Notably, the single-shot dSNE, which is mostly an implementation for tSNE of the previously existed landmark point method, tends to ignore the differences across the sites and only respects the reference dataset. An alternative solution-an average of the gradients weighted by the quality of their respective local tSNE-may even further bias results toward good local groupings but unacceptable overall embedding. In contrast, our multi-shot approach provides enough information propagation to warrant better embeddings. Yet, all that is being exchanged is pertinent to the already public reference data. We conclude that dSNE is a valuable quality control tool for virtual consortia working with private data in decentralized analysis setups.

WORST