#  dame-flame: A Python Library Providing Fast Interpretable Matching for
Causal Inference

Neha R. Gupta<EMAIL_ADDRESS>
Vittorio<EMAIL_ADDRESS>
Chia-Rui Chang<EMAIL_ADDRESS>
Tianyu Wang<EMAIL_ADDRESS>
Marco Morucci<EMAIL_ADDRESS>
Pritam Dey<EMAIL_ADDRESS>
Thomas J. Howell<EMAIL_ADDRESS>
Xian Sun<EMAIL_ADDRESS>
Angikar Ghosal<EMAIL_ADDRESS>
Sudeepa Roy<EMAIL_ADDRESS>
Cynthia Rudin<EMAIL_ADDRESS>
Alexander Volfovsky<EMAIL_ADDRESS>Primary developer and
maintainer of dame-flame

###### Abstract

dame-flame is a Python package for performing _matching_ for _observational
causal inference_ on datasets containing discrete covariates. This package
implements the _Dynamic Almost Matching Exactly (DAME)_ and _Fast, Large-Scale
Almost Matching Exactly (FLAME)_ algorithms, which match treatment and control
units on subsets of the covariates. The resulting matched groups are
interpretable, because the matches are made on covariates (rather than, for
instance, propensity scores), and high-quality, because machine learning is
used to determine which covariates are important to match on. DAME solves an
optimization problem that matches units on as many covariates as possible,
prioritizing matches on important covariates. FLAME approximates the solution
found by DAME via a much faster backward feature selection procedure. The
package provides several adjustable parameters to adapt the algorithms to
specific applications, and can calculate treatment effects after matching.
Descriptions of these parameters, details on estimating treatment effects, and
further examples, can be found in the documentation at https://almost-
matching-exactly.github.io/DAME-FLAME-Python-Package/.

Keywords: matching, causal inference, observational data, interpretability,
Python

##  1 Introduction

We introduce the dame-flame Python package for almost exact matching of
treatment and control units in discrete observational data for causal
analysis.111While handling of continuous observational data is outside the
scope of dame-flame, continuous data can be discretized before applying dame-
flame. This is discussed in the documentation. The advantage of matching is
that it reduces bias of treatment effect estimates, and permits interpretable
analyses that are easier to troubleshoot than other types of analysis for
observational causal studies. dame-flame produces high-quality matches for a
given unit by matching it exactly on as many relevant covariates as possible.
The choice of covariates to match on is aided by machine learning on a holdout
training set.

Many other matching methods either produce low-quality matches (with poor
treatment effect estimates), uninterpretable matches (e.g. in which matches
can group units with highly dissimilar covariates values), or matches that are
manually defined by an analyst. One of the most widely used algorithms is
nearest neighbor propensity score matching, provided by the R package MatchIt
(Stuart et al., 2011). Propensity score matching reduces a dataset to one
dimension, and matches can contain units even at extreme ends of a covariate
space; such matches are uninterpretable. Another common matching algorithm is
Coarsened Exact Matching (CEM), popularly available in the R package cem
(Iacus et al., 2009). CEM requires the user to manually coarsen variables,
which requires substantial human intervention to choose carefully and can be
dominated by irrelevant covariates (Dieng et al., 2019). dame-flame has none
of these disadvantages, providing interpretable matched groups, and using
machine learning on a holdout training set to learn distance metrics.

dame-flame is efficient, owing to fast bit-vector computations and a
combination of backwards feature elimination (for FLAME) and a type of
downwards closure property (for DAME) for feature elimination. The combination
of FLAME (at earlier iterations) and DAME (at later iterations) permits
scalability in matching units on many dimensions in earlier iterations, and
flexibility to almost-exactly match the remaining units on important
covariates in later iterations.

##  2 Algorithm Summary

The dame-flame package implements the algorithms _Dynamic Almost Matching
Exactly (DAME)_ (Dieng et al., 2019) and _Fast Large-Scale Almost Matching
Exactly (FLAME)_ (Wang et al., 2019). Both algorithms iterate over two nested
loops, shown in Figure 1.

The algorithms receive an input dataset with discrete observational
covariates, a binary treatment indicator, and a continuous or discrete outcome
column. It also receives either

Figure 1: Flow diagram: DAME/FLAME algorithms

a separate holdout training dataset with the same covariates as the input
dataset, or it can partition the input dataset to create the holdout set. The
matching set contains the units we aim to match and the holdout training set
is used to determine which variables are important. Users have a variety of
options for handling missing covariate data. They can (1) exclude rows with
missing values from the procedure, (2) impute the missing data via MICE
(Buuren and Groothuis-Oudshoorn, 2010), or (3) specify that matches should not
occur on missing values, without dropping or imputing them.

After handling missing data, the algorithms begin by matching any units that
can be matched exactly on all covariates, where both a treatment and control
unit are contained in each matched group. The algorithms then execute the
outer loop: updating sets of covariates to match on, referred to as active
covariate sets, until a stopping criterion is reached. In each iteration, the
algorithms execute the inner loop, examining each covariate set to select the
best one to match on. Units that have identical values in all of the
covariates that are part of the chosen covariate set form a matched group,
given that at least one treatment unit and at least one control unit are
present in the group.

To determine the best covariate set, FLAME selects the covariates yielding the
highest _match quality_ $\mathtt{MQ}$, defined as
$\mathtt{MQ}=\mathtt{C}\cdot\mathtt{BF}-\mathtt{PE}$, where $\mathtt{C}$ is a
user specified hyperparameter. DAME selects the covariate set minimizing
$\mathtt{PE}$. Here, $\mathtt{PE}$ denotes the _predictive error_. This term
estimates how important the dropped covariates are for predicting the outcome,
based on a user-chosen machine learning algorithm run on the holdout training
set. The _balancing factor_ , $\mathtt{BF}$, measures the proportion of
treatment and control units that are matched on a covariate set. DAME iterates
efficiently over all covariate sets, prioritizing matching on large covariate
sets if they can be used to effectively predict the outcome on the holdout
training set. FLAME approximates this solution for scalability: it
consistently matches on a smaller set of covariates than in the previous
iteration, while still ensuring that each covariate set can be used to
effectively predict the outcome on the holdout training set.

The outer loop has a number of possible stopping criteria. It must stop when
all units are placed in matched groups or all covariate sets have been
dropped. Additionally, users can enforce stopping based off other criteria,
e.g., (1) when there are too few unmatched (treatment or control) units, (2)
after a certain number of iterations, (3) when predictive error rises too
much, or (4) when the balancing factor for a given round is not high enough.

Users can also choose whether units are matched with replacement; that is,
whether a previously matched unit can be matched in a subsequent iteration of
the algorithm.

##  3 Package Usage

### Basic Functionality.

To run DAME or FLAME, the user provides the input data as a data frame or
file, which must contain an outcome column, and a treatment column. Each
algorithm by default produces a table consisting of the units that were
matched to at least one other unit. For each unit that was matched, the table
indicates which of the covariates were used for matching, and the covariate
values that each unit was matched on. The covariates that were not used to
match the unit are denoted with “*” as their values. Users can find matched
groups of individual units and can estimate treatment effects, including the
_average treatment effect_ (ATE) and _conditional average treatment effect_
(CATE) of a selected unit after either algorithm has been run.

Comparison to Other Matching Packages. Many of the features previously
described are unavailable in other matching packages. Table 1 compares the
characteristics of dame-flame against popular alternatives. Most matching
packages are implemented in R. R’s cem package only supports ATT treatment
effects (Iacus et al., 2009). The MatchIt package does not support any
treatment effect calculations, recommending users compute these using a
separate R package, Zelig (Stuart et al., 2011). Users of any propensity score
matching algorithm can adjust matched group sizes only by entering a ratio of
treatment to control units, forcing all matched groups to be of the same size.
Python’s PyMatch and DoWhy offer Propensity Score matching, but DoWhy doesn’t
emphasize matched groups, favoring treatment effects and other output (Sharma
et al., 2019), (Miroglio et al., 2017). A further advantage of dame-flame is
the higher quality of the matched groups generated by DAME and FLAME relative
to propensity score matching, as shown by Dieng et al. (2019).

Additional Advanced Controls of DAME-FLAME. Users are offered several ways to
control the matching procedure, including early stopping critera, and the
option of placing units in matched groups with and without replacement. The
default parameters for the package were chosen for their versatility and
speed, so the algorithm can be relevant and easy-to-use for a range of users.
Full descriptions of algorithm parameters are provided in the respective
package documentations. As part of a suite of printing options, users can
request that the predictive error and/or balancing factor at each iteration be
returned.

Language: Package |  built-in treatment-effect estimations |  flexible matched group sizes |  missing data handling options |  provide matched groups
---|---|---|---|---
Python: dame-flame |  ATE, ATT, CATE |  |  |
Python: DoWhy |  ATT, ATE, ATC |  |  |
Python: PyMatch |  |  |  |
R: cem |  ATT |  |  |
R: MatchIt Propensity Score |  |  |  |
Table 1: Features of Matching Packages

##  4 Installation, Documentation, and Project Management

dame-flame is available for download on PyPi and on GitHub222 Code is publicly
available at https://github.com/almost-matching-exactly/DAME-FLAME-Python-
Package/.. On GitHub, users may report bugs, ask questions, and make pull
requests. The documentation covers the API, installation instructions, a
quick-start tutorial, several examples, and a contributing guide. The package
is designed for Python 3.6 and above. dame-flame depends on scikit-learn
version 0.23.2 and above, pandas version 0.11.0 and above, and numpy version
1.16.5 and above. The API consists of standard Pythonic design established by
scikit-learn, such as the methods fit, where a user provides a holdout
training dataset, and predict, where a user provides a matching dataset, and
the matches of the algorithm are computed.

There is an accompanying experimental R package for the FLAME algorithm, whose
development version can also be found on GitHub.333https://github.com/almost-
matching-exactly/R-FLAME Testing was done to ensure that R-FLAME and dame-
flame yield consistent results. Furthermore, dame-flame offers Continuous
Integration, and the extensive test suite has a code coverage of 98%.

##  5 Conclusion

The dame-flame package offers efficient, easy-to-use implementations of the
DAME and FLAME algorithms, allowing users to perform fast, interpretable
matching for causal inference for observational data with discrete covariates.
The package is easily accessible and accompanied by detailed documentation,
with concrete examples. The package is written in a highly modular manner,
facilitating the introduction of new features and variations of the DAME and
FLAME algorithms.

Acknowledgments: This work was supported in part by awards NIH R01EB025021,
NSF IIS-1703431, and the Duke University Energy Initiative Energy Research
Seed Fund.

## References

  * Buuren and Groothuis-Oudshoorn (2010) Stef van Buuren and Karin Groothuis-Oudshoorn.  mice: Multivariate imputation by chained equations in R.  _Journal of Statistical Software_ , pages 1–68, 2010.
  * Dieng et al. (2019) Awa Dieng, Yameng Liu, Sudeepa Roy, Cynthia Rudin, and Alexander Volfovsky.  Interpretable Almost-Exact Matching for Causal Inference.  In _The 22nd International Conference on Artificial Intelligence and Statistics, AISTATS 2019, 16-18 April 2019, Naha, Okinawa, Japan_ , pages 2445–2453, 2019.
  * Iacus et al. (2009) Stefano Iacus, Gary King, and Giuseppe Porro.  cem: Software for coarsened exact matching.  _Journal of Statistical Software, Articles_ , 30(9):1–27, 2009.  ISSN 1548-7660.  doi: 10.18637/jss.v030.i09.  URL https://www.jstatsoft.org/v030/i09.
  * Miroglio et al. (2017) Ben Miroglio et al.  pymatch.  https://github.com/benmiroglio/pymatch, 2017.
  * Sharma et al. (2019) Amit Sharma, Emre Kiciman, et al.  DoWhy: A Python package for causal inference.  https://github.com/microsoft/dowhy, 2019.
  * Stuart et al. (2011) Elizabeth A Stuart, Gary King, Kosuke Imai, and Daniel Ho.  Matchit: nonparametric preprocessing for parametric causal inference.  _Journal of Statistical Software_ , 2011.
  * Wang et al. (2019) Tianyu Wang, Marco Morucci, M. Usaid Awan, Yameng Liu, Sudeepa Roy, Cynthia Rudin, and Alexander Volfovsky.  FLAME: A Fast Large-scale Almost Matching Exactly Approach to Causal Inference.  _CoRR_ , abs/1707.06315v7, 2019.  URL http://arxiv.org/abs/1707.06315.


# CECILIA: Comprehensive Secure Machine Learning Framework

Ali Burak Ünal 1,2, Nico Pfeifer 1, Mete Akgün 1,2

###### Abstract

Since machine learning algorithms have proven their success in many different
applications, there is also a big interest in privacy preserving machine
learning methods for building models on sensitive data. Moreover, the increase
in the number of data sources and the high computational power required by
those algorithms force individuals to outsource the training and/or the
inference of a machine learning model to the clouds providing such services.
To address this, we propose a secure three party computation framework,
CECILIA, offering privacy preserving building blocks to enable complex
operations privately. In addition to the adapted and common operations like
addition and multiplication, it offers multiplexer, most significant bit and
modulus conversion. The first two are novel in terms of methodology and the
last one is novel in terms of both functionality and methodology. CECILIA also
has two complex novel methods, which are the exact exponential of a public
base raised to the power of a secret value and the inverse square root of a
secret Gram matrix. We use CECILIA to realize the private inference on pre-
trained recurrent kernel networks, which require more complex operations than
most other deep neural networks, on the structural classification of proteins
as the first study ever accomplishing the privacy preserving inference on
recurrent kernel networks. In addition to the successful private computation
of basic building blocks, the results demonstrate that we perform the exact
and fully private exponential computation, which is done by approximation in
the literature so far. Moreover, they also show that we compute the exact
inverse square root of a secret Gram matrix up to a certain privacy level,
which has not been addressed in the literature at all. We also analyze the
scalability of CECILIA to various settings on a synthetic dataset. The
framework shows a great promise to make other machine learning algorithms as
well as further computations privately computable by the building blocks of
the framework.

## Introduction

In recent years, machine learning algorithms, especially deep learning
algorithms, have become more sophisticated and are now more data-driven to
achieve better performance. At the same time, the number of sources generating
data with sensitive information and the amount of that data have increased
dramatically. By its very nature, the privacy of this data must be protected
during collaborative model training and/or inference. Likewise, the privacy of
such a model must be kept private when the owner of the trained model deploys
it as a prediction service.

To preserve the privacy of both the data and the model, there are several
privacy techniques utilized in the literature. One of these techniques is
differential privacy (DiP) introduced by Dwork, Roth et al. (2014). In DiP,
the input data, the model parameters, and/or the output are perturbed to
protect the privacy by adding noise in a certain budget that adjusts the
privacy level (Abadi et al. 2016; Chen et al. 2020). By nature, noise addition
sacrifices the performance of the model and exactness of the result to some
extent. A cryptographic technique called homomorphic encryption (HE) addresses
the missing points of DiP in addition to the basic privacy requirements. HE
protects the privacy of the data and/or the model by encrypting the data
and/or the parameters of the model with different key schemes. Thanks to this
mechanism, various operations such as addition and multiplication can be
performed in the encrypted domain. There are several attempts to implement
machine learning algorithms using HE (Bakshi and Last 2020; Hesamifard,
Takabi, and Ghasemi 2017; Gilad-Bachrach et al. 2016; Lu et al. 2021).
However, its huge runtime and the limited number of practical operations that
can be realized with HE makes these computations impractical. Secure multi-
party computation (MPC), on the other hand, addresses the missing points of HE
and the fundamental basic privacy requirements. The data and/or the model
parameters are shared among several parties in such a way that none of the
parties can learn about the data and/or the model parameters on their own.
Then, these parties perform the desired computation privately. To address
various machine learning algorithms, there are several MPC frameworks in the
literature (Knott et al. 2021; Wagh, Gupta, and Chandran 2019; Mohassel and
Zhang 2017; Damgård et al. 2012; Rathee et al. 2020; Wagh et al. 2020; Patra
et al. 2021) some of which also utilizes HE (Mishra et al. 2020; Huang et al.
2022). Although they have some efficient and secure basic functions for
several types of deep neural networks, from which we have also benefited, they
lack some complex functions completely such as the inverse square root of a
Gram matrix or the exact computation of some functions such as the exponential
to realize the exact inference on some deep neural networks like recurrent
kernel networks (Chen, Jacob, and Mairal 2019).

In this study, we introduce a general purpose efficient secure multi-party
computation framework, CECILIA, based on $3$ computing parties. To the best of
our knowledge, the novelty of CECILIA over the other frameworks stems the
privacy preserving exact computation of the exponential of a known base raised
to the power of a secret shared value, the computation of the inverse square
root of a secret shared Gram matrix with a high probability, a new privacy
preserving modulus conversion from $2^{n-1}$ to $2^{n}$, efficient private
computation of the most significant bit and novel methodology of performing
multiplexer in a privacy preserving way. By using CECILIA, we perform privacy
preserving inference on a pre-trained recurrent kernel network, which has more
complex operations than most of the other deep neural networks and cannot be
realized using any other framework in the literature. As a summary, CECILIA

  * •

provides novel and/or efficient private computations for the most significant
bit, multiplexer and modulus conversion to address the functions used in
machine learning algorithms.

  * •

introduces the exact privacy preserving computation of the exponential of a
known base raised to the power of a secret shared value for the first time in
the literature.

  * •

enables the privacy preserving inverse square root of a secret shared Gram
matrix with a high probability as the first attempt in the literature.

  * •

realizes an efficient privacy preserving inference on recurrent kernel
networks for the first time.

## Preliminaries

### Security Model

To prove the security of CECILIA, we use the honest-but-curious security
model, in other words we analyze the security of the framework when there is a
semi-honest adversary who corrupts a party and follows the protocols honestly,
but at the same time tries to infer information during the execution of these
protocols. In our analyses of the building blocks of CECILIA, we consider a
scenario where a semi-honest adversary corrupts either one of the proxies or
the helper party and tries to infer the input and/or output values of the
executed function. We take the same approach to analyze the security of
privacy preserving RKNs and focus on if the semi-honest adversary corrupting a
party can break the privacy of the test sample of the data owner and/or the
model parameters of the model owner during the prediction of a test sample
using the outsourced model.

### Square-and-multiply Algorithm for Exponential

One of the approaches to computing the exponential of a base raised to the
power of a value that can be represented in a binary form is the square-and-
multiply algorithm. Let $b$ be the base and $a$ be the power, whose binary
representation is $\langle a\rangle$ of length n. Then, one can compute
$b^{a}$ as follows:

$b^{a}=\prod_{i=0}^{n-1}b^{\langle a\rangle_{i}\cdot 2^{i}}$ (1)

where $\langle a\rangle_{i}$ represents the bit value of $a$ at position $i$,
assuming that the indexing starts from the least significant bit, that is
$\langle a\rangle_{0}$ corresponds to the least significant bit and $\langle
a\rangle_{n-1}$ corresponds to the most significant bit.

### Recurrent Kernel Networks

Chen, Jacob, and Mairal (2019) gave a kernel perspective of RNNs by showing
that the computation of the specific construction of RNNs, which they call
recurrent kernel networks (RKNs), mimics the substring kernel allowing
mismatches and the local alignment kernel, which are widely used on sequence
data (EL-Manzalawy, Dobbs, and Honavar 2008; Nojoomi and Koehl 2017; Leslie,
Eskin, and Noble 2001). RKNs use a well-designed kernel formulation and enable
the parameters of the model to be optimized via backpropagation. As a result,
they outperform the traditional substring kernel and the local alignment
kernel as well as the LSTMs (Hochreiter, Heusel, and Obermayer 2007). Such
performance promises to have a large impact on the results of tasks with small
to medium size data.

In RKNs, small motifs, called anchor points, are used as templates to measure
similarities among sequences and each character of an anchor point is encoded
via an encoding similar to one-hot encoding. Let there be $q$ anchor points of
length $k$, each of whose characters is encoded using a vector of size $d$,
and a sequence $x$ of length $s$, whose characters are encoded using one-hot
encoded vectors of length $d$. For the $t$-th character of the input sequence,
RKNs first compute the similarity of this character to each character of all
anchor points by the following:

$K(x_{t},z_{j}^{i})=e^{\alpha(\langle x_{t},z_{j}^{i}\rangle-1)}$ (2)

where $\alpha$ is the similarity measure parameter, $x_{t}$ is the
$d$-dimensional vector representing the $t$-th character of the sequence $x$
and $z_{j}^{i}$ is the $d$-dimensional vector representing the $i$-th
character of the $j$-th anchor point. Once the similarity of the $t$-th
character of the sequence to the $j$-th character of all $q$ anchor points is
computed and represented as a vector $b_{j}[t]$, the computation continues as
follows:

$c_{j}[t]=\lambda c_{j}[t-1]+c_{j-1}[t-1]\otimes b_{j}[t]$ (3)

where $c_{j}[t]$ is the initial mapping of the sequence up to its $t$-th
character into a $q$-dimensional vector based on anchor points of length $j$
for $j\in\\{1,\ldots,k\\}$ and $t\in\\{1,\ldots,s\\}$, $\lambda$ is a scalar
value to downgrade the effect of the previous time points and $\otimes$ is the
Hadamard product. As the base of this recursive equation, $c_{0}[t]$ is a
vector of $1$s and $c_{j}[0]$ is a vector of $0$s if $j\neq 0$.

To obtain the final mapping of the sequence, RKN multiplies $c_{j}[s]$ by the
inverse square root matrix of the Gram matrix of the anchor points up to their
$j$-th characters. Afterwards, in the substring kernel allowing mismatches,
they use a linear layer, which is a dot product between $c_{k}[s]$ and the
weight vector of the classifier $w$, to obtain the prediction of the input
sequence. Figure 1 depicts these computations.

(a)

(b)

Figure 1:  The arithmetic circuit of the computation in (a) a single neuron of
RKN at position $t$ and $k$-mer level $k$ and (b) the linear classifier layer
of RKN after the last position of the input sequence $x$. Black lines are a
scalar value. Orange lines are $d$-dimensional vectors and blue lines are
$q$-dimensional vectors. Purple ones are $(q\times q)$-dimensional matrix

## Framework

CECILIA is based on three computing parties and 2-out-of-2 additive secret
sharing. Two of these parties, $P_{0}$ and $P_{1}$, are called proxy and the
external entities such as the model owner and the data sources interact with
them. The third one, $P_{2}$, is the helper party helping the proxies compute
the desired function without breaking the privacy. It provides the proxies
with the shares of the purposefully designed values. It also performs
calculations on the data masked by the proxies in the real number domain and
returns the share of the results of these calculations. In the rest of this
section, we describe the number format we use to represent values. We then
introduce the methods of CECILIA.

### Notations

To ensure the privacy in CECILIA, we use 2-out-of-2 additive secret sharing
over three different rings $\mathbb{Z}_{L}$, $\mathbb{Z}_{K}$ and
$\mathbb{Z}_{V}$ where $L=2^{n}$, $K=2^{n-1}$, $V=67$ and $n=64$. We denote
two shares of $x$ over $\mathbb{Z}_{L}$, $\mathbb{Z}_{K}$ and $\mathbb{Z}_{V}$
by ($\langle x\rangle_{0}$, $\langle x\rangle_{1}$), ($\langle
x\rangle_{0}^{K}$, $\langle x\rangle_{1}^{K}$) and ($\langle
x\rangle_{0}^{V}$, $\langle x\rangle_{1}^{V}$), respectively. If a value $x$
is shared over the ring $\mathbb{Z}_{V}$, every bit of $x$ is additively
shared in $\mathbb{Z}_{V}$. This means $x$ is shared as a vector of $n$ shares
where each share takes a value between $0$ and $(V-1)$. We also use boolean
sharing of a single bit denoted by ($\langle x\rangle_{0}^{B}$, $\langle
x\rangle_{1}^{B}$).

### Number Format

Since machine learning algorithms can require both positive and negative real
numbers, CECILIA should be able to work with these numbers. To achieve this,
we follow the number format proposed by Mohassel and Zhang (2017) in SecureML.
In this number format, the most significant bit indicates the sign. $0$ and
$1$ indicate positive and negative values, respectively. Moreover, a certain
number of least significant bits of this representation are for the fractional
part of the value and the rest expresses the integer part.

### Addition (ADD)

In CECILIA, the proxies add the shares of two secret shared values they have
to obtain the share of the addition of these values without any communication
and privacy leakage.

### Multiplication (MUL)

The multiplication operation, which we adapted from SecureML, uses the pre-
computed multiplication triples (Beaver 1991) and requires truncation because
of the special number format. For more details, refer to the Supplement and
see (Mohassel and Zhang 2017; Wagh, Gupta, and Chandran 2019).

1 Algorithm  _$\mathsf{MOC}$(__)_

input : $P_{0}$ and $P_{1}$ hold $\langle x\rangle_{0}^{K}$ and $\langle
x\rangle_{1}^{K}$, respectively

output : $P_{0}$ and $P_{1}$ get $\langle x\rangle_{0}$ and $\langle
x\rangle_{1}$, respectively

2 $P_{0}$ and $P_{1}$ hold a common random bit $u^{\prime}$

3 $P_{2}$ picks a random numbers $r\in\mathbb{Z}_{K}$ and generates $\langle
r\rangle_{0}^{K}$, $\langle r\rangle_{1}^{K}$, $\\{\langle
r[j]\rangle_{0}^{V}\\}_{j\in[n]}$ and $\\{\langle
r[j]\rangle_{1}^{V}\\}_{j\in[n]}$.

4 $P_{2}$ computes $w=\text{isWrap}(\langle r\rangle_{0}^{K},\langle
r\rangle_{1}^{K},K)$ and divides $w$ into two boolean shares $w_{0}^{B}$ and
$w_{1}^{B}$

5 $P_{2}$ sends $\langle r\rangle_{i}^{K}$, $\\{\langle
r[j]\rangle_{i}^{V}\\}_{j\in[n]}$ and $w_{i}^{B}$ to $P_{i}$, for each
$i\in\\{0,1\\}$

6 For each $i\in\\{0,1\\}$, $P_{i}$ executes Steps $7$-$8$

7 $\langle y\rangle_{i}^{K}=\langle x\rangle_{i}^{K}+\langle r\rangle_{i}^{K}$

8 $P_{i}$ reconstructs $y$ by exchanging shares with $P_{1-i}$

9 $u_{i}^{B}=\mathsf{PC}(\\{\langle
r[j]\rangle_{i}^{V}\\}_{j\in[n]},y,u^{\prime})$

10 $P_{0}$ computes $u_{i}^{B}=u_{i}^{B}\oplus u^{\prime}$

11 For each $i\in\\{0,1\\}$, $P_{i}$ computes $c_{i}^{B}=w_{i}^{B}\oplus
u_{i}^{B}$

12 $P_{0}$ computes $\langle y\rangle_{0}=\langle
y\rangle_{0}^{K}+\text{isWrap}(\langle y\rangle_{0}^{K},\langle
y\rangle_{1}^{K},K)\cdot K$

13 $P_{1}$ sets $\langle y\rangle_{1}=\langle y\rangle_{1}^{K}$

14 For each $i\in\\{0,1\\}$, $P_{i}$ computes $\langle x\rangle_{i}=\langle
y\rangle_{i}-(\langle r\rangle_{i}^{K}+c_{i}^{B}\cdot K)$

15

Algorithm 1 Modulus Conversion ($\mathsf{MOC}$)

### Modulus Conversion (MOC)

CECILIA offers the functionality $\mathsf{MOC}$ described in Algorithm 1. It
converts shares over $\mathbb{Z}_{K}$ to fresh shares over $\mathbb{Z}_{L}$
where $L=2K$. Even though other frameworks in the literature have functions
with a similar name, none of those performs this specific modulus conversion.
To the best of our knowledge, CECILIA is the first framework offering this
specific modulus conversion. Assuming that $P_{0}$ and $P_{1}$ have the shares
$\langle x\rangle_{0}^{K}$ and $\langle x\rangle_{1}^{K}$, respectively, the
first step for $P_{0}$ and $P_{1}$ is to mask their shares by using the shares
of the random value $r\in\mathbb{Z}_{K}$ sent by $P_{2}$. Afterwards, they
reconstruct $(x+r)\in\mathbb{Z}_{K}$ by first computing $\langle
y\rangle_{i}^{K}=\langle x\rangle_{i}^{K}+\langle r\rangle_{i}^{K}$ for
$i\in\\{0,1\\}$ and then sending these values to each other. Along with the
shares of $r\in\mathbb{Z}_{K}$, $P_{2}$ also sends the information in boolean
shares telling whether the summation of the shares of $r$ wraps so that
$P_{0}$ and $P_{1}$ can convert $r$ from the ring $\mathbb{Z}_{K}$ to the ring
$\mathbb{Z}_{L}$. Once they reconstruct $y\in\mathbb{Z}_{K}$, $P_{0}$ and
$P_{1}$ can change the ring of $y$ to $\mathbb{Z}_{L}$ by adding $K$ to one of
the shares of $y$ if $\langle y\rangle_{0}^{K}+\langle y\rangle_{1}^{K}$
wraps. After conversion, the important detail regarding $y\in\mathbb{Z}_{L}$
is to fix the value of it. In case $(x+r)\in\mathbb{Z}_{K}$ wraps, which we
identify by using the private compare ($\mathsf{PC}$) method (Wagh, Gupta, and
Chandran 2019), depending on the boolean share of the outcome of
$\mathsf{PC}$, $P_{0}$ or $P_{1}$ or both add $K$ to their shares. If both
add, this means that there is no addition to the value of
$y\in\mathbb{Z}_{L}$. At the end, $P_{i}$ subtracts $r_{i}\in\mathbb{Z}_{L}$
from $y_{i}\in\mathbb{Z}_{L}$ and obtains $x_{i}\in\mathbb{Z}_{L}$ for
$i\in\\{0,1\\}$.

### Most Significant Bit (MSB)

One of the biggest improvements that CECILIA introduces is in the private
determination of the most significant bit of a secret shared value $x$ via
$\mathsf{MSB}$. We deftly integrated $\mathsf{MOC}$ and PC into $\mathsf{MSB}$
to reduce communication rounds. Given the shares of $x$, $\mathsf{MSB}$ first
extracts the least significant $(n-1)$-bits via $mod\ K$. Then, it converts
the ring of this value from $\mathbb{Z}_{K}$ to $\mathbb{Z}_{L}$.
$\mathsf{MSB}$ then subtracts it from $x$, which results in either $0$ or $K$
in $\mathbb{Z}_{L}$. Finally, it secretly maps this result to $0$ or $1$ and
obtains the most significant bit of $x$ privately.

1 Algorithm  _$\mathsf{MSB}$(__)_

input : $P_{0}$ and $P_{1}$ hold $\langle x\rangle_{0}$ and $\langle
x\rangle_{1}$, respectively

output : $P_{0}$ and $P_{1}$ get $\langle z\rangle_{0}$ and $\langle
z\rangle_{1}$, respectively, where $z$ is equal to $0$ if the most significant
bit of $x$ is $0$ and $1$ otherwise.

2 $P_{0}$ and $P_{1}$ hold a common random bit $f$ and $g$. $P_{0}$ and
$P_{1}$ additionally hold $\ell$ common random values
$s_{j}\in\mathbb{Z}_{V}^{*}$ for all $j\in[\ell]$, a random permutation $\pi$
for $\ell$ elements and $\ell$ common random values
$u_{j}\in\mathbb{Z}_{V}^{*}$.

3 $P_{2}$ picks a random number $r\in K$ and generates $\langle
r\rangle_{0}^{K}$, $\langle r\rangle_{1}^{K}$, $\\{\langle
r[j]\rangle_{0}^{V}\\}_{j\in[\ell]}$ and $\\{\langle
r[j]\rangle_{1}^{V}\\}_{j\in[\ell]}$.

4 $P_{2}$ computes $w=\text{isWrap}(\langle r\rangle_{0}^{K},\langle
r\rangle_{1}^{K},K)$

5 $P_{2}$ sends $\langle r\rangle_{i}^{K}$ and $\\{\langle
r[j]\rangle_{i}^{V}\\}_{j\in[\ell]}$ to $P_{i}$, for each $i\in\\{0,1\\}$

6 For each $i\in\\{0,1\\}$, $P_{i}$ executes Steps $7$-$26$

7 $\langle d\rangle_{i}^{K}=\langle x\rangle_{i}\mod K$

8 $\langle y\rangle_{i}^{K}=\langle d\rangle_{i}^{K}+\langle r\rangle_{i}^{K}$

9 $y=\mathsf{Reconst(\langle y\rangle_{i}^{K})}$

10 $\langle y\rangle_{i}=\langle y\rangle_{i}^{K}+i\cdot\text{isWrap}(\langle
y\rangle_{0}^{K},\langle y\rangle_{1}^{K},K)\cdot K$

11 $\langle a[0]\rangle_{i}=ifK-\langle x\rangle_{i}+\langle
y\rangle_{i}-\langle r\rangle_{i}^{K}$

12 $\langle a[1]\rangle_{i}=i(1-f)K-\langle x\rangle_{i}+\langle
y\rangle_{i}-\langle r\rangle_{i}^{K}$

13

14 Let $t=y+1\bmod 2^{\ell}$

15 for  _$j=\ell;\ j >0;\ j=j-1$_ do

16 if  _$g=0$_ then

17 $\left\langle w_{j}\right\rangle_{i}^{V}=\langle
r[j]\rangle_{i}^{V}+iy[j]-2y[j]\langle r[j]\rangle_{i}^{V}$

18 $\left\langle c_{j}\right\rangle_{i}^{V}=iy[j]-\langle
r[j]\rangle_{i}^{V}+j+\sum_{k=j+1}^{\ell}\left\langle
w_{k}\right\rangle_{i}^{V}$

19

20 else if  _$g=1$ AND $r\neq 2^{\ell}-1$_ then

21 $\left\langle w_{j}\right\rangle_{i}^{V}=\langle
r[j]\rangle_{i}^{V}+it[j]-2t[j]\langle r[j]\rangle_{i}^{V}$

22 $\left\langle c_{j}\right\rangle_{i}^{V}=-it[j]+\langle
r[j]\rangle_{i}^{V}+i+\sum_{k=j+1}^{\ell}\left\langle
w_{k}\right\rangle_{i}^{V}$;

23 else

24 if  _$i\neq 1$_ then

25 $\left\langle c_{j}\right\rangle_{i}^{V}=(1-i)\left(u_{j}+1\right)-iu_{j}$

26

27 else

28 $\left\langle c_{j}\right\rangle_{i}^{V}=(-1)^{j}\cdot u_{j}$

29

30

31

32

33 Send $\left\\{\left\langle
b_{j}\right\rangle_{i}^{V}\right\\}_{j}=\pi\left(\left\\{s_{j}\left\langle
c_{j}\right\rangle_{i}^{V}\right\\}_{j}\right)$ and $\langle a\rangle_{i}$ to
$P_{2}$

34 For all $j\in{[\ell]}$, $P_{2}$ computes $d_{j}=\mathsf{Reconst}(\langle
d_{j}\rangle_{0}^{V},\langle d_{j}\rangle_{1}^{V})$ and sets $g^{\prime}=1$
iff $\exists j\in[\ell]$ such that $d_{j}=0$.

35 $P_{2}$ reconstructs $a[j]$ where $j\in\\{0,1\\}$ and computes
$a[j]=(a[j]-(g^{\prime}\oplus w)K)/K$

36 $P_{2}$ creates two fresh shares of $a[j]$ where $j\in\\{0,1\\}$ and sends
them to $P_{0}$ and $P_{1}$

37 For each $i\in\\{0,1\\}$, $P_{i}$ executes Step $31$

38 $\langle z\rangle_{i}=\langle a[f\oplus g]\rangle_{i}$

39

Algorithm 2 Most Significant Bit ($\mathsf{MSB}$)

### Comparison of Two Secret Shared Values (CMP)

CECILIA also provides $\mathsf{CMP}$ which privately compares two secret
shared values $x$ and $y$, and outputs $0$ if $x\geq y$, $1$ otherwise in
secret shared form. $\mathsf{CMP}$ utilizes $\mathsf{MSB}$ to determine the
most significant bit of $x-y$. Then, it outputs $0$ if $\mathsf{MSB}$ gives
$0$, $1$ otherwise.

### Multiplexer (MUX)

CECILIA also provides the functionality $\mathsf{MUX}$ described in Algorithm
3. It performs the selection of one of two secret shared values based on the
secret shared bit value using the randomized encoding (RE) of multiplication
(Applebaum 2017). In $\mathsf{MUX}$, the proxies compute Equation 4 using
their shares ($\langle x\rangle_{0}$, $\langle y\rangle_{0}$, $\langle
b\rangle_{0}$) and ($\langle x\rangle_{1}$, $\langle y\rangle_{1}$, $\langle
b\rangle_{1}$), respectively, and obtain the fresh shares of $z=x-b(x-y)$. As
shown in Equation 4, the proxies need to multiply two values owned by
different parties in the computation of $\langle b\rangle_{0}(\langle
x\rangle_{1}-\langle y\rangle_{1})$ and $\langle b\rangle_{1}(\langle
x\rangle_{0}-\langle y\rangle_{0}$. They outsource these multiplications to
the helper via the RE. To the best of our knowledge, it is the first time that
the RE is used in such a functionality.

$\begin{split}z={}&x-b(x-y)\\\ ={}&\langle x\rangle_{0}+\langle
x\rangle_{1}-\langle b\rangle_{0}(\langle x\rangle_{0}-\langle
y\rangle_{0})-\langle b\rangle_{1}(\langle x\rangle_{1}-\langle
y\rangle_{1})\\\ &-\langle b\rangle_{0}(\langle x\rangle_{1}-\langle
y\rangle_{1})-\langle b\rangle_{1}(\langle x\rangle_{0}-\langle
y\rangle_{0})\end{split}$ (4)

1 Algorithm  _$\mathsf{MUX}$(__)_

input : $P_{0}$ and $P_{1}$ hold ($\langle x\rangle_{0}$, $\langle
y\rangle_{0}$, $\langle b\rangle_{0}$) and ($\langle x\rangle_{1}$, $\langle
y\rangle_{1}$, $\langle b\rangle_{1}$), respectively.

output : $P_{0}$ and $P_{1}$ get $\langle z\rangle_{0}$ and $\langle
z\rangle_{1}$, respectively, where $z=x-b(x-y)$.

2 $P_{0}$ and $P_{1}$ hold four common random values $r_{i}$ where
$i\in\\{0,1,2,3\\}$

3 $P_{0}$ computes $M_{1}=\langle x\rangle_{0}-\langle b\rangle_{0}(\langle
x\rangle_{0}-\langle y\rangle_{0})+r_{1}\langle b\rangle_{0}+r_{2}(\langle
x\rangle_{0}-\langle y\rangle_{0})+r_{2}r_{3}$, $M_{2}=\langle
b\rangle_{0}+r_{0}$, $M_{3}=\langle x\rangle_{0}-\langle y\rangle_{0}+r_{3}$

4 $P_{0}$ sends $M_{2}$ and $M_{3}$ to $P_{2}$

5 $P_{1}$ computes $M_{4}=\langle x\rangle_{1}-\langle b\rangle_{1}(\langle
x\rangle_{1}-\langle y\rangle_{1})+r_{0}(\langle x\rangle_{1}-\langle
y\rangle_{1})+r_{0}r_{1}+r_{3}\langle b\rangle_{1}$, $M_{5}=(\langle
x\rangle_{1}-\langle y\rangle_{1})+r_{1}$, $M_{6}=\langle b\rangle_{1}+r_{2}$

6 $P_{1}$ sends $M_{5}$ and $M_{6}$ to $P_{2}$

7 $P_{2}$ computes $M_{2}M_{5}+M_{3}M_{6}$ $=z$

8 $P_{2}$ divides $z$ into two shares $(\langle z\rangle_{0}+\langle
z\rangle_{1})$ and sends $\langle z\rangle_{0}$ and $\langle z\rangle_{1}$ to
$P_{0}$ and $P_{1}$, respectively

9 $P_{0}$ computes $\langle z\rangle_{0}=M_{1}-\langle z\rangle_{0}$

10 $P_{1}$ computes $\langle z\rangle_{1}=M_{4}-\langle z\rangle_{1}$

11

Algorithm 3 Multiplexer ($\mathsf{MUX}$)

### Dot Product (DP)

Since the dot product between two vectors is a widely used operation in
machine learning algorithms, CECILIA provides $\mathsf{DP}$. It essentially
uses $\mathsf{MUL}$ to compute the element-wise multiplication of the vectors
and $\mathsf{ADD}$ to sum the elements of the resulting vector to obtain the
dot product result.

### Exponential Computation (EXP)

One of the novel methods of CECILIA is functionality $\mathsf{EXP}$ given in
4. It computes the exponential of a publicly known base raised to the power of
a given secret shared value. For this purpose, we have been inspired by the
square-and-multiply algorithm. To the best of our knowledge, CECILIA is the
first study that extends the core idea of the square-and-multiply algorithm to
cover not only the positive numbers, but also the negative numbers as well as
their decimal parts in a multi-party scenario and performs the exact
exponential computation. As an overview, the proxies first obtain the most
significant bit of the secret shared power and use this to select the set
containing either the power itself and the contribution of each bit of a
positive power or the absolute of the power and the contribution of each bit
of a negative power. Then, the proxies determine the value of each bit of the
power in a secret sharing form and use them to select between the previously
selected contributions of the bits and a vector of $1$s. The last step is to
multiply these selected contribution of the bits of the power to the
exponential in a binary-tree-like-structure. In total, $\mathsf{EXP}$ requires
two $\mathsf{MSB}$, two $\mathsf{MUX}$ and $log_{2}(n)$-many $\mathsf{MUL}$
calls.

1 Algorithm  _$\mathsf{EXP}$(__)_

input : $P_{0}$ and $P_{1}$ hold $\langle x\rangle_{0}$ and $\langle
x\rangle_{1}$, respectively, and publicly known base $b$

output : $P_{0}$ and $P_{1}$ get $\langle z\rangle_{0}$ and $\langle
z\rangle_{1}$, respectively

2

3 For $i\in\\{0,1\\}$, $P_{i}$ executes Steps 3-15 with the help of $P_{2}$

4

5 $\langle s\rangle_{i}=\mathsf{MSB}(\langle x\rangle_{i})$

6 $\langle|x|\rangle_{i}=0-\langle x\rangle_{i}$

7 for  _$j=n;\ j >0;\ j=j-1$_ do

8 $\langle cP[n-j]\rangle_{i}=i*b^{2^{j-f}}$

9 $\langle cN[n-j]\rangle_{i}=i*(1/b^{2^{j-f}}$)

10

11 $(\langle\hat{x}\rangle_{i},\langle cONE\rangle_{i})=\mathsf{MUX}((\langle
x\rangle_{i},\langle cP\rangle_{i}),(\langle|x|\rangle_{i},\langle
cN\rangle_{i}),\langle s\rangle_{i})$

12 for  _$j=0;\ j <n;\ j=j+1$_ do

13 $\langle\hat{x}^{n-j}\rangle_{i}=(\langle\hat{x}\rangle_{i}<<j)$

14

15 $\langle
M\rangle_{i}=\mathsf{MSB}((\langle\hat{x}^{64}\rangle_{i},\langle\hat{x}^{63}\rangle_{i},\ldots,\langle\hat{x}^{1}\rangle_{i}))$

16 $\langle cACT\rangle_{i}=\mathsf{MUX}(\langle cONE\rangle_{i},\langle
cZERO\rangle_{i},\langle M\rangle_{i})$

17 $\langle z\rangle_{i}=\langle cACT\rangle_{i}$

18 for  _$j=0;\ j <log_{2}(n);\ j=j+1$_ do

19 $\langle z\rangle_{i}=\mathsf{MUL}(\langle z[0:len(z)/2]\rangle_{i},\langle
z[len(z)/2:len(z)]\rangle_{i})$

20

Algorithm 4 Exponential computation ($\mathsf{EXP}$)

### Inverse Square Root of Gram Matrix (INVSQRT)

As a unique feature, CECILIA provides a more sophisticated and special
operation, namely the inverse square root of a secret shared Gram matrix,
whose pseudocode is given in Algorithm 5. To the best of our knowledge,
CECILIA is the first study in the literature offering such a functionality. In
$\mathsf{INVSQRT}$, we follow the idea of computing the reciprocal of the
square root of the eigenvalues of the Gram matrix and constructing the desired
matrix using the diagonal matrix of these reciprocals of the square root of
the eigenvalues and the original eigenvectors. For the privacy preserving
eigenvalue decomposition, we were inspired by the approach proposed by Zhou
and Li (2016), where the eigenvalue decomposition was designed as an
outsourcing operation from a single party. We adapted this approach to the
multi-party scenario so that the proxies with a share of Gram matrix can
perform the eigenvalue decomposition without learning anything other than
their share of the outcome.

Let $i\in\\{0,1\\}$, $G_{i}$ be the share of the Gram matrix in $P_{i}$,
$\mathcal{M}$ and $(\tau,s)$ be a common orthogonal matrix and common scalars,
respectively, known by both proxies. $P_{i}$ first masks its share of the Gram
matrix by computing $G^{{}^{\prime}}_{i}=\mathcal{M}(\tau
G_{i}+sI)\mathcal{M}^{T}$, where $I$ is the identity matrix, and then sends it
to the helper party, or $P_{2}$. First, $P_{2}$ reconstructs $G^{{}^{\prime}}$
and then performs an eigenvalue decomposition on $G^{{}^{\prime}}$. It obtains
the masked versions of the original eigenvalues, represented as
$\Lambda^{{}^{\prime}}$, and eigenvectors, represented as $Q^{{}^{\prime}}$.
After splitting $Q^{{}^{\prime}}$ into two secret shares, $P_{2}$ generates a
vector of random values $\Delta$ and a random scalar $\alpha$, and sends them
to $P_{1}$ along with $Q^{{}^{\prime}}_{1}$. Once $P_{1}$ has received
$Q^{{}^{\prime}}_{1}$, $\Delta$ and $\alpha$, it computes
$Q_{1}=\mathcal{M}^{T}Q^{{}^{\prime}}_{1}$ and the vector $U=s\Delta+\alpha$,
which will act as partial unmasker. Afterwards, $P_{1}$ sends $U$ to $P_{0}$.
The final task of $P_{2}$ is to mask the masked eigenvalues by computing
$\Lambda^{\prime\prime}=\Lambda^{\prime}\odot\Delta+\alpha$ and to send it
along with $Q^{{}^{\prime}}_{0}$ to $P_{0}$. Once $P_{0}$ receives
$Q^{{}^{\prime}}_{0}$, $U$ and $\Lambda^{\prime\prime}$, it computes
$Q_{0}=\mathcal{M}^{T}Q^{{}^{\prime}}_{0}$ and
$\Lambda^{\prime\prime\prime}=(\Lambda^{\prime\prime}-U)/\tau$. To obtain the
shares of the reciprocal of the square root of the original eigenvalues, i.e.
$\Lambda^{-1/2}$, $P_{0}$ and $P_{1}$ perform private multiplication with
inputs $((\Lambda^{\prime\prime\prime})^{-1/2},0)$ and $(0,\Delta^{1/2})$,
respectively. The final step for the proxies to obtain the share of the
inverse square root of $G$ is to reconstruct it by computing
$G_{i}^{-1/2}=Q_{i}\cdot diag(\Lambda^{-1/2})\cdot Q_{i}^{T}$ for
$i\in\\{0,1\\}$.

Note that $\mathsf{INVSQRT}$ is a private method with a high probability,
since $\mathcal{M}$, $\tau$, $s$, $\Delta$ and $\alpha$ must be in a problem
specific range not to cause overflow in the masking procedures and loss of
masked values’ connection to the real numbers. This restriction affects the
feature of complete randomness of the shares in MPC. Based on the resulting
values from the computation, the parties can reduce the range of secret
values, which in this case are the values in $G$, $Q$ and $\Lambda$, to a
range slightly smaller than the original range. However, to perform private
inference on RKNs, the inverse square root of the required Gram matrices can
be outsourced along with the other parameters of the model and the randomness
regarding the aforementioned values is completely preserved. The reason of
including $\mathsf{INVSQRT}$ here is to initiate the effort of realizing such
a computation, which is required if one wishes to train an RKN or other
algorithms requiring this operation via MPC. Main contribution of this paper,
i.e. the realization of private inference on RKNs via efficient and/or
previously unaddressed functions such as the exact exponential, is still valid
even if we use the outsourced inverse square root of the Gram matrices. We
discuss further details in Supplement.

1 Algorithm  _$\mathsf{INVSQRT}$(__)_

input : $P_{0}$ and $P_{1}$ hold the share of Gram matrix $\langle
G\rangle_{0}$ and $\langle G\rangle_{1}$, respectively, and common random
scalars $s$ and $\tau$, and a common random matrix $\mathcal{M}$

output : $P_{0}$ and $P_{1}$ get $\langle G_{invsqrt}\rangle_{0}$ and $\langle
G_{invsqrt}\rangle_{1}$, respectively

2

3 $P_{i}$ locally computes $\langle
G^{\prime}\rangle_{i}=\mathcal{M}(\tau\langle G\rangle_{i}+sI)\mathcal{M}^{T}$
and sends $\langle G^{\prime}\rangle_{i}$ to $P_{2}$ for $i\in\\{0,1\\}$

4

5 $P_{2}$ computes $G^{\prime}=\langle G^{\prime}\rangle_{0}+\langle
G^{\prime}\rangle_{1}$ to reconstruct the masked Gram matrix

6

7 $P_{2}$ computes $(\Lambda^{\prime},Q^{\prime})=localEIG(G^{\prime})$ to
locally perform eigenvalue decomposition

8

9 $P_{2}$ splits the masked eigenvectors $Q^{\prime}$ into two shares $\langle
Q^{\prime}\rangle_{0}$ and $\langle Q^{\prime}\rangle_{1}$

10

11 $P_{2}$ generates a vector of random values $\Delta$ and a random scalar
$\alpha$, and sends $\Delta$, $\alpha$ and $\langle Q^{\prime}\rangle_{1}$ to
$P_{1}$

12

13 $P_{2}$ locally computes
$\Lambda^{\prime\prime}=(\Lambda^{\prime}\odot\Delta+\alpha)$ to mask the
masked eigenvalues and sends $\Lambda^{\prime\prime}$ and $\langle
Q^{\prime}\rangle_{0}$ to $P_{0}$

14

15 $P_{1}$ computes $\langle Q\rangle_{1}=\mathcal{M}^{T}\langle
Q^{\prime}\rangle_{1}$ and $U=s\Delta+\alpha$

16

17 $P_{1}$ sends $U$ to $P_{0}$

18

19 $P_{0}$ computes $\langle Q\rangle_{0}=\mathcal{M}^{T}\langle
Q^{\prime}\rangle_{0}$ and
$\Lambda^{\prime\prime\prime}=(\Lambda^{\prime\prime}-U)/\tau$

20

21 $P_{0}$ and $P_{1}$ call $\mathsf{MUL}$ with inputs
$(((\Lambda^{\prime\prime\prime})^{-1/2},0)$ and $(0,\Delta^{1/2}))$,
respectively, to obtain $\langle\Lambda\rangle_{i}$

22

23 $P_{i}$ computes $\langle
G_{invsqrt}\rangle_{i}=\mathsf{MUL}(\mathsf{MUL}(\langle
Q\rangle_{i},diag(\langle\Lambda\rangle_{i})),\langle Q\rangle_{i})$

Algorithm 5 Inverse Square Root of a Gram Matrix ($\mathsf{INVSQRT}$)

## Privacy Preserving RKNs (ppRKN)

In this section, we explain the procedure of the privacy preserving inference
on RKNs (ppRKN) via CECILIA.

### Outsourcing

The first step for the private inference is to outsource the model parameters
from the model owner and the test sample from the data owner to the proxies.
In the outsourcing of the model parameters, which are the anchor point matrix,
the linear classifier weights and the inverse square root of the Gram matrices
of the anchor points if one wants to have the full randomness, the model owner
splits them into two arithmetic shares and sends them to the proxies in such a
way that each proxy has a single share of each parameter. To outsource the
test samples, the data owner proceeds similarly after using one-hot encoding
to convert a sequence into a vector of numbers. It divides this vector into
two shares and sends them to the proxies.

### Private Inference

After outsourcing the model parameters and the test sample, we use the
building blocks of CECILIA to realize the private inference on a pre-trained
RKN. Let $t\in\\{1,\ldots,s\\}$ be the index of the characters in the
sequence. The first step is to compute the similarity of the one-hot encoded
$t$-th character of the sequence to each character of each anchor point via
Equation 2. Such a similarity calculation involves the dot product of two
secret shared vectors, the subtraction of a plaintext scalar value from a
secret shared value, the multiplication of a plaintext scalar value by a
secret shared value and the exponential of a known base raised to the power of
a secret shared value. The output of this process corresponds to the component
$b_{j}[t]$ in Equation 3. Once the similarity computation is complete, the
proxies proceed with the private elemenent-wise product between $b_{j}[t]$ and
$c_{j-1}[t-1]$, which is the initial mapping of the sequence up to the
$(t-1)$-th character based on the anchor points of length $(j-1)$. Then the
proxies add the result of the element-wise product with the downgraded
$c_{j}[t-1]$ with a plaintext scalar value $\lambda$. At the end of this
computation, the proxies obtain the secret shared initial mapping of the
sequence $c_{j}[t]$ up to the $t$-th character to $q$-dimensional space based
on each anchor point of length $j\in\\{1,\ldots,k$}.

After computing the initial mapping of the sequence up to its total length,
that is obtaining $c_{k}[s]$, the proxies either compute the inverse square
root of the Gram matrices of the anchor points up to each length of the anchor
points via $\mathsf{INVSQRT}$, which is $K_{Z_{j}Z_{j}}^{-1/2}$ where $Z_{j}$
is the anchor points up to the $j$-th character, or use the inverse square
root of the Gram matrices directly if they are outsourced. To obtain the final
mapping, they multiply these matrices by the corresponding initial mapping
vectors of the sequence. Afterwards, the proxies perform the private dot
product computation of two secret shared vectors, which are the weights of the
classifier and the mapping of the sequence. In the end, they obtain the
prediction of the ppRKN for the given sequence in a secret shared form. These
shares can then be sent back to the owner of the data enabling the
reconstruction of the prediction of the ppRKN.

## Security Analysis

The security analyses of the methods of CECILIA and the protocol ppRKN are
given in the Supplement. For the security analysis of the adapted functions
from other studies, we kindly refer the reader to the corresponding papers.

Protocol |  RC |  Protocol |  RC
---|---|---|---
MUL |  2 |  MUX |  2
MOC |  4 |  DP |  2
MSB |  4 |  EXP |  24
CMP |  4 |  INVSQRT |  15
Table 1: Round complexities (RC) of CECILIA

## Complexity Analysis of the Framework

We analyze the communication round complexities (RC) of the methods in
CECILIA. Table 1 gives RCs required by each method in CECILIA to perform the
designated operation. Note that we also give the analysis of the methods that
we adapted from the other studies to give the comprehensive view of CECILIA.

## Results

### Dataset

To replicate the experiments of the RKN and show the correctness of ppRKN, we
utilized the same dataset as Chen, Jacob, and Mairal (2019), namely Structural
Classification of Proteins (SCOP) version 1.67 (Murzin et al. 1995). It
contains $85$ fold recognition tasks with positive or negative label, and
protein sequences of varying lengths.

We also created two synthetic test sets for the scalability analysis. One of
them contains $5$ sequences of length $128$. We used this test set to evaluate
the effects of the hyperparameters $q$ and $k$ on the execution time of
CECILIA. The second test set consists of $5$ sequences of different lengths to
evaluate the effect of sequence length on execution time.

### Experimental Setup

To conduct experiments, we used Amazon EC2 t2.xlarge instances. For the LAN
setting, we selected all instances from the Frankfurt region, and for the WAN
setting, we selected additional instances from London and Paris. We represent
numbers with $64$ bits whose least significant $20$ bits are used for the
fractional part. When computing the final mapping of the sequence, we chose to
use $\mathsf{INVSQRT}$ to compute the inverse square root of the Gram
matrices. Since the runtime of the experiments using the outsourced inverse
square root of the Gram matrices is significantly shorter, we wanted to
provide the worst case runtime results by using $\mathsf{INVSQRT}$.

### Correctness Analysis

We conducted the experiments to analyze the correctness of the ppRKN on the
LAN. We performed the predictions of test samples of a task of SCOP to verify
that we can get the same predictions as with RKN in plaintext. We selected the
first task and trained a different RKN model for each combination of the
parameters $q\in\\{16,32,64,128\\}$ and $k\in\\{5,7,10\\}$. We randomly
selected $5$ sequences from the test samples utilized by the RKN in this task
and input these sequences into ppRKN. Then, ppRKN privately performed the
prediction of these test samples on the secret shared model. At the end of
this process, we received the shares of the predictions from the proxies and
performed the reconstruction to obtain these predictions in plaintext. When we
compared the predictions of ppRKN with the predictions of plaintext RKN for
these randomly selected test samples, the largest absolute difference between
the corresponding predictions is less than $2\times 10^{-5}$. Such close
predictions of ppRKN suggest that CECILIA can correctly perform almost the
same predictions as RKN without sacrificing the privacy of the test samples or
the model.

There are two reasons for the small deviation in prediction results. The first
is the truncation operation at the end of the multiplication to preserve the
format of the resulting number (Mohassel and Zhang 2017). Although this causes
only a small error for a single multiplication, it can add up for thousands of
multiplications. Second, the utilized number format has inherently limited
precision and this could cause a small error between the representation of the
value in our format and the actual value itself. Such an error could also
accumulate during the calculation and lead to a relatively high error in the
end. Depending on the problem, one can address these problems by setting the
number of bits for the fractional part to a higher number to increase the
precision and minimize the overall error.

(a)

(b)

(c)

Figure 2: The results of the execution time analysis of CECILIA on both WAN
and LAN settings for varying (a) number of anchor points for a fixed k-mer
length and sequence length, (b) length of k-mers for a fixed number of anchor
points and sequence length, and (c) length of sequences for a fixed number of
anchor points and length of k-mer.

### Execution Time Analysis

We also examined the effects of the parameters of the RKN, namely the number
of anchor points, the length of the k-mers and the length of the sequence, on
the execution time of CECILIA on both LAN and WAN. For this purpose, we used
two synthetic datasets. In the analysis of a parameter, we vary that
parameters and fixed the rest. Figure 2 summarizes the results and
demonstrates the linear trend in the execution time of CECILIA for different
parameters.

## Conclusion

In this study, we propose a new comprehensive secure framework, called
CECILIA, based on three party computation to enable privacy preserving
computation of complex functions utilized in machine learning algorithms. Some
of these functions have not been addressed before such as the exact
exponential computation and the inverse square root of a Gram matrix and some
of them had been solved either inefficiently or using a different methodology.
We used CECILIA to realize the privacy preserving inference on a pre-trained
RKN. To the best of our knowledge, CECILIA is the first study performing such
an inference on RKNs. We demonstrate the correctness of CECILIA on the LAN
using the tasks of the original RKN paper. The experiments demonstrate that
CECILIA can give almost the exact same prediction result as one could obtain
from the RKN without privacy. This demonstrates the correctness of the
proposed methods. Moreover, we evaluated the performance of CECILIA in terms
of execution time on two synthetic datasets. We carried out these analyses on
both WAN and LAN, and showed that CECILIA scales linearly with the length of
the k-mers and the length of the input sequence, and almost linearly with the
number of anchor points. Note that CECILIA is not just limited to the privacy
preserving inference on RKNs but also capable of realizing other machine
learning algorithms requiring the provided building blocks. In future work,
building on CECILIA, we plan to address the privacy preserving training of the
RKN and other state-of-the-art machine learning algorithms. We also aim to
improve the security of $\mathsf{INVSQRT}$ so that we do not lose the complete
randomness of MPC in the computation of the inverse square root of a secret
shared Gram matrix.

## References

  * Abadi et al. (2016) Abadi, M.; Chu, A.; Goodfellow, I.; McMahan, H. B.; Mironov, I.; Talwar, K.; and Zhang, L. 2016.  Deep learning with differential privacy.  In _Proceedings of the 2016 ACM SIGSAC conference on computer and communications security_ , 308–318.
  * Applebaum (2017) Applebaum, B. 2017.  Garbled circuits as randomized encodings of functions: a primer.  In _Tutorials on the Foundations of Cryptography_ , 1–44. Springer.
  * Bakshi and Last (2020) Bakshi, M.; and Last, M. 2020.  CryptoRNN-Privacy-Preserving Recurrent Neural Networks Using Homomorphic Encryption.  In _International Symposium on Cyber Security Cryptography and Machine Learning_ , 245–253. Springer.
  * Beaver (1991) Beaver, D. 1991.  Efficient Multiparty Protocols Using Circuit Randomization.  In _Advances in Cryptology - CRYPTO ’91, 11th Annual International Cryptology Conference, Santa Barbara, California, USA, August 11-15, 1991, Proceedings_ , 420–432.
  * Chen, Jacob, and Mairal (2019) Chen, D.; Jacob, L.; and Mairal, J. 2019.  Recurrent Kernel Networks.  _arXiv preprint arXiv:1906.03200_.
  * Chen et al. (2020) Chen, S.; Fu, A.; Shen, J.; Yu, S.; Wang, H.; and Sun, H. 2020.  RNN-DP: A new differential privacy scheme base on Recurrent Neural Network for Dynamic trajectory privacy protection.  _Journal of Network and Computer Applications_ , 168: 102736.
  * Damgård et al. (2012) Damgård, I.; Pastro, V.; Smart, N.; and Zakarias, S. 2012.  Multiparty computation from somewhat homomorphic encryption.  In _Annual Cryptology Conference_ , 643–662. Springer.
  * Dwork, Roth et al. (2014) Dwork, C.; Roth, A.; et al. 2014.  The algorithmic foundations of differential privacy.  _Found. Trends Theor. Comput. Sci._ , 9(3-4): 211–407.
  * EL-Manzalawy, Dobbs, and Honavar (2008) EL-Manzalawy, Y.; Dobbs, D.; and Honavar, V. 2008.  Predicting linear B-cell epitopes using string kernels.  _Journal of Molecular Recognition: An Interdisciplinary Journal_ , 21(4): 243–255.
  * Gilad-Bachrach et al. (2016) Gilad-Bachrach, R.; Dowlin, N.; Laine, K.; Lauter, K.; Naehrig, M.; and Wernsing, J. 2016.  Cryptonets: Applying neural networks to encrypted data with high throughput and accuracy.  In _International conference on machine learning_ , 201–210. PMLR.
  * Hesamifard, Takabi, and Ghasemi (2017) Hesamifard, E.; Takabi, H.; and Ghasemi, M. 2017.  Cryptodl: Deep neural networks over encrypted data.  _arXiv preprint arXiv:1711.05189_.
  * Hochreiter, Heusel, and Obermayer (2007) Hochreiter, S.; Heusel, M.; and Obermayer, K. 2007.  Fast model-based protein homology detection without alignment.  _Bioinformatics_ , 23(14): 1728–1736.
  * Huang et al. (2022) Huang, Z.; Lu, W.-j.; Hong, C.; and Ding, J. 2022.  Cheetah: Lean and Fast Secure Two-Party Deep Neural Network Inference.  _IACR Cryptol. ePrint Arch._ , 2022: 207.
  * Knott et al. (2021) Knott, B.; Venkataraman, S.; Hannun, A.; Sengupta, S.; Ibrahim, M.; and van der Maaten, L. 2021.  CrypTen: Secure multi-party computation meets machine learning.  _arXiv preprint arXiv:2109.00984_.
  * Leslie, Eskin, and Noble (2001) Leslie, C.; Eskin, E.; and Noble, W. S. 2001.  The spectrum kernel: A string kernel for SVM protein classification.  In _Biocomputing 2002_ , 564–575. World Scientific.
  * Lu et al. (2021) Lu, W.-j.; Huang, Z.; Hong, C.; Ma, Y.; and Qu, H. 2021.  PEGASUS: bridging polynomial and non-polynomial evaluations in homomorphic encryption.  In _2021 IEEE Symposium on Security and Privacy (SP)_ , 1057–1073. IEEE.
  * Mishra et al. (2020) Mishra, P.; Lehmkuhl, R.; Srinivasan, A.; Zheng, W.; and Popa, R. A. 2020.  Delphi: A cryptographic inference service for neural networks.  In _29th USENIX Security Symposium (USENIX Security 20)_ , 2505–2522.
  * Mohassel and Zhang (2017) Mohassel, P.; and Zhang, Y. 2017.  SecureML: A System for Scalable Privacy-Preserving Machine Learning.  In _2017 IEEE symposium on security and privacy (SP)_ , 19–38. IEEE.
  * Murzin et al. (1995) Murzin, A. G.; Brenner, S. E.; Hubbard, T.; and Chothia, C. 1995.  SCOP: a structural classification of proteins database for the investigation of sequences and structures.  _Journal of molecular biology_ , 247(4): 536–540.
  * Nojoomi and Koehl (2017) Nojoomi, S.; and Koehl, P. 2017.  A weighted string kernel for protein fold recognition.  _BMC bioinformatics_ , 18(1): 1–14.
  * Patra et al. (2021) Patra, A.; Schneider, T.; Suresh, A.; and Yalame, H. 2021.  $\\{$ABY2. 0$\\}$: Improved $\\{$Mixed-Protocol$\\}$ Secure $\\{$Two-Party$\\}$ Computation.  In _30th USENIX Security Symposium (USENIX Security 21)_ , 2165–2182.
  * Rathee et al. (2020) Rathee, D.; Rathee, M.; Kumar, N.; Chandran, N.; Gupta, D.; Rastogi, A.; and Sharma, R. 2020.  CrypTFlow2: Practical 2-party secure inference.  In _Proceedings of the 2020 ACM SIGSAC Conference on Computer and Communications Security_ , 325–342.
  * Wagh, Gupta, and Chandran (2019) Wagh, S.; Gupta, D.; and Chandran, N. 2019.  SecureNN: 3-Party Secure Computation for Neural Network Training.  _Proc. Priv. Enhancing Technol._ , 2019(3): 26–49.
  * Wagh et al. (2020) Wagh, S.; Tople, S.; Benhamouda, F.; Kushilevitz, E.; Mittal, P.; and Rabin, T. 2020\.  Falcon: Honest-majority maliciously secure framework for private deep learning.  _arXiv preprint arXiv:2004.02229_.
  * Zhou and Li (2016) Zhou, L.; and Li, C. 2016.  Outsourcing Eigen-Decomposition and Singular Value Decomposition of Large Matrix to a Public Cloud.  _IEEE Access_ , 4: 869–879.

## Supplement

### Number Format

To illustrate the number format, let $n$ be the number of bits to represent
numbers, $f$ be the number of bits allocated for the fractional part and
$\mathbb{S}$ be the set of values that can be represented in this number
format, one can convert $x\in\mathbb{R}$ to $\hat{x}\in\mathbb{S}$ as follows:

$\hat{x}=\left\\{\begin{array}[]{ll}\lfloor x*2^{f}\rfloor&\quad x\geq 0\\\
2^{n}-\lfloor|x*2^{f}|\rfloor&\quad x<0\end{array}\right.$ (5)

For example, $x=3.42$ is represented as $\hat{x}=112066$, which, if you omit
the leading zeros, is $1101101011100001010001$ in binary. The two most
significant bits, the $21$st and $22$nd bits, are used to represent $3$ and
the remainder represents $0.42$.

The choice of $f$ depends on the task for which CECILIA is used. If the
numbers in the task are mostly small or they need to be as accurate as
possible, the precision of the representation of the numbers is crucial. In
such a case, a higher value for $f$ is required to allow for more decimal
places of the original value, which in turn requires sacrifices in the upper
limit that can be calculated by the building blocks of CECILIA. However, if
the numbers that appear during the process are large, $f$ must be set to lower
values to allow the integer part to represent larger numbers with fewer
decimal places.

In our experiments with CECILIA, we represent the numbers with $64$ bits and
set $f$ to $20$ to have a relatively high precision in representing the
numbers and a sufficient range for the exponential. Although we have used the
$64$ bit representation in the explanations of CECILIA’s functions, CECILIA
can also work with a different number of bit representations. The reason for
our choice is to take advantage of the natural modular operation of most
modern CPUs, which operate with $64$ bits. In this way, we avoid performing
modular operations after each arithmetic operation performed on the shares.
However, one can still use a larger or smaller number of bits on $64$ bit CPUs
by performing the modular operations manually.

### Discussion on MSB

We give the detailed pseudocode of $\mathsf{MSB}$ in Algorithm 2 in the main
paper. Since it is optimized to reduce the communication round by integrating
$\mathsf{MOC}$, in this section we give the compact version of the pseudocode
of it in Algorithm 6 in the Supplement. We also give the security analysis of
$\mathsf{MSB}$ based on this compact version to make it easy to follow.

1 Algorithm  _$\mathsf{MSB}$(__)_

input : $P_{0}$ and $P_{1}$ hold ($\langle x\rangle_{0}$) and ($\langle
x\rangle_{1}$), respectively

output : $P_{0}$ and $P_{1}$ get $\langle z\rangle_{0}$ and $\langle
z\rangle_{1}$, respectively, where $z$ is equal to $0$ if the most significant
bit of $x$ is $0$ and $1$ otherwise.

2 $P_{0}$ and $P_{1}$ hold a common random bit $f$

3 For each $i\in\\{0,1\\}$, $P_{i}$ executes Steps $4$-$9$.

4 $\langle d\rangle_{i}^{K}=\langle x\rangle_{i}\mod K$

5 $\langle d\rangle_{i}=\mathsf{MOC}(\langle d\rangle_{i}^{K})$.

6 $\langle z\rangle_{i}=\langle x\rangle_{i}-\langle d\rangle_{i}$.

7 $\langle a[0]\rangle_{i}=ifK-\langle z\rangle_{i}$

8 $\langle a[1]\rangle_{i}=i(1-f)K-\langle z\rangle_{i}$

9 $P_{i}$ sends $\langle a\rangle_{i}$ to $P_{2}$

10 $P_{2}$ reconstructs $a[j]$ where $j\in\\{0,1\\}$ and computes
$a[j]=a[j]/K$

11 $P_{2}$ creates two fresh shares of $a[j]$ where $j\in\\{0,1\\}$ sends them
to $P_{0}$ and $P_{1}$

12 For each $i\in\\{0,1\\}$, $P_{i}$ executes Step $13$

13 $\langle z\rangle_{i}=\langle a[f]\rangle_{i}$

14

Algorithm 6 Compact and unoptimized version of Most Significant Bit
($\mathsf{MSB}$)

###  $\mathsf{MUL}$ and Details of Truncation

Let $i\in\\{0,1\\}$ be and $\langle c\rangle_{i}=\langle
a\rangle_{i}\cdot\langle b\rangle_{i}$ be the shares of the multiplication
triple and $\langle x\rangle_{i}$ and $\langle y\rangle_{i}$ be the shares of
$x$ and $y$, respectively, in $P_{i}$. To compute $\langle z\rangle=\langle
x\rangle\cdot\langle y\rangle$, $P_{i}$ first computes $\langle
e\rangle_{i}=\langle x\rangle_{i}-\langle a\rangle_{i}$ and $\langle
f\rangle_{i}=\langle y\rangle_{i}-\langle b\rangle_{i}$, and then sends them
to $P_{1-i}$. Afterwards, $P_{i}$ reconstructs $e$ and $f$ and calculates
$\langle z\rangle_{i}=i\cdot e\cdot f+f\cdot\langle a\rangle_{i}+e\cdot\langle
b\rangle_{i}+\langle c\rangle_{i}$. In the last step to obtain the shares of
the multiplication, $P_{0}$ shifts $z_{0}$ to the right by $f$. Slightly
differently, $P_{1}$ first shifts $(-1*z_{1})$ to the right by $f$ and then
multiplies the result by $-1$.

### Details of Security Analysis

In this part, we give the security analysis of $\mathsf{MOC}$, $\mathsf{MSB}$,
$\mathsf{CMP}$, $\mathsf{MUX}$, $\mathsf{EXP}$ and $\mathsf{INVSQRT}$. In the
security analysis of some of these methods, we take the security of the
private compare functionality $\mathcal{F}_{\mathsf{PC}}$ (Wagh, Gupta, and
Chandran 2019) as base. We kindly refer the readers to Wagh, Gupta, and
Chandran (2019) for the security analysis of $\mathcal{F}_{\mathsf{PC}}$.

######  Lemma 1.

The protocol $\mathsf{MOC}$ in Algorithm 1 in the main paper securely realizes
the functionality $\mathcal{F}_{MOC}$ in $\mathcal{F}_{\mathsf{PC}}$ hybrid
model.

###### Proof.

First, we prove the correctness of our protocol by showing $(\langle
x\rangle_{0}^{K}+\langle x\rangle_{1}^{K})\mod K=(\langle x\rangle_{0}+\langle
x\rangle_{1})\mod L$. In the protocol, $y=(x+r)\mod K$ and
$\mathsf{isWrap}(x,r,K)=r\overset{?}{>}y$, that is, $\mathsf{isWrap}(x,r,K)=1$
if $r>y$, $0$ otherwise. At the beginning, $P_{0},P_{1}$ and $P_{2}$ call
$\mathcal{F}_{PC}$ to compute $c=r\overset{?}{>}y$ and $P_{0}$ and $P_{1}$
obtain the boolean shares $c_{0}$ and $c_{1}$, respectively. Besides, $P_{2}$
sends also the boolean shares $w_{0}$ and $w_{1}$ of
$w=\mathsf{isWrap}(\langle r\rangle_{0},\langle r\rangle_{1},K)$ to $P_{0}$
and $P_{1}$, respectively. If $\mathsf{isWrap}(\langle y\rangle_{0},\langle
y\rangle_{1},K)$ is $1$ then $P_{0}$ adds $K$ to $\langle y\rangle_{0}$ to
change the ring of $y$ from $K$ to $L$. To convert $r$ from ring $K$ to ring
$L$, $P_{0}$ and $P_{1}$ add $K$ to their shares of $r$ based on their boolean
shares $w_{0}$ and $w_{1}$, respectively. If $w_{0}=1$, then $P_{0}$ adds $K$
to its $r_{1}$ and $P_{1}$ does the similar with its shares. Later, we need to
fix is the summation of $x$ and $r$, that is the value $y$. In case of
$x+r\geq K$, we cannot fix the summation value $y$ in ring $L$ by simply
converting it from ring $K$ to ring $L$. This summation should be $x+r$ in
ring $L$ rather than $(x+r)$ $mod$ $K$. To handle this problem, $P_{0}$ and
$P_{1}$ add $K$ to their shares of $y$ based on their shares $c_{0}$ and
$c_{1}$. As a result, we convert the values $y$ and $r$ to ring $L$ and fix
the value of $y$ if necessary. The final step to obtain $x_{i}$ for party
$P_{i}$ is simply subtract $r_{i}$ from $y_{i}$ where $i\in\\{0,1\\}$.

Next, we prove the security of our protocol. $P_{2}$ involves this protocol in
execution of $\mathcal{F}_{PC}$. We give the proof $\mathcal{F}_{PC}$ above.
At the end of the execution of $\mathcal{F}_{PC}$, $P_{2}$ learns
$u^{\prime}$. However, $u^{\prime}=u\oplus(x>r)$ and $P_{2}$ does not know
$u$. Thus $u^{\prime}$ is uniformly distributed and can be perfectly simulated
with randomly generated values. $P_{i}$ where $i\in\\{0,1\\}$ sees fresh
shares of $\langle r\rangle_{i}^{K}$, $\\{\langle
r[j]\rangle_{i}^{V}\\}_{j\in[n]}$, $w_{i}^{B}$ and $u_{i}^{B}$. These values
can be perfectly simulated with randomly generated values. ∎

######  Lemma 2.

The protocol $\mathsf{MSB}$ in Algorithm 2 in the main paper securely realizes
the functionality $\mathcal{F}_{MSB}$ in $\mathcal{F}_{\mathsf{MOC}}$ hybrid
model.

###### Proof.

First, we prove the correctness of our protocol. Assume that we have $n$-bit
number $u$. $v=u-(u\mod 2^{n-1})$ is either $0$ or $2^{n-1}$. In our protocol,
$\langle z\rangle_{i}$ is the output of $P_{i}$ where $i\in\\{0,1\\}$. We have
to prove that $\mathsf{Reconstruct(\langle z\rangle_{i})}$ is equal to the
most significant bit (MSB) of $x$. $P_{i}$ where $i\in\\{0,1\\}$ computes
$d_{i}^{K}=x_{i}\mod K$ which is a share of $d$ over $K$. $P_{i}$ computes
$d_{i}$ which is a share of $d$ over $L$ by invoking $\mathsf{MOC}$. Note that
$z=x-\mathsf{Reconstruct(\langle d\rangle_{i})}$ and all bits of $z$ are $0$
except the MSB of $z$, which is equal to the MSB of $x$. Now, we have to map
$z$ to $1$ if it is equal to $K$ or $0$ if it is equal to $0$. $P_{0}$ sends
the $z_{0}$ and $z_{0}+K$ in random order to $P_{2}$ and $P_{1}$ sends the
$z_{1}$ to $P_{2}$. $P_{2}$ reconstructs two different values, divides these
values by $K$, creates two additive shares of them, and sends these shares to
$P_{0}$ and $P_{1}$. Since $P_{0}$ and $P_{1}$ know the order of the real MSB
value, they correctly select the shares of its mapped value.

Second, we prove the security of our protocol. $P_{i}$ where $i\in\\{0,1\\}$
sees $\langle d\rangle_{i}$, which is a fresh share of $d$, and $\langle
a[0]\rangle_{i}$ and $\langle a[1]\rangle_{i}$, one of which is a fresh share
of the MSB of $x$ and the other is a fresh share of the complement of the MSB
of $x$. Thus the view of $P_{i}$ can be perfectly simulated with randomly
generated values. ∎

######  Lemma 3.

The protocol $\mathsf{CMP}$ in the main paper securely realizes the
functionality $\mathcal{F}_{CMP}$ in $\mathcal{F}_{\mathsf{MSB}}$ hybrid
model.

###### Proof.

First, we prove the correctness of our protocol. Assume that we have $x$ and
$y$. We first compute $z=x-y$. If $z$ is negative, which corresponds to $1$ in
the most significant bit of $z$, it means that $y>x$. In this case,
$\mathsf{CMP}$ outputs $1$. If $z$ is non-negative, which corresponds to $0$
in the most significant bit of $z$, then it indicates that $x\geq y$. In this
case, the output of $\mathsf{CMP}$ is $0$. Since the output of $\mathsf{CMP}$
exactly matches the output of $\mathsf{MSB}$ and we have already proved the
correctness of $\mathsf{MSB}$, we can conclude that $\mathsf{CMP}$ works
correctly.

Second, we prove the security of our protocol. Since $\langle
z\rangle_{i}=\langle x\rangle_{i}-\langle y\rangle_{i}$ is computed locally by
$P_{i}$ for $i\in\\{0,1\\}$, it does not reveal any information about $x$ and
$y$. Afterward, $\mathsf{MSB}$ is called on $\langle z\rangle_{i}$ to
determine the most significant bit of $z$ in secret shared form. Considering
that the security of $\mathsf{MSB}$ is proven, we can conclude that
$\mathsf{CMP}$ compares two secret shared values without compromising their
privacy. ∎

######  Lemma 4.

The protocol $\mathsf{MUX}$ in Algorithm 3 in the main paper securely realizes
the functionality $\mathcal{F}_{\mathsf{MUX}}$.

###### Proof.

We first prove the correctness of our protocol. $\langle z\rangle_{i}$ is the
output of $P_{i}$ where $i\in\\{0,1\\}$. We need to prove that
$\mathsf{Reconstruct(\langle z\rangle_{i})}=(1-b)x+by$.

$\begin{split}\langle z\rangle_{0}+\langle z\rangle_{1}={}&\langle
x\rangle_{0}-\langle b\rangle_{0}(\langle x\rangle_{0}-\langle
y\rangle_{0})+r_{1}\langle b\rangle_{0}+r_{2}(\langle x\rangle_{0}-\langle
y\rangle_{0})\\\ &+r_{2}r_{3}+\langle x\rangle_{1}-\langle
b\rangle_{1}(\langle x\rangle_{1}-\langle y\rangle_{1})+r_{0}(\langle
x\rangle_{1}-\langle y\rangle_{1})\\\ &+r_{0}r_{1}+r_{3}\langle
b\rangle_{1}-\langle b\rangle_{0}\langle x\rangle_{1}+\langle
b\rangle_{0}\langle y\rangle_{1}-\langle b\rangle_{0}r_{1}-r_{0}\langle
x\rangle_{1}\\\ &+r_{0}\langle y\rangle_{1}-r_{0}r_{1}-\langle
x\rangle_{0}\langle b\rangle_{1}-\langle x\rangle_{0}r_{2}+\langle
y\rangle_{0}\langle b\rangle_{1}\\\ &+\langle y\rangle_{0}r_{2}-r_{3}\langle
b\rangle_{1}-r_{3}r_{2}\\\ ={}&(1-\langle b\rangle_{0}-\langle
b\rangle_{1})(\langle x\rangle_{0}+\langle x\rangle_{1})+(\langle
b\rangle_{0}+\langle b\rangle_{1})(\langle y\rangle_{0}+\langle
y\rangle_{1})\\\ ={}&(1-b)x+by\end{split}$ (6)

Next we prove the security of our protocol. $P_{2}$ gets $M_{2},M_{3},M_{5}$
and $M_{6}$. All these values are uniformly random values because they are
generated using uniformly random values $r_{0},r_{1},r_{2},r_{4}$. $P_{2}$
computes $M_{2}M_{5}+M_{3}M_{6}$. The computed value is still uniformly random
because it contains uniformly random values $r_{0},r_{1},r_{2},r_{4}$. As a
result, any value learned by $P_{2}$ is perfectly simulated. For each
$i\in\\{0,1\\}$, $P_{i}$ learns a fresh share of the output. Thus $P_{i}$
cannot associate the share of the output with the shares of the inputs and any
value learned by $P_{i}$ is perfectly simulatable. ∎

######  Lemma 5.

The protocol $\mathsf{EXP}$ securely computes the exponential of a publicly
known base raised to the power of a secret shared value.

###### Proof.

We begin the proof by showing the correctness of the method. Let $x$ be the
power whose representation in our number format is $\langle x\rangle$ and $b$
be the publicly known base. $P_{0}$ or $P_{1}$ computes
$C_{p}=\\{\ldots,b^{8},b^{4},b^{2},b,b^{1/2},b^{1/4},b^{1/8},\ldots\\}$ and
$C_{n}=\\{\ldots,b^{-8},b^{-4},b^{-2},b^{-1},b^{-1/2},b^{-1/4},b^{-1/8},\ldots\\}$
and the other generates a corresponding set of $0$s for $C_{p}$ and $C_{n}$.
These values in $C_{p}$ and $C_{n}$ correspond to $b^{2^{i}}$ and $b^{-1\cdot
2^{i}}$, respectively, for $i\in\\{(n-f),\ldots,2,1,0,-1,-2,\ldots,-f\\}$
assuming that only the corresponding bit value of the power $x$ is $1$. They
choose one of these sets based on the sign of $x$ and let $C$ be the selected
set. Afterwards, they must choose between $c_{j}\in C$ and $1$ depending on
$\langle x\rangle_{j}$ where $j\in\\{0,1,\ldots,n\\}$. For this selection,
they use the most significant bit operation on all cases where each bit of $x$
is at the most significant bit position. This is done by shifting the shares
of $x$ to the left. Once they have the correct set of contributions, they
basically multiply all of those contributions to obtain the result of the
exponential. This proves the correctness of $\mathsf{EXP}$.

Corruption of a proxy: At the beginning, since the adversary corrupting a
proxy knows only one share of the power $x$, that is either $x_{0}$ or
$x_{1}$, it cannot infer any information about the other share. The first step
of the exponential is to compute the possible contribution of every bit of
positive and negative power. This is publicly known. The following step is to
select between these contribution depending on the result of $\mathsf{MSB}(x)$
by using $\mathsf{MUX}$. Since both $\mathsf{MSB}$ and $\mathsf{MUX}$ are
secure, the adversary can neither infer anything about $x_{1-j}$ nor relate
the share of the result it obtains to $x$ in general. In the next step, they
obtain each bit of $x$ in secret shared form by using $\mathsf{MSB}$ and bit
shifting on the shares of $x$. Considering the proven security of
$\mathsf{MSB}$ and the shifting being simply local multiplication of each
share by $2$, there is no information that the adversary could obtain.
Afterwards, the proxies select the correct contributions by employing
$\mathsf{MUX}$. Since $\mathsf{MUX}$ gives the fresh share of what is
selected, the adversary cannot associate the inputs to the output. The last
step is to multiply these selected contributions via $\mathsf{MUL}$, which is
also proven to be secure. Therefore, we can conclude that $\mathsf{EXP}$ is
secure against a semi-honest adversary corrupting a proxy.

Corruption of the helper: Since the task of the helper party in the
computation of the exponential of a secret shared power is either to provide
multiplication triples or to perform the required computation on the masked
data, there is nothing that the adversary corrupting the helper party could
learn about $x$. Therefore, it is fair to state that $\mathsf{EXP}$ is secure
against a semi-honest adversary corrupting the helper. ∎

######  Lemma 6.

The protocol $\mathsf{INVSQRT}$ privately computes the inverse square root of
a secret shared Gram matrix with a high probability.

###### Proof.

We first demonstrate the correctness of $\mathsf{INVSQRT}$. Once the Gram
matrix $G$ whose eigenvalues and eigenvectors are $\Lambda$ and $Q$,
respectively, is masked by $G^{{}^{\prime}}=\mathcal{M}(\tau
G_{i}+sI)\mathcal{M}^{T}$, we conduct eigenvalue decomposition on
$G^{{}^{\prime}}$ resulting in the masked eigenvalues and eigenvectors of $G$,
which are $\Lambda^{{}^{\prime}}=\tau\Lambda+s$ and
$Q^{{}^{\prime}}=\mathcal{M}Q$, respectively. The recovery of $Q$ can be
achieved by $Q=\mathcal{M}^{T}Q^{{}^{\prime}}$ since $M^{T}=M^{-1}$ due to the
orthogonality. Since the masked eigenvalues $\Lambda^{{}^{\prime}}$ are, once
again, masked by computing
$\Lambda^{{}^{\prime\prime}}=\Delta\Lambda^{{}^{\prime}}+\alpha$, we have to
compute the partial unmasker $U=s\Delta+\alpha$ by using $\Delta$, $\alpha$
and $s$, and subtract it from $\Lambda^{{}^{\prime\prime}}$ and divide the
result by $\tau$ to obtain $\Lambda^{{}^{\prime\prime\prime}}=\Delta\Lambda$.
Then, the multiplication of $(\Lambda^{{}^{\prime\prime\prime}})^{-1/2}$ and
$\Delta^{1/2}$ gives us the reciprocal of the square root of the original
eigenvalues, i.e. $\Lambda^{-1/2}$. Finally, we construct the inverse square
root of $G$ by computing $Q\cdot diag(\Lambda^{-1/2})\cdot Q^{T}$. This
concludes the correctness of $\mathsf{INVSQRT}$.

Corruption of a proxy: In case of a corruption of $P_{1}$ by a semi-honest
adversary, there is no gain for the adversary, since it only receives the
share of the masked eigenvectors and the masks utilized by the helper party to
mask the eigenvalues. At the end, $P_{1}$ gets the share of the reciprocal of
the square root of the eigenvalues as a secret shared result of
$\mathsf{MUL}$. Therefore, $\mathsf{INVSQRT}$ is secure against an adversary
corrupting $P_{1}$. On the other hand, if an adversary corrupts $P_{0}$, then
it can reduce the possible eigenvalues of the Gram matrix into a specific set
considering that $\mathcal{M}$, $\tau$, $s$, $\Delta$ and $\alpha$ are from a
specific range. Even in the worst case that the adversary deduces all the
eigenvalues, which is extremely unlikely and negligible, it cannot reconstruct
and learn the Gram matrix since $P_{0}$ receives only the share of the
eigenvectors, which is indistinguishable from a random matrix for $P_{0}$.
Thus, $\mathsf{INVSQRT}$ is also secure against an adversary corrupting
$P_{0}$.

Corruption of the helper: Since the values in $\mathcal{M}$, $\tau$ and $s$
are from a specific range, the helper party can utilize the linear systems
$\Lambda^{{}^{\prime}}=\tau\Lambda+s$ and $Q^{\prime}=MQ$ to narrow down the
possible values of $Q$ and $\Lambda$, and $G$ implicitly. Unless a masked
value in the linear systems is at one of its edge cases, it is not possible
for the adversary to infer the actual value. It is plausible to assume that
not every value is at the edge case in these linear systems. Furthermore,
considering that obtaining the eigenvalues of $G$ has a relatively high
chance, which is still extremely low, compared to obtaining the eigenvectors
of $G$, the adversary cannot form a linear system based on the inferred
eigenvalues resulting in a unique solution for the eigenvectors, eventually a
Gram matrix. Considering all these infeasibilities of obtaining the
eigenvalues and eigenvectors in addition to the computational infeasibility
and lack of validation method of a brute-force approach, the adversary cannot
obtain $G$ by utilizing $\Lambda^{\prime}$ and $Q^{\prime}$. Therefore,
$\mathsf{INVSQRT}$ is secure against a semi-honest adversary corrupting the
helper party with a high probability. ∎

######  Lemma 7.

In case of the utilization of the outsourced inverse square root of the Gram
matrices of the anchor points, the protocol ppRKN securely realizes the
inference of a test sample via the pre-trained RKN model. If
$\mathsf{INVSQRT}$ is utilized, then the protocol ppRKN privately performs the
inference of a test sample via the pre-trained RKN model with a high
probability.

###### Proof.

In addition to the experimental proof of the correctness of ppRKN, based on
the hybrid model, we can state that ppRKN correctly performs the private
inference since it utilizes the methods of CECILIA which are already shown to
perform their corresponding tasks correctly.

Corruption of a proxy by a semi-honest adversary: In case of utilization of
the outsourced inverse square root of the Gram matrices, the private inference
of RKN via ppRKN involves a series of $\mathsf{MUL}$, $\mathsf{ADD}$,
$\mathsf{EXP}$, $\mathsf{MSB}$, $\mathsf{MUX}$ and $\mathsf{CMP}$ operations
on secret shared data, and addition and multiplication between secret shared
data and a public scalar. Since we proved or gave the references to the proof
of the security of those methods, we can conclude that ppRKN is secure via the
hybrid model. Instead of using the outsourced inverse square root of the Gram
matrices, if we employ $\mathsf{INVSQRT}$ to compute those matrices, the
adversary corrupting $P_{0}$ can narrow down the set of possible eigenvalues
of the matrices to a relatively small set of values. Even in the worst case,
the adversary cannot obtain the Gram matrix since it has only the share of the
eigenvectors of the Gram matrix. For an adversary compromising $P_{1}$, there
is no information gain at all. Therefore, ppRKN performs the fully private
inference in case of the utilization of the outsourced inverse square root of
the Gram matrices and realizes the private inference with a high probability
when $\mathsf{INVSQRT}$ is used.

Corruption of the helper party by a semi-honest adversary: In case of
utilization of the outsourced inverse square root of the Gram matrices,
considering the same reason we stated in the scenario where the adversary
corrupts a proxy, we can conclude that ppRKN is secure via the hybrid model.
Instead of using the outsourced inverse square root of the Gram matrices, if
we employ $\mathsf{INVSQRT}$ to compute those matrices, the adversary
corrupting the helper party can narrow down the set of possible values of the
eigenvalues and the eigenvectors of the Gram matrix to a relatively small set
of values. As we discuss in the security analysis of $\mathsf{INVSQRT}$, such
a reduction on the set is not enough for the adversary to reconstruct the Gram
matrix with a high probability. Therefore, ppRKN is either fully secure or
secure with a high probability depending on whether $\mathsf{INVSQRT}$ is used
or not, respectively. ∎

### Discussion on Range of Random Values in $\mathsf{INVSQRT}$

As we discussed in the main paper, $\mathsf{INVSQRT}$ requires the utilized
random values to be from a specific range so that the masked values do not
overflow and lose its connection to the real numbers. Let
$\mathbb{E}_{\tau},\mathbb{E}_{s},\mathbb{E}_{\mathcal{M}},\mathbb{E}_{\Delta}$
and $\mathbb{E}_{\alpha}$ be the set of values from which we can select
$\tau,s,\mathcal{M}_{ij},\Delta$ and $\alpha$ in a certain setting,
respectively. And, let $\mathbb{E}_{\Lambda}$ be the set of eigenvalues of
$G$. Let us also assume that all these sets contain positive values. Then, we
have the following constraints on these random values:

  * •

$\mathbf{\Lambda}$: Since the anchor points used in the private inference are
$L_{2}$ normalized vectors, the dot product of two $L_{2}$ normalized vectors
cannot exceed $1$. Furthermore, the diagonal of the Gram matrix of these
anchor points are all $1$. Therefore, considering the fact that
$\sum_{i}^{q}\Lambda_{i}=Tr(G)$, we can confidently have the following
constraint:

$max(\mathbb{E}_{\Lambda})\leq q$ (7)

  * •

$\mathbf{\mathcal{M}}$: In the masking of $G$, we compute
$G^{{}^{\prime}}=\mathcal{M}(\tau G_{i}+sI)\mathcal{M}^{T}$. Considering the
boundary case where every entry of $G$ is $1$, we can rewrite the computation
of $G_{kl}$ as follows:

$G^{{}^{\prime}}_{kl}=\sum_{u}^{q}\sum_{v}^{q}\mathcal{M}_{kv}\big{(}\delta(u,v)(\tau+s)+(1-\delta(u,v)\tau\big{)}\mathcal{M}_{lu}$

In this computation, we have $q^{2}$-many addition operations each of which
contains $2$ multiplication operations. Assume that $\mathcal{M}_{ij}$ is the
largest entry of $\mathcal{M}$, that is
$\mathcal{M}_{ij}=max(\mathbb{E}_{\mathcal{M}})$. At some point, we have to
compute $z=\mathcal{M}_{ij}(\tau+s)\mathcal{M}_{ij}$. Due to the
multiplication and the number format that we used, $z$ has to be smaller than
$2^{n-2f-1}$. Since $G^{{}^{\prime}}_{kl}\leq 2^{n-f-1}$ not to lose the
connection to the real numbers, in the worst case that every entry of
$\mathcal{M}$ equals to $\mathcal{M}_{ij}$, we have the constraint
$q^{2}(\mathcal{M}_{ij}(\tau+s)\mathcal{M}_{ij})<2^{n-f-1}$. Considering all
these constraints, we have the following limitation on the maximum entry of
$\mathcal{M}$:

$max(\mathbb{E}_{\mathcal{M}})^{2}<\begin{cases}\frac{2^{n-2f-1}}{\tau+s}\vspace{0.2cm}\\\
\frac{2^{n-f-1}}{q^{2}(\tau+s)}\\\ \end{cases}$ (8)

The lowest of these two constraints determines the boundary of the values of
$\mathcal{M}$’s entries.

  * •

$\mathbf{\Delta}$: In $\mathsf{INVSQRT}$, $P_{2}$ masks the masked eigenvalues
by computing
$\Lambda^{{}^{\prime\prime}}=\Delta\odot\Lambda^{{}^{\prime}}+\alpha$ or
$\Lambda^{{}^{\prime\prime}}=\Delta\odot(\tau\Lambda+s)+\alpha$ if we write
$\Lambda^{{}^{\prime}}$ more explicitly. In the computation of
$\Lambda^{{}^{\prime\prime}}$, the result of an elementwise multiplication
cannot be larger than $2^{n-2f-1}$ and the overall result cannot exceed
$2^{n-f-1}$. For the range of $\Delta_{i}$ where $\Delta_{i}$ is the maximum
element of $\Delta$, we have the following piecewise limitation:

$max(\mathbb{E}_{\Delta})<\begin{cases}\frac{2^{n-2f-1}}{\tau
q+s}\vspace{0.2cm}\\\ \frac{2^{n-f-1}-\alpha}{\tau q+s}\\\ \end{cases}$ (9)

  * •

$\mathbf{\tau}$: For the upper bound on $\tau$, let us rewrite Equations 8 and
9:

$max(\mathbb{E}_{\tau})<\begin{cases}\frac{2^{n-2f-1}}{\mathcal{M}_{ij}^{2}}-s\vspace{0.2cm}\\\
\frac{2^{n-f-1}}{q^{2}\mathcal{M}_{ij}^{2}}-s\vspace{0.2cm}\\\
\frac{2^{n-2f-1}-s\Delta}{\Delta q}\vspace{0.2cm}\\\
\frac{2^{n-f-1}-\alpha-s\Delta}{\Delta q}\\\ \end{cases}$ (10)

  * •

$\mathbf{s}$: In the initial masking of the eigenvalues of $G$, we compute
$\tau G_{i}+sI$ ignoring the orthogonal matrix multiplications. This leads to
the masking $\Lambda^{{}^{\prime}}=\tau\Lambda+s$. In this masking, we have
limited range for the values. Due to its nature of such a limitation, we have
edge cases which one can uniquely map to corresponding $\tau$, $\Lambda$ and
$s$ values. In order to minimize this, $s$ has to be large enough to cover the
total gap between the largest and the third largest possible result of the
multiplication between the possible values of $\tau$ and $\Lambda$. Note that
the largest possible gap is
$min(max(\mathbb{E}_{\tau}),max(\mathbb{E}_{\Lambda}))$ and it is mostly fair
to assume that $max(\mathbb{E}_{\tau})>max(\mathbb{E}_{\Lambda})$. Then, we
can approximate the lower bound of $min(\mathbb{E}_{s})$ as follows:

$min(\mathbb{E}_{s})>2\cdot max(\mathbb{E}_{\Lambda})$

When we set $max(\mathbb{E}_{\Lambda})$ to $q$ for the worst case scenario, we
end up with the following constraint:

$min(\mathbb{E}_{s})>2q$ (11)

With the assumption that $\mathbb{E}_{\tau}$ and $\mathbb{E}_{\Lambda}$
contains positive values, this lower bound on $s$ reduces the number of cases
that can be uniquely mapped to
$min(max(\mathbb{E}_{\tau}),max(\mathbb{E}_{\Lambda}))+1$. For the upper bound
on $s$, let us rewrite the constraints given in Equation 10, we obtain the
following upper bound on $s$:

$max(\mathbb{E}_{s})<\begin{cases}\frac{2^{n-2f-1}}{\mathcal{M}_{ij}^{2}}-\tau\vspace{0.2cm}\\\
\frac{2^{n-f-1}}{q^{2}\mathcal{M}_{ij}^{2}}-\tau\vspace{0.2cm}\\\
\frac{2^{n-2f-1}}{\Delta}-\tau q\vspace{0.2cm}\\\
\frac{2^{n-f-1}-\alpha}{\Delta}-\tau q\\\ \end{cases}$ (12)

  * •

$\mathbf{\alpha}$: Considering the masking that we described in the discussion
of $\Delta$ and the argument on the choice of lower bound of $s$, we have the
following lower and upper bounds for $\alpha$:

$\begin{split}min(\mathbb{E}_{\alpha})>2\cdot
min\big{(}max(\mathbb{E}_{\Delta}),&max(\mathbb{E}_{\tau})\cdot
max(\mathbb{E}_{\Lambda})\\\ &+max(\mathbb{E}_{s})\big{)}\end{split}$ (13)
$max(\mathbb{E}_{\alpha})<2^{n-f-1}-2^{n-2f-1}$ (14)

In order to demonstrate the resulting ranges, let us image a scenario where we
set the upper limits of the random values as
$max(\mathbb{E}_{\Lambda})=2^{5}$, $max(\mathbb{E}_{\tau}=2^{5})$,
$max(\mathbb{E}_{s}=2^{10})$, $max(\mathbb{E}_{\mathcal{M}}=2^{5})$,
$max(\mathbb{E}_{\Delta}=2^{10})$ and $max(\mathbb{E}_{\alpha}=2^{30})$. Let
us also set $f=20$ and $n=64$. Considering the number format that we use, the
size of these set are $max(\mathbb{E}_{\Lambda})=2^{25}$,
$max(\mathbb{E}_{\tau}=2^{25})$, $max(\mathbb{E}_{s}=2^{30})$,
$max(\mathbb{E}_{\mathcal{M}}=2^{25})$, $max(\mathbb{E}_{\Delta}=2^{30})$ and
$max(\mathbb{E}_{\alpha}=2^{50})$.

These ranges are smaller than the original range of the values in CECILIA. As
we also stated in the main paper, $\mathsf{INVSQRT}$ does not provide the full
randomness of MPC. It sacrifices the randomness to some extent for the sake of
the computability of the inverse square root of a secret shared Gram matrix.
$\mathsf{INVSQRT}$ is an attempt towards such a computation via MPC without
sacrificing the randomness of it at all. For the time being, one can choose
the option of utilizing outsource inverse square root of Gram matrices if the
full randomness is required. In that case, the complete randomness of MPC is
valid for every operation performed in the privacy preserving inference on
RKNs.