# Obfuscation of Discrete data

Saswata Naha111
Sayantan Roy222
Arkaprava Sanki333
Diptanil Santra444Indian Statistical Institute Kolkata

(June 2021)

###### Abstract

Data obfuscation deals with the problem of masking a data-set in such a way
that the utility of the data is maximized while minimizing the risk of the
disclosure of sensitive information. To protect data we address some ways that
may as well retain its statistical uses to some extent. One such way is to
mask a data with additive noise and revert to certain desired parameters of
the original distribution from the knowledge of the noise distribution and
masked data. In this project, we discuss the estimation of any desired
quantile and range of a quantitative data set masked with additive noise.

## Introduction

Privacy protection and data security have recently received a substantial
amount of attention due to the increasing need to protect various sensitive
information like economic data and medical data. In case of continuous data, a
very famous example is publication of marks of students. An Institute would
like to publish the performance of its students without disclosing the actual
marks. Whereas, for discrete case, data on economic class study can be used,
where it is of national interest to know about the economic status of the
country, not disclosing the sensitive information of income of any individual.
We are interested in the obfuscation of discrete case only and the desired
estimations from such a data.

## Basic Problem

Our main goal through this project is to construct a way of obfuscation of
discrete data using some additive noise in such a way that the probability of
estimating the actual data given the obfuscated data is quite low. However,
given the nature of the noise added to the original data, we can estimate
certain population statistics with low error.

Suppose, $X$ be the data vector containing the individual values of the
population. We want to publish the distribution of $X$,
$p_{i}=\mathbf{P}(X=x_{i})\>\forall\>x_{i}\in$ range of $X$. For additive
model, we add noise vector $Y$ with
$q_{j}=\mathbf{P}(Y=y_{j})\>\forall\>y_{j}\in$ range of $Y$ depending on the
range of $X$ and get the masked data $Z=X+Y$. Our aim is to adjust $Z$
accordingly to publish the obfuscated data $Z^{\prime}$ such that
$\mathbf{P}(X\mid Z^{\prime})\approx 0$.

As for the next part of the project we want to estimate certain statistical
properties of the true data viz. Median or other quantiles and Range. For
this, we would like to estimate the actual probabilities $p_{0},p_{1},\cdots$.

Note that,
$\displaystyle\mathbf{P}(Z=z_{k})=\displaystyle\sum_{i=0}^{k}\mathbf{P}(X=x_{i})\mathbf{P}(Y=z_{k}-x_{i})=\sum_{i=0}^{k}p_{i}q_{j}$.
Estimating $\mathbf{P}(Z=z_{k})$ from the obfuscated data, we can try to solve
for $p_{i}$’s from the given data of $q_{j}$’s.

## Dataset

We have generated some poisson datasets and added some noise to them. Then
truncated that in such a way that the probability of estimating the actual
data from the obfuscated data is quite low.

At first, 10 lakhs data from poisson with parameter generated from exponential
(2.5).

As we have mentioned before that economic data are more sensitive in nature
and instead of publishing the actual data we want to obfuscate it before
publishing it. This kind of data usually take positive values and they are
positively skewed. We have generated the data in such a way that it can be
compared with an economic data. Here in the data set, the variable in the
X-axis denotes different economic class. The whole population is divided into
12 economic classes where $0^{th}$ class denote the lowest income group and
$12^{th}$ class denote the highest income group. The variable in Y-axis
denotes the frequency of each class. Then we obfuscate it in the following
way:

We add an additive noise from unif(0,10).

Then we collect the household size data set from the website of “Census of
India” and have used the data for the whole Indian region.

The data is shown in the diagram below -

Here X-axis denotes the household size and Y-axis denotes the frequency
corresponding to each household size.

Then we similarly add an additive noise from unif(1,10).

## Goodness of Obfuscation

In this section, we want to show that using this obfuscation procedure,
probability of getting back the actual data from the published obfuscated data
is quite low.

Let X,Y,Z be the values of the actual data, added noise and values of
obfuscated data respectively. Now, we have to calculate $P(X|Z)$. We know,

$\mathbf{P}(X=x|Z=z)=\frac{\mathbf{P}(X=x,Z=z)}{\mathbf{P}(Z=z)}$

.

To calculate the probabilities in the numerator and denominator separately,
first we construct a matrix with the values of actual data (X) in the rows and
the corresponding error values (Y) in the columns. So the $(i,j)^{th}$ element
of the matrix represents number of individuals in class i has taken error
values j. So the size of the matrix is $12\times 10$. For our calculation, as
the values of X have been fixed as x, the row wise sum of the matrix is fixed
for all rows. Also note that, the sum of the off-diagonals of the matrix
represents the values of Z. Fixing these values according to our data, we get
the total number of possible cases, from which we can calculate the exact
probability.

Next, for the probability in the denominator we take all possible cases for X
and Y such that the data in our case is valid in real situation. The
computations are very cumbersome∗, however, the actual conditional probability
turns out to be of the order $\mathcal{O}(10^{-10^{6}})$. So we can conclude
that, since this probability in itself is very small, if we further mask the
obfuscated data, the probability of getting the true data conditioned on that
will be even smaller. So this method of obfuscation can be said to be
satisfactory.

[* Cumbersome in the sense that, we have made an algorithm to count the number
of matrices which satisfy the required conditions. In this way we got the
proportion of matrices satisfying the conditions which will be close to the
required probability.]

## Different approaches to estimate quantiles

Quantiles are robust statistics and hence is not much affected by the presence
of outlying data-points. To do the estimation of quantiles we first estimate
the c.d.f. of X from the obfuscated data using different approaches.

### First Approach

First we try the most common method to estimate quantiles that is, using least
square method from the complete obfuscated data. Moreover, we know that the
least square method gives an unbiased estimate of parameter. So we shall use
this in the following way:

We know that for a least square model,

$\tilde{Y}=X\tilde{\beta}+\tilde{\epsilon},$ where,
$\tilde{\epsilon}\sim\mathcal{N}(\tilde{0},\sigma^{2}I)$,

We get the solution as $\hat{\beta}=(X^{\prime}X)^{-1}X^{\prime}\tilde{Y}$.

However, as the components of $\tilde{\beta}$ are probabilities and add up to
1, implementing the constraint of $\tilde{\beta}^{\prime}\tilde{1}=1$, we get
the solution as

$\hat{\beta}=(X^{\prime}X)^{-1}(X^{\prime}\tilde{Y}-\frac{\tilde{1}^{\prime}(X^{\prime}X)^{-1}X^{\prime}\tilde{Y}-1}{\tilde{1}^{\prime}(X^{\prime}X)^{-1}\tilde{1}}\tilde{1})$.

The result for our data set is like the following -

[Necessary R codes are given in Appendix section.]

But in this approach,as can be seen in the estimated solution, some of the
components come out to be negative and hence cannot be probabilities. So, we
need a method to further implement the non-negativity criteria of the solution
or need to change the published data in some suitable way.

First we try to modify the solution method keeping the whole obfuscated data.
We suspect that the inverse function in numerical methods can affect the
solution to a certain effect. So, we then try using QR decomposition to check
if there is any improvement. But some of the probabilities are found to be
still negative. So, we conclude that, we need to rectify the approach itself.

### Second Approach

Since the first approach failed, we shall try MLE though it is not unbiased.
But we know that MLE maximizes the likelihood (i.e, probability of the value
of the parameter is the maximum according to the sample) and consistent. Also
MLE belongs to the parameter space. Therefore, our next attempt is using MLE
with truncation. Suppose, $\tilde{p}$ and $\tilde{r}$ be the probability
distributions of the actual data and the obfuscated data respectively.
Maximizing the likelihood w.r.t $\tilde{r}$ we get the MLE estimates as
$\displaystyle r_{i}=\frac{n_{i}}{\Sigma n_{j}}\>\forall\>i$.

So, solving w.r.t $\tilde{p}$ we get-

$\displaystyle p_{0}=\frac{n_{0}}{\sum n_{i}}*11$, $\displaystyle
p_{j}=\frac{n_{j}-n_{j-1}}{\sum n_{i}}*11$ for j = 1,2,…,11 and $\displaystyle
p_{12}=1-\sum p_{j}$.

Similarly as the above, we can also estimate $\tilde{p}$ as -

$\displaystyle p_{12}=\frac{n_{22}}{\sum n_{i}}*11$, $\displaystyle
p_{j}=\frac{n_{j+11}-n_{j+10}}{\sum n_{i}}*11$ for j = 1,2,…,11 and
$\displaystyle p_{0}=1-\sum p_{j}$

Using R, both of the estimates are calculated, the results are shown below -

[Necessary R codes are given in Appendix section.]

However, in this case as well, some of the estimates turn negative and one way
to deal in such situation is combine the above to estimates to get a better
estimate.

A general technique in this combination is depending on data. Take that part
of first estimate where the height of the columns of obfuscated data are
increasing and that part of second estimate where the height of columns are
decreasing.

For our data set, the combined estimate is as following -

[Necessary R codes are given in Appendix section.]

Now, there is two issues regarding this method - First of all, this method is
on non-truncated obfuscated data, where there is a data leak in the extreme
values. Secondly, the technique of the combination as mentioned above, failed
to work if there are ups and downs in frequency all over the range.

To avoid the first issue, we shall work on the truncated obfuscated data and
for the second issue, we have to somehow extend our MLE approach to
Constrained MLE.

### Third Approach

For the constrained MLE we need to include the constraints $p_{i}\geq 0\hskip
2.84526pt\forall\hskip 2.84526pti$ i.e $r_{1}\leq r_{2}\leq r_{3}\leq...\leq
r_{n}$. For that we are transforming the unconstrained vector
($u_{1},u_{2},...,u_{n}$) to constrained vector ($r_{1},r_{2},...,r_{n}$) by
using nested logistic transformation -

$\displaystyle r_{k}=T_{k}(u)=\frac{1}{1+\sum_{i=k}^{n}e^{-u_{i}}}\hskip
2.84526pt\forall\hskip 2.84526ptk=1,2,...,n$

Now by replacing ($r_{1},r_{2},...,r_{n}$) in our likelihood function and
differentiating it w.r.t. each $u_{i}$, we get -

$\displaystyle\frac{\delta l}{\delta
u_{k}}=-n_{0}\frac{\sum_{j=1}^{k}\frac{e^{-u_{k}}}{(1+\sum_{i=j}^{n}e^{-u_{i}})^{2}}}{1-\sum_{j=1}^{n}\frac{1}{1+\sum_{i=j}^{n}e^{-u_{i}}}}+\sum_{j=1}^{k}\frac{n_{j}.e^{-u_{k}}}{1+\sum_{i=j}^{n}e^{-u_{i}}}=0\hskip
2.84526pt\forall\hskip 2.84526ptk=1,2,...,n$

To solve these equations we tried Sequential Quadratic Programming method for
which we need to find the gradient and hessian matrix by differentiating the
system of equations w.r.t $u_{i}$’s once and twice respectively. But it was
very difficult to calculate. So we have tried another method.

### Fourth approach

We have seen that obtaining the probabilities analytically was difficult and
problematic. So we shall now try numerical method. In this method also we
shall try to use constrained MLE by calculating the $p_{i}$’s in iterative
procedure. The steps are given below:

Step 1: Define the likelihood function as per the obfuscation procedure.

Step 2: Initialize $p_{i}$’s such that $p_{i}\geq 0\hskip
2.84526pt\forall\hskip 2.84526pti=1,2,...,n$ and $\sum p_{i}=1$

Step 3: To determine the values of $p_{i}$’s we want the precision of the
values of $p_{i}$’s upto third digit (say). In each step of the iteration we
fix the sum and so we update $p_{i}$ by $sum\times\frac{1}{1000}$ and
$p_{i+1}$ by $sum\times\frac{999}{1000}$ and check whether the value of the
likelihood function is greater or not. In this way we shall update $p_{i}$ and
$p_{i+1}$ by $sum\times\frac{j}{1000}$ and $p_{i+1}$ by
$sum\times\frac{(1000-j)}{1000}\hskip 2.84526pt\forall\hskip
2.84526ptj=0,1,2,...,1000$ and compare with the likelihood value. In this way
we shall get the values of $p_{i}$ and $p_{i+1}$ for which likelihood function
will be maximum. We do the same thing for different values of i and repeat the
whole process k times. Where k is sufficiently large.

Correctness of the above process: In the process we are keeping the sum fixed
and vary the ‘j’ to see for which proportion of ‘sum’ should be $p_{i}$. So
the sum of the $p_{i}$’s are always be 1 and they all are positive. [Thus we
avoided the problems occurred in the previous methods.] In each step we are
comparing the value of likelihood function in current step with that of the
previous step. And as we are taking k to be sufficiently large (i.e we are
repeating the iterative process sufficient no. of times), the process will
return the correct constrained MLE values.

### Implementation

1\. We first implement this procedure to the generated data set which we have
used in previous approaches.

Here the log likelihood function is $L=n[1]\times log(x[1])+n[2]\times
log(x[1]+x[2])+\cdots+n[11]\times log(x[1]+x[2]+x[3]+\cdots+x[11])+n[12]\times
log(x[2]+x[3]+\cdots+x[12])+\cdots+n[23]\times log(x[13])$

where, $x[i]$ be the probaility of $i^{th}$ class and $n[i]$ be the observed
count in that class.

[Necessary R codes are given in Appendix section.]

The results are as following -

The vector ‘p’ gives the estimate of $p_{i}$’s using the iterative method and
‘table(u)/sum(table(u))’ shows the proportions of values from actual data. We
can observe that the estimation is quite well i.e the values in the two
vectors are close enough. That’s why quantiles can be estimated with high
accuracy by taking the cumulative sums of the vector p.

2\. To check the consistency we generated a data set in the same way as the
previous and implemented the procedure for that data set.

The results are shown below:

Here also the vector ‘p’ gives the estimate of $p_{i}$’s and
‘table(u)/sum(table(u))’ shows the proportions of values from actual data. We
can observe that in first few components the estimation is quite well but in
last part the difference between the estimated values and the actual values
are not that small and also there is a 0(zero) in $15^{th}$ component. To
describe this observation we go through the diagram of obfuscated data.

We can notice that there is no observation for the value 24 which is affecting
the likelihood function and consequently giving a poor estimate near to
corresponding value in actual data.

So for such cases we suggest to merge the class with any of the neighbour
class and modify the likelihood function accordingly.

3\. Now we shall try to implement the procedure for the household size data
set.

Before implementing, we need to construct the raw data set from the summarize
data. The main issue for construction is that there are some classes in the
data set that contain a class of values and also the last class is unbounded.
So we need to think about some techniques about breaking the classes in
discrete points in order to construct the raw data set that we need.

We have considered poisson distribution with parameter average household size
(i.e 4.8) as the best fit for the data. We break those classes into discrete
values according to the proportion of expected frequencies corresponding to
each value.

After breaking the data set is like following:

Here, the likelihood function is -

$L=n[1]\times log(x[1])+n[2]\times log(x[1]+x[2])+\cdots+n[10]\times
log(x[1]+x[2]+\cdots+x[10])+n[11]\times
log(x[2]+x[3]+\cdots+x[11])+\cdots+n[22]\times
log(x[13]+x[14]+\cdots+x[22])+\cdots+n[29]\times log(x[20]+x[21]+x[22])$

where, $x[i]$ be the probaility of $i^{th}$ class and $n[i]$ be the observed
count in that class.

The results are as following:

Similarly as the generated data for this household size data set also this
procedure gives a quite good estimate of the distribution as we can see in the
result.

Here we are using the obfuscated data set for our estimation but to maintain
the data privacy at extreme values we need to do truncation. So our next aim
is to implement our procedure for the truncated data set.

We shall use our household size data set and will truncate the obfuscated data
at 22. The truncated obfuscated data set will be as following:

In this case the log likelihood function is -$L=n[1]\times
log(x[1])+n[2]\times log(x[1]+x[2])+\cdots+n[10]\times
log(x[1]+x[2]+\cdots+x[10])+n[11]\times
log(x[2]+x[3]+\cdots+x[11])+\cdots+n[21]\times
log(x[12]+x[14]+\cdots+x[21])+n[22]\times
log(10-10\times(x[1]+\cdots+x[12])-9\times x[13]-\cdots-2\times x[20]-x[21])$

where, $x[i]$ be the probaility of $i^{th}$ class and $n[i]$ be the observed
count in that class.

The result are shown below -

As we can see from the result, the first few estimated probabilities are very
close to that of actual data but the estimation gets worse at the end. This is
because the fact that we truncated the observations at 22 and hence
considering all the observations that have value more than or equal 22 as same
quantity.

Although if we focus on estimation of quantiles, we proceed to check the
cumulative sum of each vector and the result are as following -

As we can observe upto 99th percentile, the estimation is quite accurate. Only
the last quantile can’t be predicted due to the same reason of truncation.

Therefore, our procedure works fine if the data need to be published after
truncation.

Point to be noted here, in the whole procedure we described, needed the
likelihood function to be perfectly defined. Now to define the likelihood
function, it is necessary to have a bit of knowledge of the range and if the
range is published again we have compromise some data privacy. To maintain the
privacy, we suggest to publish a wider range and we have checked that our
procedure can be implemented with a knowledge of wider range also.

Implementation: Here we consider the without truncation case where the actual
range was (1,22), instead here we are publishing a wider range (1,25).

Therefore, for this case, log likelihood function is -

$L=n[1]\times log(x[1])+n[2]\times log(x[1]+x[2])+\cdots+n[10]\times
log(x[1]+x[2]+\cdots+x[10])+n[11]\times
log(x[2]+x[3]+\cdots+x[11])+\cdots+n[22]\times
log(x[13]+x[14]+\cdots+x[22])+\cdots+n[29]\times
log(x[20]+x[21]+x[22]+\cdots+x[25])$

where, $x[i]$ be the probaility of $i^{th}$ class and $n[i]$ be the observed
count in that class.

So we can see from the result that even if we declare a wider range, the
estimated probability are close to the actual value and the probabilities
corresponding to values 23-25 which are not in the actual data coming to be
zero as expected.

## Goodness of estimates

To test for the goodness of the estimates we compute the individual variance
and MSEs of each component of the estimate. Bootstrapping from the predicted
data we numerically generate multiple copies of the raw obfuscated data and
numerically calculate the variances and MSEs. The computed vectors are as seen
below:

We note that the variances are of the order $\mathcal{O}(10^{-3})$ and the
corresponding MSEs are very close to them. Hence, we can say that the
estimates are of low bias and variance and so it is consistent.

## Approaches for estimating Range

Since our actual goal is to hide the information about the raw data. So we
want to see whether the estimate of range of the actual data can be estimated
sufficiently close to the real value.

Using Law of large number555Nandi Mridul $\&$ Roy Bimal [2021] Estimation of
extremum of obfuscated data and its error of estimation:

We first try to estimate the range from the data set that we generated from
poisson distribution and added an additive noise Discrete Uniform (0,10) to
it. Suppose‘W’ is the obfuscated data. We shall add the same noise 999 times
to the obfuscated data. So we have ‘W1’(say) which is a data after obfuscation
the actual data 1000 times. If we consider 1000 as a large number, we know
that the maximum element of actual data should be close to (the maximum
element of ‘W1’-5*1000)[Using Law of large number]. But we know that the
maximum should be 12 according to our data and we observe the results to be in
between 450-500. The results in R are given below:

We have also observed that if we increase the number of noise, the resulted
maximum value is more far from that we expect.

### Another approach

We know the minimum household size is 1. So here we are trying to estimate the
maximum of the household size from the obfuscated data.

Method of finding maximum:

Step 1: Estimate the quantiles from obfuscated data by previously mentioned
method (Fourth approach to estimate quantiles).

Step 2: The maximum is estimated as the value at $100^{th}$ percentile. It may
happen that the $100^{th}$ percentile is not unique. In that case, the maximum
will be estimated from the minimum value that is at $100^{th}$ percentile.

Implementation:

Case 1: (Without Truncation)

By the procedure described above we check the $100^{th}$ percentile using our
estimated probabilities and we can see that $100^{th}$ percentile occurs first
at value 21. [The maximum value of the actual data is 22. So there is a small
error in estimation.]

Case 2: (Truncated)

We can see that in this case $100^{th}$ percentile occurs first at value 15,
which is far from that of the actual value of maximum. This result was
expected because we have truncated the obfuscated data at 22, that’s why the
values close to the end are given less priority. So, to estimate the maximum
value from the truncated obfuscated data is not as easy as the previous case.

## Conclusion

1. The method that we have used to obfuscate our data i.e. the additive white noise model is a good way of masking a discrete data set i.e. probability of estimating the actual values of the data (e.g. household size) from the obfuscated data is very low.

2. Given real life data set related to household size, we have shown the actual probability with estimated probability for different household size in the following plot:

We can notice in the above plot that the actual probability and the estimated
values of the probabilities are very close to each other (the plots almost
coincide with each other.

3. To estimate the probabilities from obfuscated data we tried Least square, MLE and constrained MLE methods. But there were some problems (e.g. estimate outside parameter space) which were solved in the fourth approach (numerical method). Also we have seen that the result obtained using iteration in this approach was very good.

4. While estimating the range of the raw data set, we first tried the method using LLN. But in that case the estimated values were very far from the actual value. Then we have tried another approach to estimate maximum value of the raw data set. In this method, we see that upto the $99^{th}$ percentile can be estimated quite accurately, however, the estimate for the $100^{th}$ percentile varies quite a lot from its real value. So there is further scope of study in this regard.

## Appendix:

R codes for First Approach :

R codes for the second approach :

R codes for Fourth approach :

Code for the generated data :

Code for the non-truncated household size data :

Code for the truncated household size data :

Code for the non-truncated household size data with wider range given:

## Acknowledgements:

We give our thanks to Prof. Bimal K. Roy, Prof. Nachiketa Chottopadhyay and
Prof. Prabal Choudhury for their kind advice and support in this project. This
project was also inspired by some previous works done by Prof. Bimal K. Roy et
al.

## References :

  1. 1.

Census of India, https://censusindia.gov.in/2011census/hh-series/hh01.html

  2. 2.

https://math.stackexchange.com/questions/2764174/how-to-solve-maximum-
likelihood-estimates-with-inequality-constraint

  3. 3.

Ghatak, Debolina $\&$ Roy, Bimal. (2018). Estimation of True Quantiles from
Quantitative Data Obfuscated with Additive Noise. Journal of Official
Statistics. 34. 671-694. 10.2478/jos-2018-0032.

  4. 4.

Nandi, Mridul $\&$ Roy, Bimal. (2021). Estimation of extremum of obfuscated
data and its error of estimation.


Tools and Methods This section contains background information about the
various components of the robotics simulation pipeline. After reading this
section, one should have enough context to understand how these components fit
into the workings of the pipeline.

Overview Webots.HPC is a pipeline, and thus owes existence to many popular
tools, software, and methods. Table tab:materials shows all of the major
materials involved in the Webots.HPC pipeline development. The following
subsections describe the involved project materials more at length.

Docker | Containerization
---|---
Singularity | Python
pip | External Python Libraries
SSH | X11/Xvfb
Webots | SUMO
Palmetto Cluster | DICE Lab Queue
Portable Batch System | Distributed Computing
Parallelization | Big Data
Table : Tools, Software, and Methods

##  1 Docker, Singularity, and Containerization

###  Containerization

At a broad level, Docker and Singularity are both a packaged set of platform
as a service (PaaS) products that implement containerization
[noauthor_docker_2021]. Containerization is the packaging of software code
with just the operating system (OS) libraries and dependencies that are
essential to run the code; the result is a singular, low-footprint executable
(called a container) that runs invariably on any infrastructure
[ibm_containerization_2021]. Containers have overtaken virtual machines in
recent years due to their portable and resource-efficient nature.

###  Docker

As mentioned, Docker is a containerization tool. Docker is the platform that
revolutionized the cloud computing field in the early 2010s, being a unified
tool that allows users to easily package an application and all required
dependencies into a container [scheepers_virtualization_2014]. Docker is
undoubtedly the de-facto standard for building containerized apps, having
acquired millions of users since its inception in March 2013
[noauthor_docker_2021]. Lastly, it is important to note that Docker is more
popular with commercial and professional users; alternative platforms are
sometimes more apt for academic usage.

###  Singularity

Singularity is the containerization software designed explicitly for High-
Performance Computing (HPC) resources (like Clemson’s Palmetto Cluster) and
scientific computing. Between Docker and Singularity, the end result is the
same; a portable, reproducible container that can run on any environment
similar to the environment the container was originally created in. A major
difference between the two platforms is that Singularity’s containers can be
created and deployed on HPC clusters, which was a previously unmet need prior
to its inception in 2015 [kurtzer_singularity_2017]. Another key difference
between the two platforms is that Singularity is much more secure, offering
features like no privileged or trusted daemons, access to host filesystem
without safety implications, a stricter privilege model, and the ability for
admins to control and limit user capabilities [kurtzer_singularity_2017].

##  2 SSH, X11, and Xvfb

###  SSH

Secure Shell (SSH) is a protocol for remote connections between computers
[ylonen_sshsecure_1996]. SSH is currently the most popular method for inter-
computer connection, as it offers several security features like encrypted
connections and password or public key-based authentication
[ylonen_sshsecure_1996]. As Webots.HPC is a resource for HPC environments,
users typically access the pipeline by SSH-ing onto said HPC environment from
their personal computer.

###  X11 and Xvfb

X is a communication protocol that provides an abstract interface for managing
a GUI. When SSH-ing into a remote computer with the intent to work with GUI-
enabled applications, X11-forwarding is required for streaming the screen
renderings from the remote computer [nye_x_1995]. When running GUI-enabled
applications in headless mode, you need to configure an Xvfb, or X virtual
framebuffer, to send renderings to [wolfe_virtual_1990]. In contrast to the X
Window 11 (X11) system, Xvfb performs all graphical operations in virtual
memory without showing any screen output [wolfe_virtual_1990].

##  3 Python, Pip, and Dependent Libraries

###  Python

Python is one of the most recognizable programming languages worldwide, and
for a good reason. While it is a slower, interpreted language, its extreme
popularity is derived from the readability of Python code, with simplified
syntax made with natural language in mind. [van_rossum_python_2007] Python was
factored into this pipeline from foresight, as we knew that the Webots-SUMO
simulation we planned to run on the pipeline was coded almost entirely in
Python. However, it is important to note that Webots also offers support for
C, C++, Java, and MATLAB with identical results.

###  Pip Package Installer

Pip is the package installer bundled with Python that allows users to install
packages from the Python Package Index and other indexes
[ansari_prerequisites_2020]. Throughout pipeline development, it was necessary
to install pip onto our singularity image, as to enable the ability to install
other prerequisite libraries.

##  4 Simulation Software

The pipeline leverages Webots as its foundational simulation software. If a
simulation includes traffic simulation, the pipeline also leverages SUMO for
additional support.

###  Webots

Webots is a cutting-edge open-source robotics simulator used for robotics and
transportation research. Webots offers native integration with SUMO and allows
us to model an actual CAV, simulating a range of sensors found on real-world
prototypes. As Webots is the primary simulation software that our pipeline
leverages, it is useful to review the way that Webots handles typical
simulation tasks. Webots stores scenes with a tree structure, where the root
node is the world and children nodes are the different items (robots,
aesthetic scenery, and sensors, etc.) in the world. Robot nodes should always
be under the root node, and also are assigned a controller. Controllers are
scripts (in C, C++, Java, MATLAB, or Python) that determine a node’s
functionality. Controllers are also the main interface through which Robot and
Sensor nodes communicate. Robot and Sensor parameters are tweaked through the
Webots User Interface; however, for sensors, the sampling period (a very
important parameter that influences simulation accuracy and performance) is
specified in the controller for an attached robot [cyberbotics_webots_2021].

###  Simulation of Urban Mobility (SUMO)

Simulation of Urban MObility (SUMO) is an open-source traffic simulation
package with net import and demand modeling components [behrisch_sumo_2011].
Since its initial development by the German Aerospace Center (DLR) in 2001,
SUMO has grown into a full-featured suite of traffic modeling tools, including
a road network platform capable of reading different input formats, demand
generation and routing utilities from various sources (origin-destination
matrices, traffic counts, etc.), a high-performance simulation usable for
situations from single junctions to entire cities, and a ”remote control”
interface called TraCI (Traffic Control Interface) to adapt a simulation
online.

###  Pairing Webots and SUMO

Independently, SUMO is capable of simulating many different types of travel,
including trains, boats, bicycles, and pedestrians, in supplement to basic
vehicular travel. However, when pairing SUMO with Webots, only automotive
travel is easy to simulate; this is due to a lack of pre-existing research
addressing alternative transportation methods. A SUMO vehicular traffic
simulation can easily be paired with Webots through the bundled SUMO
Interface, which is implemented as a child node in the simulation. Most of the
tunable SUMO parameters must be manually edited through the simulation’s
traffic configuration files (sumo.net.xml, sumo.flow.xml, sumo.edg.xml, etc.);
however, opposite of sensors, the sampling period of the SUMO Interface must
be specified in the Webots user interface. When Webots is paired with SUMO,
Webots takes on the role of the user-visible simulation front-end, offering
complex and customizable camera angles and highly detailed vehicle models.
SUMO, on the other hand, acts as the simulation back-end, sending TraCI
signals through the Webots SUMO Interface to directly control any involved
traffic. As a general analogy, think of Webots as a puppet and SUMO as the
puppeteer. However, Webots does also offer augmented functionality through
sensors. SUMO simulations can provide extensive output, but Webots sensors can
augment said output substantially. Radars, cameras, compasses, distance
sensors, light sensors, and touch sensors can all be added to Webots vehicles
to better understand the state of both the central vehicle and surrounding
traffic.

##  5 Palmetto Cluster

The Palmetto Cluster is a large-scale distribution testbed for scientific
research applications consisting of 1,700 machines maintained by the Clemson
University Cyberinfrastructure Technology Integration (CITI) group
[du_distributed_2018]. Palmetto is Clemson University’s primary high-
performance computing (HPC) resource. Through Palmetto, we are able to run
multiple simulations at the same time to collect data in parallel. Within
Palmetto, we utilized the Clemson University Data Intensive Computing
Environments (DICE) Lab queue, which contains several high-powered, big
memory-oriented computers. Table 5 contains the hardware specifications of the
11 nodes that the DICE Lab has purchased.

Make | Dell
---|---
Model | R740
Chip | Intel Xeon
Cores | 40
RAM | 744 GB
Local Scratch | 1.8 TB
Interconnect | 100g, HDR, 25GE
GPUs | 2
GPU Make | Nvidia
GPU Model | Tesla V100
Table : DICE Lab Queue Hardware Specifications (Phase 18b)

The large amounts of memory and CPU cores present on these compute nodes make
them ideal hosts for a pipeline that utilizes parallelization and distributed
computing.

##  6 Portable Batch System (PBS)

Portable Batch System (PBS) is the software that performs job scheduling on
the Palmetto Cluster, as well as the majority of other HPC clusters. The
primary task of PBS is to allocate computational tasks among available
computing resources. Within our project, we use the PBS job array
functionality to several instances of our simulation in parallel. PBS handles
all the resource allocation in the background, and end-users only have to
specify what resources they need and for how long they need said resources
[bayucan_pbs_1998].

##  7 Distributed Computing

In our project, most of the potential applications of this pipeline benefit
immensely from a distributed computing approach. For example, in the
distributed computing problem of generating a massive simulation output
dataset, the dataset would be generated much faster by running two thousand
simulations over fifty computing nodes versus running those same simulations
on one singular computer. At the conceptual level, this is distributed
computing: taking a larger problem and dividing it into many sub-tasks, each
of which can be solved by one or more computers. Lastly, these computers are
usually connected in some fashion, allowing for inter-computer communication.
In the case of most HPC resources, this is inherent
[kshemkalyani_distributed_2011].

##  8 Parallelization

Parallelization, through parallel computing, is the simultaneous use of two or
more compute resources to solve a computational problem. The following are the
four steps of a parallel computing problem [barney_introduction_2010].

  * •

There must be more than one central processing unit (CPU) present.

  * •

A problem is broken down into discrete sections that can be solved
concurrently.

  * •

Each of these parts is then broken down further into a series of instructions.

  * •

Instructions from each part execute simultaneously on separate CPUs.

The parallel environment of Webots.HPC is twofold: first, at the pipeline
level, batches of simulation instances are spread across that available
hardware resources. Then, at the simulation level, Webots engages multi-
threaded processing for each simulation instance. The effectiveness of Webots’
multi-threading is evaluated in Section 5: Performance Evaluation.

##  9 Big Data

Webots.HPC was developed with big data applications in mind. ”Big” data can be
defined as a dataset too large and/or too complex to be processed using
traditional storage and analytic architecture [lantz_potentials_2015]. HPC
resources are very apt for dealing with big data due to their immense
processing capabilities and specialized storage configurations. Additionally,
simulation output datasets can aggregate exponentially; a simulation with a 10
MB output dataset, after being run 100,000 times in sequence, would then swell
to a 1 TB size.

##  10 Chapter Conclusion

Overall, all of these tools, software, and methods came together to produce
Webots.HPC. It is important to note that certain components are not fixed,
meaning that a different commodity could have been used in its place. In
Figure LABEL:fig:arch_diagram, Singularity, X11, and Xvfb are the only fixed
components. The pipeline process behind Webots.HPC can be generalized further
to platform another simulation software in place of the Webots and SUMO
bundle. In addition, Python and Pip could be replaced with another programming
language supported by this new simulation software. Even PBS is not a fixed
component of the pipeline, as some HPC resources have alternative job
scheduling programs installed [goos_job_2001].