Privacy-Preserving Logistic Regression Training with A Faster Gradient Variant

## Abstract

Abstract Training logistic regression over encrypted data has been a compelling approach in addressing security concerns for several years. In this paper, we introduce an efficient gradient variant, called q ⁢ u ⁢ a ⁢ d ⁢ r ⁢ a ⁢ t ⁢ i ⁢ c 𝑞 𝑢 𝑎 𝑑 𝑟 𝑎 𝑡 𝑖 𝑐 quadratic italic_q italic_u italic_a italic_d italic_r italic_a italic_t italic_i italic_c g ⁢ r ⁢ a ⁢ d ⁢ i ⁢ e ⁢ n ⁢ t 𝑔 𝑟 𝑎 𝑑 𝑖 𝑒 𝑛 𝑡 gradient italic_g italic_r italic_a italic_d italic_i italic_e italic_n italic_t , for privacy-preserving logistic regression training. We enhance Nesterov’s Accelerated Gradient (NAG), Adaptive Gradient Algorithm (Adagrad) and Adam algorithms by incorporating their quadratic gradients and evaluate these improved algorithms on various datasets. Experimental results demonstrate that the enhanced algorithms achieve significantly improved convergence speed compared to traditional first-order gradient methods. Moreover, we applied the enhanced NAG method to implement homomorphic logistic regression training, achieving comparable results within just 4 iterations.
There is a good chance that the quadratic gradient approach could integrate first-order gradient descent/ascent algorithms with the second-order Newton-Raphson methods, and that it could be applied to a wide range of numerical optimization problems.

## 1 Introduction

Given a person’s healthcare data related to a certain disease, we can train a logistic regression (LR) model capable of telling whether or not this person is likely to develop this disease. However, such personal health information is highly private to individuals. The privacy concern, therefore, becomes a major obstacle for individuals to share their biomedical data.
The most secure solution is to encrypt the data into ciphertexts first by Homomorphic Encryption (HE) and then securely outsource the ciphertexts to the cloud, without allowing the cloud to access the data directly.
iDASH is an annual competition that aims to call for implementing interesting cryptographic schemes in a biological context. Since 2014, iDASH has included the theme of genomics and biomedical privacy. The third track of the 2017 iDASH competition and the second track of the 2018 iDASH competition were both to develop HE-based solutions for building an LR model over encrypted data.

Several studies on logistic regression models are based on homomorphic encryption.
Kim et al. () discussed the problem of performing LR training in an encrypted environment. They employed full-batch gradient descent during the training process and utilized the least-squares method to approximate the sigmoid function.
In the 2017 iDASH competition, Bonte & Vercauteren (), Kim et al. (), Chen et al. (), and Crawford et al. () all addressed the same problem explored by Kim et al. (). In the 2018 iDASH competition, Kim et al. () and Blatt et al. () further developed the problem, focusing on efficient packing and a semi-parallel algorithm. There are other related works (Kim et al., ; Bergamaschi et al., ; Ogilvie et al., ) focusing on various aspects but the papers most relevant to this work are (Bonte & Vercauteren, ) and (Kim et al., ). Bonte & Vercauteren () developed a practical algorithm called the simplified fixed Hessian (SFH) method. Our study extends their work and adopts the ciphertext packing technique proposed by Kim et al. () for efficient homomorphic computation.

Our specific contributions in this paper are as follows:

- 1.
We propose a new gradient variant, quadratic gradient, which can combine the first-order gradient algorithms and the second-order Newton-Raphson method (aka Newton’s method).
- 2.
We develop three enhanced gradient algorithms by equipping the original ones with quadratic gradient. The resulting algorithms all converge and, in most cases, demonstrate strong performance in terms of convergence speed.
- 3.
We implement privacy-preserving logistic regression training using the enhanced NAG method, to our best knowledge, which seems to be a great choice without compromising much on computation and storage.

## 2 Preliminaries

We use the square brackets “$[\ ]$” to denote the index of a vector or matrix element in what follows. For example, for a vector $\bm{v}\in\mathbb{R}^{(n)}$ and a matrix $M\in\mathbb{R}^{m\times n}$, $\bm{v}[i]$ or $\bm{v}_{[i]}$ means the $i$-th element of vector $\bm{v}$ and $M[i][j]$ or $M_{[i][j]}$ the $j$-th element in the $i$-th row of $M$.

### 2.1 Fully Homomorphic Encryption

Fully Homomorphic Encryption (FHE) is a type of cryptographic scheme that can be used to compute an arbitrary number of additions and multiplications directly on the encrypted data. It was not until 2009 that Gentry constructed the first FHE scheme via a bootstrapping operation (Gentry, ). FHE schemes themselves are computationally time-consuming; the choice of dataset encoding matters likewise to the efficiency. In addition to these two limits, how to manage the magnitude of plaintext (Jäschke & Armknecht, ) also contributes to the slowdown. Cheon et al. () proposed a method to construct an HE scheme with a rescaling procedure which could eliminate this technical bottleneck effectively. We adopt their open-source implementation HEAAN while implementing our homomorphic LR algorithms. In addition, it is inevitable to pack a vector of multiple plaintexts into a single ciphertext for yielding a better amortized time of homomorphic computation.
HEAAN supports a parallel technique (aka SIMD) to pack multiple complex numbers in a single polynomial and provides rotation operation on plaintext slots.
The underlying HE scheme in HEAAN is well described in (Kim et al., ; ; Han et al., ).

### 2.2 Database Encoding Method

Kim et al. () proposed an efficient and promising database-encoding method by using SIMD technique, which could make full use of the computation and storage resources.
Suppose a database contains a training dataset $Z$ comprising $n$ samples and $(1+d)$ covariates. The dataset matrix $Z$ is packed into a single ciphertext in a row-by-row manner.

When employing this encoding scheme, we can manipulate the data matrix $Z$ through HE operations on the ciphertext $Enc[Z]$, utilizing only three HE operations - rotation, addition, and multiplication. For instance, if we wish to isolate the first column of $Enc[Z]$ and exclude the other columns, we can create a constant matrix $F$ with ones in the first column and zeros elsewhere. Multiplying $Enc[Z]$ by $F$ will yield the desired ciphertext.

Han et al. () mentioned several basic but important operations used by Kim et al. () in their implementation, such as a procedure named “SumColVec” to compute the summation of the columns of a matrix. With these fundamental operations, more intricate computations, like computing gradients in logistic regression models, become achievable.

### 2.3 Logistic Regression Model

Logistic regression is widely used in binary classification tasks to infer whether a binary-valued variable belongs to a certain class or not. LR can be generalized from linear regression (Murphy, ) by mapping the whole real line $(\bm{\beta}^{\top}\mathbf{x})$ to $(0,1)$ via the sigmoid function $\sigma(z)=1/(1+\exp(-z))$, where the vector $\bm{\beta}\in\mathbb{R}^{(1+d)}$ is the main parameter of LR and the vector $\mathbf{x}=(1,x_{1},\ldots,x_{d})\in\mathbb{R}^{(1+d)}$ the input covariate. Thus logistic regression can be formulated with the class label $y\in\{\pm 1\}$ as follows:

$$ $\displaystyle\Pr(y=+1|\mathbf{x},\bm{\beta})$ $\displaystyle=\sigma(\bm{\beta}^{\top}\mathbf{x})$ $\displaystyle=\frac{1}{1+e^{-\bm{\beta}^{\top}\mathbf{x}}},$ $\displaystyle\Pr(y=-1|\mathbf{x},\bm{\beta})$ $\displaystyle=1-\sigma(\bm{\beta}^{\top}\mathbf{x})$ $\displaystyle=\frac{1}{1+e^{+\bm{\beta}^{\top}\mathbf{x}}}.$ $$

LR sets a threshold (usually $0.5$) and compares its output with it to decide the resulting class label.

The logistic regression problem can be transformed into an optimization problem that seeks a parameter $\bm{\beta}$ to maximize $L(\bm{\beta})=\prod_{i=1}^{n}\Pr(y_{i}|\mathbf{x}_{i},\bm{\beta})$ or its log-likelihood function $l(\bm{\beta})$ for convenience in the calculation:

$$ $\displaystyle l(\bm{\beta})=\ln L(\bm{\beta})=-\sum_{i=1}^{n}\ln(1+e^{-y_{i} \bm{\beta}^{\top}\mathbf{x}_{i}}),$ $$

where $n$ is the number of examples in the training dataset. LR does not have a closed form of maximizing $l(\bm{\beta})$ and two main methods are adopted to estimate the parameters of an LR model: (a) gradient descent method via the gradient; and (b) Newton’s method by the Hessian matrix. The gradient and Hessian of the log-likelihood function $l(\bm{\beta})$ are given by, respectively:

$$ $\displaystyle\nabla_{\bm{\beta}}l(\bm{\beta})$ $\displaystyle=\sum_{i}(1-\sigma(y_{i}\bm{\beta}^{\top}\mathbf{x}_{i}))y_{i} \mathbf{x}_{i},$ $\displaystyle\nabla_{\bm{\beta}}^{2}l(\bm{\beta})$ $\displaystyle=\sum_{i}(y_{i}\mathbf{x}_{i})(\sigma(y_{i}\bm{\beta}^{\top} \mathbf{x}_{i})-1)\sigma(y_{i}\bm{\beta}^{\top}\mathbf{x}_{i})(y_{i}\mathbf{x} _{i})$ $\displaystyle=X^{\top}SX,$ $$

where $S$ is a diagonal matrix with entries $S_{ii}=(\sigma(y_{i}\bm{\beta}^{\top}\mathbf{x}_{i})-1)\sigma(y_{i}\bm{\beta}^
{\top}\mathbf{x}_{i})$ and $X$ the dataset.

The log-likelihood function $l(\bm{\beta})$ of LR has at most a unique global maximum (Allison, ), where its gradient is zero. Newton’s method is a second-order technique to numerically find the roots of a real-valued differentiable function, and thus can be used to solve the $\bm{\beta}$ in $\nabla_{\bm{\beta}}l(\bm{\beta})=0$ for LR.

## 3 Technical Details

It is quite time-consuming to compute the Hessian matrix and its inverse in Newton’s method for each iteration. One way to limit this downside is to replace the varying Hessian with a fixed matrix $\bar{H}$. This novel technique is called the fixed Hessian Newton’s method. Böhning & Lindsay () have shown that the convergence of Newton’s method is guaranteed as long as $\bar{H}\leq\nabla_{\bm{\beta}}^{2}l(\bm{\beta})$, where $\bar{H}$ is a symmetric negative-definite matrix independent of $\bm{\beta}$ and “$\leq$” denotes the Loewner ordering in the sense that the difference $\nabla_{\bm{\beta}}^{2}l(\bm{\beta})-\bar{H}$ is non-negative definite. With such a fixed Hessian matrix $\bar{H}$, the iteration for Newton’s method can be simplified to:

$$ $\bm{\beta}_{t+1}=\bm{\beta}_{t}-\bar{H}^{-1}\nabla_{\bm{\beta}}l(\bm{\beta}).$ $$

Böhning and Lindsay also suggest the fixed matrix $\bar{H}=-\frac{1}{4}X^{\top}X$ is a good lower bound for the Hessian of the log-likelihood function $l(\bm{\beta})$ in LR.

### 3.1 Simplified Fixed Hessian

Bonte & Vercauteren () simplify this lower bound $\bar{H}$ further due to the need for inverting the fixed Hessian in the encrypted domain. They replace the matrix $\bar{H}$ with a diagonal matrix $B$ whose diagonal elements are simply the sums of each row in $\bar{H}$. They also suggest a specific order of calculation to optimize the computation of $B$ more efficiently.
Their new approximation $B$ of the fixed Hessian is:

$$ $B=\left[\begin{array}[]{cccc}\sum_{i=0}^{d}\bar{h}_{0i}&0&\ldots&0\\ 0&\sum_{i=0}^{d}\bar{h}_{1i}&\ldots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\ldots&\sum_{i=0}^{d}\bar{h}_{di}\\ \end{array}\right],$ $$

where $\bar{h}_{ki}$ is the element of $\bar{H}$.
This diagonal matrix $B$ is in a very simple form and can be obtained from $\bar{H}$ without much difficulty. The inverse of $B$ can even be approximated in the encrypted form by computing the inverse of each diagonal element of $B$ using an iterative Newton’s method with an appropriate initial value.
Their simplified fixed Hessian method can be formulated as follows:

$$ $\displaystyle\bm{\beta}_{t+1}$ $\displaystyle=\bm{\beta}_{t}-B^{-1}\cdot\nabla_{\bm{\beta}}l(\bm{\beta}),$ $\displaystyle=\bm{\beta}_{t}-\left[\begin{array}[]{cccc}b_{00}&0&\ldots&0\\ 0&b_{11}&\ldots&0\\ \vdots&\vdots&\ddots&\vdots\\ 0&0&\ldots&b_{dd}\\ \end{array}\right]\cdot\left[\begin{array}[]{c}\nabla_{0}\\ \nabla_{1}\\ \vdots\\ \nabla_{d}\\ \end{array}\right]=\bm{\beta}_{t}-\left[\begin{array}[]{c}b_{00}\cdot\nabla_{0 }\\ b_{11}\cdot\nabla_{1}\\ \vdots\\ b_{dd}\cdot\nabla_{d}\\ \end{array}\right],$ $$

where $b_{ii}$ is the reciprocal of $\sum_{i=0}^{d}\bar{h}_{0i}$ and $\nabla_{i}$ is the element of $\nabla_{\bm{\beta}}l(\bm{\beta})$.

Consider a special situation: if all the elements $b_{00},\ldots,b_{dd}$ had the same value $-\eta$ with $\eta>0$, the iterative formula of the SFH method could be given as:

$$ $\displaystyle\bm{\beta}_{t+1}$ $\displaystyle=\bm{\beta}_{t}-(-\eta)\cdot\left[\begin{array}[]{c}\nabla_{0}\\ \nabla_{1}\\ \vdots\\ \nabla_{d}\\ \end{array}\right]=\bm{\beta}_{t}+\eta\cdot\nabla_{\bm{\beta}}l(\bm{\beta}),$ $$

which is the same as the formula of the naive gradient $ascent$ method. Such a coincidence not only helps generate the idea behind this work but also leads us to believe that there is a connection between the Hessian matrix and the learning rate of the gradient (descent) method.

We regard $B^{-1}\cdot\nabla_{i}$ as a novel enhanced gradient variant and allocate a distinct learning rate to it. As long as we ensure that this new learning rate decreases from a positive floating-point number greater than 1 (such as 2) to 1 in a bounded number of iteration steps, the fixed Hessian Newton’s method guarantees the algorithm will converge eventually.

The SFH method proposed by Bonte & Vercauteren () has two limitations: (a) in the construction of the simplified fixed Hessian, all entries in the symmetric matrix $\bar{H}$ need to be non-positive. For machine learning applications, datasets are typically normalized in advance to the range [0,1], satisfying the convergence condition of the SFH method. However, in other cases, such as numerical optimization, this condition may not always hold; and (b) the simplified fixed Hessian matrix $B$, as well as the fixed Hessian $\bar{H}=-\frac{1}{4}X^{\top}X$ , can still be singular, especially when the dataset is a high-dimensional sparse matrix, such as the MNIST datasets.
We extend SFH by removing these limitations to generalize this simplified fixed Hessian to be invertible in any case and propose a faster gradient variant, which we term $quadratic$ $gradient$.

### 3.2 Quadratic Gradient Definition

Suppose that a differentiable scalar-valued function $F(\mathbf{x})$ has its gradient $\bm{g}$ and Hessian matrix $H$, with any matrix $\bar{H}\leq H$ in the Loewner ordering for a $maximization$ problem as follows:

$$ $\displaystyle\bm{g}=\left[\begin{array}[]{c}g_{0}\\ g_{1}\\ \vdots\\ g_{d}\\ \end{array}\right],$ $\displaystyle H=\left[\begin{array}[]{cccc}\nabla_{00}^{2}&\nabla_{01}^{2}& \ldots&\nabla_{0d}^{2}\\ \nabla_{10}^{2}&\nabla_{11}^{2}&\ldots&\nabla_{1d}^{2}\\ \vdots&\vdots&\ddots&\vdots\\ \nabla_{d0}^{2}&\nabla_{d1}^{2}&\ldots&\nabla_{dd}^{2}\\ \end{array}\right],\quad\bar{H}$ $\displaystyle=\left[\begin{array}[]{cccc}\bar{h}_{00}&\bar{h}_{01}&\ldots&\bar {h}_{0d}\\ \bar{h}_{10}&\bar{h}_{11}&\ldots&\bar{h}_{1d}\\ \vdots&\vdots&\ddots&\vdots\\ \bar{h}_{d0}&\bar{h}_{d1}&\ldots&\bar{h}_{dd}\\ \end{array}\right],$ $$

where $\nabla_{ij}^{2}=\nabla_{ji}^{2}=\frac{\partial^{2}F}{\partial x_{i}\partial x_
{j}}$. We construct a new diagnoal Hessian matrix $\tilde{B}$ with each diagnoal element $\tilde{B}_{kk}$ being $-\epsilon-\sum_{i=0}^{d}|\bar{h}_{ki}|$,

$$ $\displaystyle\tilde{B}=\left[\begin{array}[]{ccc}-\epsilon-\sum_{i=0}^{d}|\bar {h}_{0i}|&0&\ldots\\ 0&-\epsilon-\sum_{i=0}^{d}|\bar{h}_{1i}|&\ldots\\ \vdots&\vdots&\ddots\\ 0&0&\ldots\\ \end{array}\right],$ $$

where $\epsilon$ is a small positive constant to avoid division by zero (usually set to $1e-8$).

As long as $\tilde{B}$ satisfies the convergence condition of the aforementioned fixed Hessian method, namely $\tilde{B}\leq H$, we can use this approximation $\tilde{B}$ of the Hessian matrix as a lower bound. Given that we have already assumed $\bar{H}\leq H$, it suffices to demonstrate that $\tilde{B}\leq\bar{H}$. We establish $\tilde{B}\leq\bar{H}$ in a manner similar to that employed by Bonte & Vercauteren ().

###### Lemma 3.1 .

Let $A\in\mathbb{R}^{n\times n}$ be a symmetric matrix, and let $B$ be the diagonal matrix whose diagonal entries $B_{kk}=-\epsilon-\sum_{i}|A_{ki}|$ for $k=1,\ldots,n$, then $B\leq A$.

###### Proof.

By definition of the Loewner ordering, we have to prove the difference matrix $C=A-B$ is non-negative definite, which means that all the eigenvalues of $C$ need to be non-negative. By construction of $C$ we have that $C_{ij}=A_{ij}+\epsilon+\sum_{k=1}^{n}|A_{ik}|$ for $i=j$ and $C_{ij}=A_{ij}$ for $i\neq j$. By means of Gerschgorin’s circle theorem, we can bound every eigenvalue $\lambda$ of $C$ in the sense that $|\lambda-C_{ii}|\leq\sum_{i\neq j}|C_{ij}|$ for some index $i\in\{1,2,\ldots,n\}$. We conclude that $\lambda\geq A_{ii}+\epsilon+|A_{ii}|\geq\epsilon>0$ for all eigenvalues $\lambda$ and thus that $B\leq A$.
∎

###### Definition 3.2 ( Quadratic Gradient ) .

Given the previously defined $\tilde{B}$, we define the quadratic gradient as $G=\bar{B}\cdot\bm{g}$ with a new learning rate $N_{t}$, where $\bar{B}$ is a diagonal matrix with diagonal entries $\bar{B}_{kk}=1/|\tilde{B}_{kk}|$. The learning rate $N_{t}$ should always be no less than 1 and is designed to decrease to 1 over a limited number of iterations. It is important to note that $G$ remains a column vector of the same dimension as the gradient $\bm{g}$. To maximize the function $F(\mathbf{x})$, we can employ the iterative formula: $\mathbf{x}_{t+1}=\mathbf{x}_{t}+N_{t}\cdot G$, similar to the naive gradient approach.

Minimizing the function $F(x)$ is equivalent to maximizing the function $-F(x)$. In this context, we need to construct $\tilde{B}$ using either a good $lower$ bound $\bar{H}$ for the Hessian $-H$ of $-F(x)$ or a good $upper$ bound $\bar{H}$ for the Hessian $H$ of $F(x)$. It is worth noting that $\bar{H}$ can be taken as the Hessian $H$ itself, thereby eliminating the need to find either lower or upper bounds.

Böhning & Lindsay () did not propose a systematic method for determining or constructing a constant Hessian approximation for their fixed Hessian method. This may be due to the absence of such matrices for most objective functions. We note that our quadratic gradient approach does not require a constant (fixed) Hessian replacement and instead allows for the use of the original Hessian itself to construct our varying diagonal matrix $\tilde{B}$ for the quadratic gradient. Under this less restrictive condition, we present a methodical approach to evaluate whether a fixed Hessian matrix can be employed for efficient computation. Specifically, we first construct $\tilde{B}$ from the Hessian of the given objective function and then check if, for each diagonal element $\tilde{B}_{kk}$, there is a constant maximum for $\bar{B}_{kk}=1/|\tilde{B}_{kk}|$ (or a constant minimum $|\tilde{B}_{kk}|$). If each diagonal element $|\tilde{B}_{kk}|$ achieves its minimum $|\tilde{M}_{kk}|$, then a constant diagonal matrix with the corresponding diagonal elements $|\tilde{M}_{kk}|$ would satisfy the convergence condition for the SFH method. Finally, we should evaluate the performance of such a fixed Hessian replacement, as it might not always provide a good bound. This search method can be applied to any optimization function to identify a fixed Hessian. However, the absence of a constant Hessian replacement via this method does not necessarily imply that such an approximation does not exist.

Figure: Algorithm 1 The Enhanced Nesterov’s Accelerated Gradient Algorithm

###### Lemma 3.1 .

###### Proof.

###### Definition 3.2 ( Quadratic Gradient ) .

### 3.3 Quadratic Gradient Algorithms

Quadratic gradient can be used to enhance various first-order gradient algorithms:

- 1.
NAG is a different variant of the momentum method to give the momentum term much more prescience. The iterative formulas of the gradient $ascent$ method for NAG are as follows:
$\displaystyle V_{t+1}$
$\displaystyle=\bm{\beta}_{t}+\eta_{t}\cdot\nabla J(\bm{\beta}_{t}),$
(3)
$\displaystyle\bm{\beta}_{t+1}$
$\displaystyle=(1-\gamma_{t})\cdot V_{t+1}+\gamma_{t}\cdot V_{t},$
(4)
where $V_{t+1}$ is the intermediate variable used for updating the final weight $\bm{\beta}_{t+1}$ and $\gamma_{t}\in(0,1)$ is a smoothing parameter of moving average to evaluate the gradient at an approximate future position (Kim et al., ). The enhanced NAG is to replace equation 3 with $V_{t+1}=\bm{\beta}_{t}+N_{t}\cdot G,$ where $N_{t}$ is the new learning rate for the enhanced algorithm usually setted to $1+\eta_{t}$. Our enhanced NAG is described in Algorithm 1 .
- 2.
Adagrad is a gradient-based algorithm suitable for dealing with sparse data. The updated operations of Adagrad and its quadratic-gradient version, for every parameter $\bm{\beta}_{[i]}$ at each iteration step $t$, are as follows, respectively:
$\displaystyle\bm{\beta}_{[i]}^{(t+1)}$
$\displaystyle=\bm{\beta}_{[i]}^{(t)}-\frac{\eta_{t}}{\epsilon+\sqrt{\sum_{k=1}
^{t}\bm{g}_{[i]}^{(t)}\cdot\bm{g}_{[i]}^{(t)}}}\cdot\bm{g}_{[i]}^{(t)},$
$\displaystyle\bm{\beta}_{[i]}^{(t+1)}$
$\displaystyle=\bm{\beta}_{[i]}^{(t)}-\frac{N_{t}}{\epsilon+\sqrt{\sum_{k=1}^{t
}G_{[i]}^{(t)}\cdot G_{[i]}^{(t)}}}\cdot G_{[i]}^{(t)}.$
- 3.
Adam is an optimization algorithm that combines the benefits of Adagrad, which adapts learning rates for each parameter, with momentum, allowing it to converge faster and handle sparse gradients more effectively. Similar to how Adagrad transitions to its enhanced version, the update formula for enhanced Adam, also employs a new learning rate $N_{t}$ and replaces the gradient with its quadratic form.

Figure: (a) The iDASH dataset
Refer to caption: x1.png

Figure: (a) The iDASH dataset
Refer to caption: x9.png

Figure: (a) The iDASH dataset
Refer to caption: x17.png

Performance Evaluation
We evaluate the performance of various algorithms in a non-encrypted context using Python on a desktop computer with an Intel Core G640 CPU at 1.60 GHz and 7.3 GB of RAM. Given that our focus is on the convergence speed of the algorithms during training, we use the loss function $l(\bm{\beta})$, specifically maximum likelihood estimation (MLE), as the only indicator. We evaluate six algorithms—NAG, Adagrad, Adam, and their quadratic-gradient variants (denoted as Enhanced NAG, Enhanced Adagrad, and Enhanced Adam, respectively)—using the datasets adopted by Kim et al. (): the iDASH genomic dataset (iDASH), the Myocardial Infarction dataset from Edinburgh (edin), the Low Birth Weight Study (lbw), Nhanes III (nhanes3), the Prostate Cancer Study (pcs), and the Umaru Impact Study datasets (uis). The genomic dataset from the third task of the 2017 iDASH competition comprises 1,579 records, each featuring 103 binary genotypes and a binary phenotype indicating whether the patient has cancer. The remaining five datasets each contain a single binary dependent variable. We also evaluate these algorithms on two large datasets from Han et al. (): a real financial dataset with 422,108 samples and 200 features, and the restructured public MNIST dataset, which includes 11,982 samples from the training set with 196 features. For a fair comparison with the baseline work by Kim et al. (), the enhanced NAG algorithm employs the learning rate $1+\frac{10}{1+t}$, consistent with Kim et al.’s choice of $\frac{10}{1+t}$ for their learning rate. The enhanced Adagrad algorithm utilizes $N_{t}=1+0.01$ as the learning rate, while the enhanced Adam algorithm adopts the following learning rate settings: $\alpha=1+0.001$, $\beta_{1}=0.9$, and $\beta_{2}=0.999$. Refer to (Kingma & Ba, ) for further details regarding these parameters.
In our experiments, we consistently use $\bar{H}=-\frac{1}{4}X^{\top}X$ to construct $\tilde{B}$ for the function $l(\bm{\beta})$.

From the empirical results presented in Figures 1, 2, and 3, we conclude that all the enhanced algorithms converge, and that the enhanced NAG algorithm consistently outperforms the original NAG method in terms of convergence rate.
Although the enhanced Adagrad and Adam methods do not achieve faster convergence than their original counterparts on the iDASH genomic dataset, they demonstrate clear advantages in all other cases. The reason for this discrepancy is beyond the scope of this paper and will be addressed as part of future work. One possible explanation for the superior performance of the enhanced algorithms is that the quadratic gradient integrates curvature information into first-order gradient methods.

An important observation is that the enhanced algorithms demonstrate better performance with learning rates between $1$ and $2$ compared to other values. When the learning rate exceeds $3$, the algorithm is nearly guaranteed not to converge. For quadratic gradient algorithms, we recommend employing an exponentially decaying learning rate, such as $1+A\cdot\gamma^{t}$, where $t$ denotes the iteration number, $A$ is a positive constant typically set to no less than $1$, and $\gamma$ is a positive number less than $1$ that controls the rate of decay.

Gradient And Quadratic Gradient
We executed the raw gradient ascent algorithm and the raw quadratic gradient algorithm on the lbw dataset using various learning rates. Figure 4 presents the detailed results of this experiment. It was precisely these results that directly inspired the authors to develop the concept of the quadratic gradient. Similar to the gradient, the quadratic gradient, which incorporates curvature, also demonstrates smooth progression. In particular, as the learning rate is gradually increased, both the gradient and the quadratic gradient reveal corresponding gradual changes in performance, rather than sudden jumps.

Figure: (a) The raw first-order gradient ascent algorithm with iterative formulas: $\bm{\beta}_{t+1}=\bm{\beta}_{t}+\mathrm{lr}\cdot\bm{g}$
Refer to caption: x25.png

## 4 Secure Training

We adopt the enhanced NAG method to implement secure logistic regression training based on HE.
The difficulty in applying the quadratic gradient is to invert the diagonal matrix $\tilde{B}$ in order to obtain $\bar{B}$. We leave the computation of matrix $\bar{B}$ to data owner and let the data owner upload the ciphertext encrypting the $\bar{B}$ to the cloud. Since data owner has to prepare the dataset and normalize it, it would also be practicable for the data owner to calculate the $\bar{B}$ owing to no leaking of sensitive data information.

Privacy-preserving logistic regression training based on HE techniques faces a difficult dilemma that no homomorphic schemes are capable of directly calculating the sigmoid function in the LR model. A common solution is to replace the sigmoid function with a polynomial approximation by using the widely adopted least-squares method. We can call a function named “ polyfit($\cdot$) ” in the Python package Numpy to fit the polynomial in a least-square sense. We adopt the degree 5 polynomial approximation $g(x)$ developed by Kim et al. (), utilizing the least squares approach to approximate the sigmoid function over the interval $[-8,8]$:
$g(x)=0.5+0.19131\cdot x-0.0045963\cdot x^{3}+0.0000412332\cdot x^{5}$ .

Given the training dataset $\text{X}\in\mathbb{R}^{n\times(1+d)}$ and training label $\text{Y}\in\mathbb{R}^{n\times 1}$, we adopt the same method that Kim et al. () used to encrypt the data matrix consisting of the training data combined with training-label information into a single ciphertext $\text{ct}_{Z}$. The weight vector $\beta^{(0)}$ consisting of zeros and the diagnoal elements of $\bar{B}$ are copied $n$ times to form two matrices. The data owner then encrypt the two matrices into two ciphertexts $\text{ct}_{\beta}^{(0)}$ and $\text{ct}_{\bar{B}}$, respectively. The ciphertexts $\text{ct}_{Z}$, $\text{ct}_{\beta}^{(0)}$ and $\text{ct}_{\bar{B}}$ are as follows:

$$ $\displaystyle\text{X}=\left[\begin{array}[]{cccc}1&x_{11}&\ldots&x_{1d}\\ 1&x_{21}&\ldots&x_{2d}\\ \vdots&\vdots&\ddots&\vdots\\ 1&x_{n1}&\ldots&x_{nd}\\ \end{array}\right],\text{Y}=\left[\begin{array}[]{c}y_{1}\\ y_{2}\\ \vdots\\ y_{n}\\ \end{array}\right],$ $\displaystyle\text{ct}_{Z}=Enc\left[\begin{array}[]{cccc}y_{1}&y_{1}x_{11}& \ldots&y_{1}x_{1d}\\ y_{2}&y_{2}x_{21}&\ldots&y_{2}x_{2d}\\ \vdots&\vdots&\ddots&\vdots\\ y_{n}&y_{n}x_{n1}&\ldots&y_{n}x_{nd}\\ \end{array}\right],$ $\displaystyle\text{ct}_{\beta}^{(0)}=Enc\left[\begin{array}[]{cccc}\beta_{0}^{ (0)}&\beta_{1}^{(0)}&\ldots&\beta_{d}^{(0)}\\ \beta_{0}^{(0)}&\beta_{1}^{(0)}&\ldots&\beta_{d}^{(0)}\\ \vdots&\vdots&\ddots&\vdots\\ \beta_{0}^{(0)}&\beta_{1}^{(0)}&\ldots&\beta_{d}^{(0)}\\ \end{array}\right],$ $\displaystyle\text{ct}_{\bar{B}}=Enc\left[\begin{array}[]{cccc}\bar{B}_{[0][0] }&\bar{B}_{[1][1]}&\ldots&\bar{B}_{[d][d]}\\ \bar{B}_{[0][0]}&\bar{B}_{[1][1]}&\ldots&\bar{B}_{[d][d]}\\ \vdots&\vdots&\ddots&\vdots\\ \bar{B}_{[0][0]}&\bar{B}_{[1][1]}&\ldots&\bar{B}_{[d][d]}\\ \end{array}\right],$ $$

where $\bar{B}_{[i][i]}$ is the diagonal element of $\bar{B}$ that is built from $-\frac{1}{4}X^{\top}X$.

The pulbic cloud takes the three ciphertexts $\text{ct}_{Z}$, $\text{ct}_{\beta}^{(0)}$ and $\text{ct}_{\bar{B}}$ and evaluates the enhanced NAG algorithm to find a decent weight vector by updating the vector $\text{ct}_{\beta}^{(0)}$. Refer to (Kim et al., ) for a detailed description about how to calculate the gradient by HE programming.

#### Limitations

In a privacy-preserving setting, when compared to the NAG method, the primary limitation of the Enhanced NAG method is that it requires one additional ciphertext multiplication to construct the quadratic gradient. In addition, the data owner needs to upload one more ciphertext encrypting the matrix $\bar{B}$. However, the enhanced algorithm converges faster, and we believe it can compensate for the mentioned shortcomings.

**Table 1: Implementation Results for iDASH datasets with 10-fold CV**
| Dataset | Sample<br>Num | Feature<br>Num | Method | Iter<br>Num | Storage<br>(GB) | Learn<br>Time<br>(min) | Accuracy<br>(%) | AUC |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| iDASH | 1579 | 18 | Ours | 4 | 0.08 | 4.43 | 61.46 | 0.696 |
| Baseline | 7 | 0.04 | 6.07 | 62.87 | 0.689 |  |  |  |

**Table 2: Implementation Results for other datasets with 5-fold CV**
| Dataset | Sample<br>Num | Feature<br>Num | Method | Iter<br>Num | Storage<br>(GB) | Learn<br>Time<br>(min) | Accuracy<br>(%) | AUC |
| --- | --- | --- | --- | --- | --- | --- | --- | --- |
| Edin | 1253 | 8 | Ours | 4 | 0.04 | 0.6 | 89.52 | 0.943 |
| Baseline | 7 | 0.02 | 3.6 | 91.04 | 0.958 |  |  |  |
| lbw | 189 | 8 | Ours | 4 | 0.04 | 0.6 | 71.35 | 0.667 |
| Baseline | 7 | 0.02 | 3.3 | 69.19 | 0.689 |  |  |  |
| nhanes3 | 15649 | 15 | Ours | 4 | 0.31 | 4.5 | 79.23 | 0.637 |
| Baseline | 7 | 0.16 | 7.3 | 79.22 | 0.717 |  |  |  |
| pcs | 379 | 9 | Ours | 4 | 0.04 | 0.6 | 63.20 | 0.733 |
| Baseline | 7 | 0.02 | 3.5 | 68.27 | 0.740 |  |  |  |
| uis | 575 | 8 | Ours | 4 | 0.04 | 0.6 | 74.43 | 0.597 |
| Baseline | 7 | 0.02 | 3.5 | 74.44 | 0.603 |  |  |  |

## 5 Experiments

#### Implementation

We implement the enhanced NAG based on HE with the library HEAAN. The C++ source code is publicly available at [https://anonymous.4open.science/r/IDASH2017-245B](https://anonymous.4open.science/r/IDASH2017-245B) . All the experiments on the ciphertexts were conducted on a public cloud with 32 vCPUs and 64 GB RAM.

For a fair comparison with the baseline (Kim et al., ), we utilized the same 10-fold cross-validation (CV) technique on the same iDASH dataset consisting of 1579 samples with 18 features and the same 5-fold CV technique on the other five datasets. Like (Kim et al., ), We consider the average accuracy and the Area Under the Curve (AUC) as the main indicators. Tables 1 and 2 display the results of the two experiments, respectively. The two tables also provide the average evaluation running time for each iteration and the storage (encrypted dataset for the baseline work and encrypted dataset and $\bar{B}$ for our method). We adopt the same packing method that Kim et al. () proposed and hence our solution has similar storage of ciphertexts to (Kim et al., ) with some extra ciphertexts to encrypt the $\bar{B}$. We chose $1+0.9^{t}$ as our learning rate configuration.

The parameters of HEAAN we set are same to (Kim et al., ): $logN=16$, $logQ=1200$, $logp=30$, $slots=32768$, which ensure the security level $\lambda=80$. We use a larger $logp=40$ to encrypt the matrix ${\bar{B}}$ for preserving the precision of ${\bar{B}}$. Refer to (Kim et al., ) for details on these parameters. Since our enhanced NAG method need one more ciphertext multiplication than the baseline work, consuming more modulus, our solution thus can only perform $4$ iterations of the enhanced NAG method. Yet despite only $4$ iterations, our enhanced NAG method still produces a comparable result.

## 6 Conclusion

In this paper, we proposed a faster gradient variant, termed quadratic gradient, and implemented the quadratic-gradient version of NAG in the encrypted domain to train a logistic regression model. The quadratic gradient introduced in this work can be directly constructed from the Hessian, effectively merging first-order gradient (descent) methods with the second-order Newton’s method. There is a promising potential for the quadratic gradient to accelerate other gradient methods, and it may serve as an alternative or even a replacement for the traditional line-search method.

## Peer Reviews

confidence: 3
contribution: 2
initial_rating: 3
presentation: 3
questions: 1. In Tables 1 and 2, the proposed enhanced NAG method generally performs worse than the baseline method for accuracy and AUC. It is also mentioned then when compared to the Kim et al. baseline, the enhanced NAG requires one additional ciphertext multiplication and so you run it for less iterations.  Is the only advantage of the enhanced NAG computation time?
2. In addition to the point above, there could be a more thorough set of experiments for the privacy-preserving setting. Restricting the comparison to 4 iterations vs. 7 does not do enough to highlight when the enhanced NAG method should be used. There is a likely trade-off in terms of convergence and compute that is not being fully explored in the experiments.
3. As far as I can tell, the proposed quadratic gradient method is a straightforward extension to the work of Bonte and Vercauteren (2018). Why do you not compare to their method in Section 5? I feel the privacy-preserving experiments are lacking further baselines.
4. In Line 185, is the gradient meant to have a subscript i here? Is is meant to be referring to $\nabla_\beta$?
5. Minor point: There are a few typos that could be corrected — L477: pulbic → public, L248: practicable -> practical, L269: setted → set
soundness: 2
strengths: - The paper tackles an important and useful problem, focusing on improving the convergence speed of privacy-preserving logistic regression.
- The paper is clearly written and well-structured.
summary: The authors propose a privacy-preserving approach to training logistic regression based on homomorphic encryption. Their approach utilises second-order information through quadratic gradients to implement various optimisers such as NAG/Adagrad and claim this helps achieve improved convergence speeds compared to first-order methods.
weaknesses: - The evaluation in the privacy-preserving setting is limited. Most of the experiments focus on the non-private setting and the private experiments compare only to a single baseline method.
- The empirical results are often unconvincing. It is not clear that the proposed quadratic gradient approach gives any clear benefit on some of the benchmark datasets in terms of convergence speed or accuracy/AUC.
- There is not enough of a comparison with prior privacy-preserving logistic regression methods to enable a clear understanding of why this method should be used compared to prior work.

---

confidence: 3
contribution: 2
initial_rating: 3
presentation: 3
questions: - Further insights (through experiments and/or theory) should be provided on the performance of the quadratic gradient approach (see weaknesses).
soundness: 2
strengths: The quadratic gradient approach provides a new combination of first-order gradient and second-order Newton-Raphson methods. Multiple experiments were considered to demonstrate the improvements in convergence speed achieved by the enhanced NAG, Adagrad, and Adam algorithms. The paper aims to address a practical challenge in encrypted logistic regression training, where efficient computation is crucial.
The quadratic gradient concept builds on several ideas from the literature, which are well reported in the paper.
summary: The paper proposes a new optimization technique for privacy-preserving logistic regression training called quadratic gradient. This method combines elements from first order gradient algorithms and the second-order Newton-Raphson method. The variant constructs a diagonal Hessian approximation that allows for the computation of an optimized gradient update without recalculating the computationally expensive Hessian matrix. As a result, the quadratic gradient facilitates faster convergence by leveraging curvature information while retaining the simplicity of first-order methods.
The authors show that quadratic-gradient versions of Nesterov Accelerated Gradient (NAG), Adagrad, and Adam optimizers converge, and, in most cases, demonstrate strong performance in terms of convergence speed, by empirical evaluation on multiple datasets. The authors apply the quadratic gradient-enhanced NAG to homomorphic encryption-based logistic regression training, showing comparable model accuracy within only a few iterations.
weaknesses: - The claim that the enhanced algorithms achieve significantly improved convergence speed compared to traditional first-order gradient
methods is not very well supported, as the performance of the algorithms varies across different datasets. The enhanced Adagrad and Adam methods do not achieve faster convergence on the iDASH, and while they do perform better on other datasets, it is hard to assess the overall extent and consistency of the improvement. Namely, the authors mention: _"Although the enhanced Adagrad and Adam methods do not achieve faster convergence than their original counterparts on the iDASH genomic dataset, they demonstrate clear advantages in all other cases. The reason for this discrepancy is beyond the scope of this paper and will be addressed as part of future work."_. It is worth providing additional analysis in this work on why certain datasets showed less improvement, as it would help identify conditions under which the method is most effective (i.e., how do the methods behave under different number of dimensions, dataset size, etc.). Providing bounds on convergence speed or computational efficiency would also strengthen the contribution.
- The paper's writing makes its focus somewhat unclear. The quadratic gradient method is proposed specifically as a technique for privacy-preserving logistic regression, but a significant portion of the paper is dedicated to showing that quadratic gradient can perform better than first-order gradient methods in the general (non-private) setting. The private setting is proposed as an application.

---

confidence: 4
contribution: 2
initial_rating: 3
presentation: 2
questions: - The empirical results presented in Section 5 require further explanation. The proposed method uses twice as much memory as the baseline, which is understandable given that it needs to store diagonal entries of the preconditioning matrix. However, it seems to be faster than the baseline in terms of the learn time. What exactly does the learn time measure and why is it the case? Given that the proposed algorithm includes an addition step compared to the baseline, it is not clear why the proposed algorithm is faster.
- What is the motivation for approximating the diagonal entries of Hessian with the row-wise sum of $\bar{H}$?
- I am not sure how the formulation presented in line 176 - 181 help better explain the idea of proposed work.
soundness: 3
strengths: - The proposed algorithm employs the fixed diagonal approximation of Hessian which reduces computational cost for calculating the Hessian during the training and can be approximated in the encrypted form.
- The proposed quadratic gradient is applied to three optimization algorithms: Nesterov’s Accelerated Gradient (NAG), AdaGrad, and Adam, and their performance is compared across various datasets.
summary: This paper presents a secure training algorithm for logistic regression models. Specifically, the proposed training algorithm uses a gradient estimator called quadratic gradient which preconditions the gradient using the estimated second-order information. To reduce the computational burden, the second-order information (i.e., Hessian) is approximated as a diagonal matrix. The training algorithm is implemented with a homomorphic encryption framework.
weaknesses: - The proposed quadratic gradient seems to be similar to the preconditioned gradient, which has long been utilized in the optimization community, and there exists a substantial body of literature on this topic. The proposed work seems to make an incremental contribution, as the quadratic gradient appears to be a rebranding of the existing preconditioned gradient. In addition, the idea of approximating the parameter-dependent Hessian matrix with a parameter-independent input covariance matrix was introduced in the literature. For example,
  1. Ida, Yasutoshi, Yasuhiro Fujiwara, and Sotetsu Iwamura. “Adaptive Learning Rate via Covariance Matrix Based Preconditioning for Deep Neural Networks.” In Proceedings of the Twenty-Sixth International Joint Conference on Artificial Intelligence, IJCAI-17
  2. Mehta, Harsh, Walid Krichene, Abhradeep Guha Thakurta, Alexey Kurakin, and Ashok Cutkosky. “Differentially Private Image Classificafrom Features.” Transactions on Machine Learning Research, November 29, 2022.

- The proposed algorithm is implemented using the existing homomorphic encryption framework. It is unclear what the significance of the use of quadratic gradient in the context of homomorphic encryption.
- The empirical results are not well explained. The performance comparison with the existing algorithm by Kim et al. shows that the proposed algorithm achieves lower accuracy than the baseline method on some datasets (for example, Edin and pcs datasets). 
- The presentation of the paper can be improved.

---

confidence: 4
contribution: 1
initial_rating: 3
presentation: 2
questions: I do not have any questions.
soundness: 2
strengths: This paper proposes a preconditioner that improves the training loss or training convergence of 3 different styles of gradient descent on 8 datasets.
summary: This paper proposes a preconditioning scheme for accelerating gradient-descent style algorithms. The paper explains how to incorporate this scheme into the logistic regression fitting with homomorphic encryption framework introduced in prior work (Kim et al.  2018a).

Experiments on a variety of datasets in very tiny plots show that the proposed preconditioner improves the training convergence and/or training loss for nesterov accelerated gradient, adagrad, and adam. However, comparisons in the homomorphic encryption framework to Kim et al. are inconclusive as often the Kim et al. version can improve accuracy at a cost of a few more minutes of training.
weaknesses: Before discussing  weaknesses, it is important to add some context to the paper. Multiplying the gradient by a positive (or negative) definite matrix during gradient descent/ascent is known as preconditioning.  This is perhaps the connection the paper was looking for on Line 178. Thus this paper is mostly about finding a satisfactory preconditioner. With this context, one can observe the following weaknesses.

1) The literature on preconditioning is huge and classical texts such as Numerical Optimization by Nocedal and Wright discuss it in depth. When the main contribution is a preconditioner, not only does the author need to start with the classical practical preconditioners mentioned there, but also conduct a literature search to identify what's new and check if this specific preconditioner has already been proposed. On top of that, since so many preconditioners and variations are available, it is always possible to find something that works well for any handful of datasets and a theoretical analysis of when it is superior is needed.

2) The proposed preconditioner is a small variation of Kim et al. As far as I can tell, it adds absolute values and an additive epsilon term that is a fairly common techniques for making matrices more positive definite (or negative definite, depending on the sign).

3) The paper does not propose anything new on the privacy side. The framework of Kim et al. is re-used.

4) The paper offloads some computation on the client, specifically the computation of the preconditioner. It raises the question of what is reasonable to offload on the client and what is not.
- One possibility is to offload nothing onto the client -- the client simply uploads the dataset to the cloud for long term storage and the cloud servers can provide various analytics functionality in the future. However, for this setting the paper needs a privacy preserving preconditioner computation.
- On the other hand, if the client is already computing the preconditioner, why wouldn't the client run gradient descent (or sdg for large datasets) to train the logistic regression themselves? There are many software packages that do this and it is a very fast computation.

5) Most of the experiments in the paper suspiciously only focus on the training loss and the paper explicitly refuses to examine other metrics, such as training classification error, testing loss, and testing classification error. Since loss and classification error are not perfectly correlated, and can exhibit strange behaviors when an optimizer is changed, it raises red flags that the paper should avoid (by running these more complete experiments).

---

confidence: 2
contribution: 3
initial_rating: 5
presentation: 3
questions: Address W1, W2, and W3
soundness: 3
strengths: - S1: The proposed enhanced NAG method is well-analyzed theoretically.
- S2: The advantages of the enhanced NAG method are demonstrated empirically in both plaintext and ciphertext settings. Under homomorphic encryption, it achieves comparable accuracy to baselines with improved efficiency in learning time and iteration count.
- S3: The paper is well-organized and easy to follow.
summary: This paper introduces an enhanced NAG (Nesterov Accelerated Gradient) method for efficient training of logistic regression, specifically designed for secure training under homomorphic encryption. The empirical study shows that the proposed method achieves comparable results in just four iterations.
weaknesses: - W1: The related work section lacks sufficient detail. In particular, the relationship to other learning algorithms under FHE, beyond just NAG, should be more clearly described.
- W2: The experimental baseline is limited to a single method. Given the existence of various homomorphic computation methods, comparisons with these alternatives are necessary to establish the novelty of the proposed approach.
- W3: Although the paper shows that the proposed method accelerates learning computation in FHE, the improvement in performance is modest, and there is a slight drop in accuracy. The paper should also discuss potential applications where these specific efficiency gains are essential.