# Fast Training Logistic Regression via Adaptive Sampling

Authors: Yunsheng Song; Xiaohan Kong; Shuoping Huang; Chao Zhang
Journal: Scientific Programming
 (2021)

---

## Abstract

Logistic regression has been widely used in artificial intelligence and machine learning due to its deep theoretical basis and good practical performance. Its training process aims to solve a large-scale optimization problem characterized by a likelihood function, where the gradient descent approach is the most commonly used. However, when the data size is large, it is very time-consuming because it computes the gradient using all the training data in every iteration. Though this difficulty can be solved by random sampling, the appropriate sampled examples size is difficult to be predetermined and the obtained could be not robust. To overcome this deficiency, we propose a novel algorithm for fast training logistic regression via adaptive sampling. The proposed method decomposes the problem of gradient estimation into several subproblems according to its dimension; then, each subproblem is solved independently by adaptive sampling. Each element of the gradient estimation is obtained by successively sampling a fixed volume training example multiple times until it satisfies its stopping criteria. The final estimation is combined with the results of all the subproblems. It is proved that the obtained gradient estimation is a robust estimation, and it could keep the objective function value decreasing in the iterative calculation. Compared with the representative algorithms using random sampling, the experimental results show that this algorithm obtains comparable classification performance with much less training time.


---

## Body

## 1. Introduction

Supervised learning is to train a learner with a labelled training set and correctly determine the outputs for the unseen instances [1]. As one of the most famous classification algorithms in supervised learning, logistic regression learner is a generalized linear regression model in which the output is discrete [2]. Logistic regression has been widely used in various kinds of applications owing to its good performance, as process tomography [3], customer churn prediction [4], spatial prediction [5], major chronic diseases and clinical risk prediction [6, 7], and so on [8–10]. The process of training a logistic regression model aims to solve an unconstrained convex optimization problem, where gradient descent is one of the most important solutions [11]. Because of the computation of the gradient using the whole training instances, it is very time-consuming to use gradient descend (GD) when the data size is large [12].To speed up GD, many improved algorithms have been developed. According to the volume of data to obtain the gradient estimation, these algorithms can be divided into two groups: stochastic gradient descent and batch gradient descent [13]. Stochastic gradient descent (SGD) uses only one randomly selected training example to compute the gradient, and this can be very efficient for large datasets [14]. So, SGD is much faster and very suitable for online learning. However, the estimated gradient obtained by SGD is difficult to be a descent direction at each iteration, so that it needs a vast number of iterations. Furthermore, SGD is difficult to be suitable for the parallel environment [15].Different from SGD, batch gradient descent (BGD) obtains the aim gradient estimation using randomly choosing a certain amount of training examples. In this way, BGD could largely reduce the error and instability of the estimation, and it also obtains an effective solution [16–18]. As the sampled examples play an important role in estimating the gradient, BGD needs to carefully choose an appropriate sample size before sampling. However, it is difficult to predetermine an appropriate sample size for different datasets. Furthermore, samples of the same size could vary in terms of their qualities, because some examples are more representative or resembling the original data than others [19].This paper presents an improved adaptive sampling (AS) algorithm for accelerating the logistic regression training process. This method firstly gives a rule for estimating the gradient by some examples, and the obtained gradient estimation can guarantee that the objective function value keeps decreasing in the iterative calculation of GD. Then, the problem of obtaining an appropriate vector that meets the rule can be decomposed into several subproblems, where each subproblem determines a component of the vector that satisfies the stopping rule. Finally, the examples are drawn successively from the training set into the sample, and it terminates as soon as each component of the estimated gradient over the obtained sample satisfies its own rule. To speed up this process, the estimated components satisfying their own stopping rules are not estimated in the subsequent iteration, and they are the corresponding components of the final estimate of the gradient. The main contributions of this paper are as follows:(1)
Giving the rule to judge whether the direction of a vector is a descent direction of the current objective function or not, it is critical for the execution efficiency of the gradient descend method.(2)
Providing an adaptive sampling method to overcome the difficulty of the predetermining sample size before sampling, this method can adaptively determine the sample size according to the character of datasets and avoid the influence of human subjective factors.(3)
Applying a strategy of divide-and-conquer to efficiently obtain the gradient estimation on the sampled examples, the aim gradient vector estimation problem is divided into several one-dimensional estimation subproblems, and each subproblem can be solved independently.(4)
Proving the obtained gradient estimation is robust using probably approximately correct theory, and this estimation could be a descent direction of the current objective function at each iteration.(5)
Designing an efficient mechanism to solve the multivariate estimation problem for large-scale data.The rest of the paper is organized as follows: Section2 reviews related methods according to their characteristic. Section 3 proposes a sampling on-demand algorithm for logistic regression and proves its effectiveness. Section 4 reports the experimental results through the comparison with existing methods. Section 5 gives the conclusion of this paper and shows some future research work.


## 2. Related Work

Related work in improving GD has been widely developed nowadays. According to the amount of data to obtain the gradient estimation, existing GD algorithms can be divided into two groups: SGD algorithm and BGD algorithm.The original SGD (OSGD) algorithm computes the gradient with only one sample from the training set. The OSGD algorithm does not consider the effect of different dimensions on its convergence, so its rate of convergence could be slow when the surface of the objection function curves steeply for different dimensions. Qian [20] proposed the Momentum algorithm to an accelerated OSGD algorithm, where the current update vector is appended with a fraction of the obtained vector in the previous iteration. To efficiently solve the sparse data, Duchi et al. [21] proposed an Adagrad method to deal with the online learning task, and it set different learning rates for different components of the vector. Besides, there are lots of algorithms to determine the adaptive learning rate for different components [22–24]. These kinds of algorithms can deal with large-scale data, they must perform a vast number of iterations before an appreciable improvement of the objective function, and it is difficult to parallelize.Different from the SGD algorithm, the BGD algorithm computes the gradient with some randomly sampled examples from the training set at each iteration. So, the BGD algorithm can reduce the variance of the estimate of the gradient, and it achieves more stable convergence. In order to sufficiently often achieve convergence to the optimal solution, the estimate of the gradient by the BGD algorithm needs to enforce descent in the objective function at every iteration. Therefore, the sample size is carefully determined. Byrd et al. [25] proposed a dynamic sample gradient (DSG) algorithm, which can dynamically determine the sample size before sampling. For the convex optimization problem, the DSG algorithm can get the optimal solution. However, the sample size determined by the DSG algorithm could increase with the increasing steps, so that the total running time of the DSG algorithm also increases. Furthermore, owing to the fact that samples of the same size could vary in terms of their qualities, this leads to the estimate of the gradient with the fixed-size sample that may not enforce descent in the objective function. On the other hand, choosing a proper learning rate is an important issue for the performance of the BGD algorithm. A smaller learning rate could cause the convergence rate to become slower, but a larger one fluctuates obviously around the optimal solution. Robbins and Sutton [26] proposed a schedule to select the appropriate learning rate during training, where the predefined schedule is conducive to reducing the learning rate. Liang et al. [27] have proposed a sampling on-demand to speed up logistic regression, but the theoretical proof about its robust result is not given, and it does not compare its classification performance with other state-of-the-art approaches. There are many similar algorithms to yield high accuracy in the solution of the optimization problem [28].


## 3. Main Content

### 3.1. Preliminary

The logistic regression classier is generated using the posterior probabilities of two labels ({0, 1}) denoted by a linear function inx∈Rm, where the sum of these two posterior probabilities is one. The form of this model is that(1)logPy=1|xPy=0|x=ω⊤z,where the weight vector ω=ω1,ω2,…,ωm+1⊤∈Rm+1, z=1,x⊤∈Rm+1, ω⊤z is the inner product between ω and z, y is the label of x, and m is the dimensional size.LetTS=x1,y1,x2,y2,…,xN,yN be a training set, each training instance xii=1,2,…,N is represented by a m-dimensional vector xi1,xi2,…,xim⊤, and its label is yi∈0,1. Therefore, the optimal weight vector is obtained by minimizing the following problem:(2)minω∈Rm+1−1N∑i=1Nyiω⊤zi−log1+expω⊤zi,where zi=zi1,zi2,…,zi,m+1⊤=1,xi1,…,xim⊤. Clearly, the objective function Jω=−1/N∑i=1Nyiω⊤zi−log1+expω⊤zi is a strictly convex function. Let ω∗ be the optimal solution of the optimization problem (2), and the classifier of logistic regression is Ifz,0.5:(3)Ifz,0.5=1iffz≥0.5,0otherwise,where fz=1/1+expω∗⊤z. With these notations, the predicted label y^ for a given test instance x0=x01,x02,…,x0m can be derived from the function Ifz0,0.5, where z0=1,x01,x02,…,x0m⊤.The GD algorithm is an iterative optimization algorithm, where it updates the current solution with the solution solved in the previous step and the gradient of the current objective function each time. Owning to its advantage of simple implementation and effectiveness, GD is widely adopted to solve the unconstraint optimization problem (2). Let ωk be the optimal solution at the kth iteration. The gradient estimation ∇Jωk is obtained as follows:(4)∇Jωk=∇J1ωk,∇J2ωk,…,∇Jm+1ωk=1N∑i=1N∇l1ωk,zi,yi,…,1N∑i=1N∇lm+1ωk,zi,yi=1N∑i=1N∇l1ωk,zi,yi,…,∇lm+1ωk,zi,yi=1N∑i=1N∇lωk,zi,yi,where ∇ljωk,zi,yi=−expωk⊤zi/1+expωk⊤zi−yizij, j=1,2,…,m+1. Problem (2) is a convex optimization problem and the obtained solution is the global optimal solution. The computation of the gradient estimation needs all the data, so that it is very time-consuming for large-scale data. To overcome the deficiency, we propose a novel algorithm for fast training logistic regression via adaptive sampling (LLR-AS).


### 3.2. LLR-AS Algorithm

In fact, the LLR-AS algorithm is also a type of minibatch gradient descent algorithm. It randomly samples a subset fromT to compute the gradient estimation ∇J^ω to replace ∇Jω computed with all the data in every iteration. According to the feature of the GD algorithm, its convergence efficiency is closely related to the quality of the gradient estimation ∇Jω. The value of Jω is difficult to be smaller all the time and the final execution time becomes longer if the objective function value cannot keep decreasing using the gradient estimation in each iteration. Because gradient is the fastest direction in which the current objective function value decreases, a similar direction between ∇Jω,∇Jω, and ∇J^ω should be considered to achieve this aim. In the following, a rule is given to obtain an estimate ∇J^ω with high quality. This results in the following inequality:(5)1−ε∇Jω22≤∇J^ω⊤∇Jωwhere X2 is a norm of the vector X and the relaxation parameter 0≤ε<1/2. According to the conclusion [29], the vector −∇J^ω must keep the objective function Jω value decreasing in the iterative calculation of GD if it satisfies inequality (5). Moreover, parameter ε controls the descent speed for the function Jω. The larger the value of ε, the smaller the directional derivative, and the lower the descent speed for the function Jω. Because the optimization problem (2) is convex, the LLR-AS algorithm can get the globally optimal solution. With the rule, our algorithm is outlined in Algorithm 1.Algorithm 1:LLR-AS algorithm.
Input: Dataset TS=x1,y1,x1,y2,…,xN,yN, the stepsize λ>0.Output: The optimal vector ω∗.(1)
Initialize:ω0∈Rm+1 and k=0repeat(2)
Obtain the vector∇J^ωk satisfying inequality (5)(3)
ωk+1=ωk−λ∇J^ωk(4)
k⟶k+1until a convergence test is satisfied;(5)
Returnω∗=ωk


### 3.3. AS Algorithm

In order to get the vector∇J^ω meeting inequality (5), a novel adaptive sampling algorithm is proposed in the following: let T=z1,y1z2,y2,…,zN,yN be the new dataset generated by the set TS, where zi=1,xi1,…,xim⊤, i=1,2,…,N. According to inequality (5), the simplest way of getting the vector ∇J^ω is to sample a subset S from the set T and obtain an estimate from ∇Jω;S=1/S∑zi,yi∈S∇lωk,zi,yi as ∇J^ω, where S is the size of subset S. However, the sample size is difficult to determine for different tasks though there exist some theoretical and empirical results in the literature such as PAC [30] learning theory and learning curves [31]. Moreover, the theoretical results are usually worst-case and learning curves are average-case, so they are not necessarily consistent with each other [19]. To tackle this difficulty, we propose an adaptive sampling algorithm. Our method obtains the aim subset S by continually sampling examples from T until the gradient estimation on the sampled subset S satisfying the stopping rule. It decides the sample size through the information of the sampled examples and solves the difficulty of predetermining the sample size. Therefore, the key issue becomes the problem that how to design the stopping rule for our adaptive sampling to satisfy inequality (5).From a statistical point of view, the estimation of gradient∇Jω satisfying inequality (5) is a m+1-dimensional vector estimating problem. However, the existing sampling procedures mainly focus on a one-dimensional estimate problem, and they cannot be directly applied to the multidimensional problem. Although there exists a close relationship between the components of the gradient ∇Jω, each component of ∇Jω can still be seen as a one-dimensional estimating subproblem. Therefore, the multivariate estimation problem (5) can be solved by solving these one-dimensional estimating subproblems. Inequality (5) is equivalent to ∑j=1m+1∇Jjω;S×∇Jjω−1−ε∇Jjω2≥0 according to the formula of vector inner product. In other words, inequality (5) must hold if each component ∇Jjω;S of ∇Jjω;S simultaneously satisfies its own inequality 1−ε∇Jjω2≤∇Jjω;S×∇Jjω, where j=1,2,…,m+1. So, the problem of seeking for the gradient estimation satisfying inequality (5) can be approximatively divided into m+1 subproblems, where each subproblem is solved by ∇Jjω;S×∇Jjω−1−ε∇Jjω2≥0 at the same time.Each component∇Jjω of ∇Jω can be considered as the exception of the population T, and ∇Jjω,S is the estimation of ∇Jjω on the subset S sampled from the population T. According to the central limit theorem, the difference of value between ∇Jjω,S and ∇Jjω continually becomes small with the enlarging sampled subset S. There exists a critical value of the sample subset S size for satisfying inequality 1−ε∇Jjω2≤∇Jjω;S×∇Jjω with a large probability. Therefore, both the early stopping rule and consecutive sampling are adopted to get the sampled subset S and estimation over it. Given the aim vector ω and the objective function Jω, each element ∇Jjω can be a constant value and its upper bound of the absolute value could be estimated by a function αjt,s=djln6m+1/δ/2t3s3j=1,2,…,m+1, where dj=2min1≤i≤Nzi,j+max1≤i≤Nzi,j and t and s are the size of sampled subset per sampling and total number of sampling, 0<δ<1/2. In the next section, we will show that such a function can achieve this goal, and the stopping condition can be finally transformed into ∇Jjω;S≥1+1/εαjt,s for each component ∇Jjω. However, it needs to recompute ∇Jω;S on the sampled subset S and test whether all the components ∇Jjω;S satisfy their own stopping rule during each round of sampling, and this computational burden could cost much more execution time to achieve this aim. Two improvements are made to solve this problem (Algorithm 2).Algorithm 2:AS algorithm.
Input: Dataset T=z1,y2,z1,y2,…,zN,yN, a weight vector ω, parameters ε∈0,1/2,δ∈0,1/2.Output: The vector ∇J^ω=∇J^1ω,∇J^2ω,/p…,∇J^m+1ω.(1)
Initialize:Λ=1,2,…,m+1, S˜0=∅, t=0(2)
Computedj=2min1≤i≤Nzi,j+max1≤i≤Nzi,j, where j=1,2,…,m+1whileΛ≥1do(3)
t=t+1(4)
Sample a random subsetSt with the size of s from the set T, S˜t=S˜t−1∪Stfor eachj∈Λdo(5)
Compute∇Jjω;S˜tif∇Jjω;S˜t≥1+1/εdjln6m+6/δ/2s3t3then(6)
Λ=Λ−j(7)
∇J^jω=∇Jjω;S˜tendendend(8)
Return∇J^ω=∇J^1ω,∇J^2ω,/p…,∇J^m+1ωLetS˜t−1 be the cumulatively sampled subset obtained by the first t−1 iterations sampling. The set S˜t=S˜t−1∪St, where St is the sampled subset at the t−th iteration. We have(6)∇Jω;S˜t=1S˜t∑zi,yi∈S˜tlω,zi,yi=1S˜t∑zi,yi∈S˜t−1lω,zi,yi+1S˜t∑zi,yi∈St−1lω,zi,yi=S˜t−1S˜t1S˜t−1∑zi,yi∈S˜t−1lω,zi,yi+1S˜t∑zi,yi∈St−1lω,zi,yi=S˜t−1S˜t−1+St∇Jω;S˜t+1S˜t−1+St∑zi,yi∈St−1lω,zi,yi.The computation of ∇Jω;S˜t can be largely reduced using the result on the set S˜t−1 at each iteration. On the other hand, we adopt an asynchronous way to get each component of ∇Jω;S˜. If one or more components satisfy their own stopping rules, then they are directly the corresponding components of the final result without considering in subsequent iteration. Algorithm 2 describes the multivariate adaptive sampling method.


### 3.4. The Effectiveness of AS Algorithm

In this subsection, we study the effectiveness of AS algorithm. To derive our main result, we need the following lemmas and theorems.Lemma 1.
LetA1,A2,…,Am+1 be m+1 random events. For any positive integer m,(7)P∩i=1m+1Ai≥∑i=1m+1PAi−m.Proof.
The lemma is proved by mathematical induction.(1)
Considerm=1 [32]. We have A1∩A2c⊆A2c and PA1∩A2c≤PA2c for m=1, where Ac is the complementary set of the set A. Meanwhile, PA1∩A2+PA2c≥PA1∩A2+PA1∩A2c=PA1. So, PA1∩A2≥PA1−PA2c=PA1−1−PA2=PA1+PA2−1.(2)
Assume that the inequality holds form=k−1; that is, P∩i=1kAi≥∑i=1kPAi−k+1.(3)
Considerm=k. We have(8)P∩j=1k+1Ai=P∩j=1kAi∩Ak+1≥P∩j=1kAi+PAk+1−1≥∑j=1kPAj−k+1+PAk+1−1≥∑j=1k+1PAj−k.
The lemma follows immediately from mathematical induction.Lemma 2.
LetS˜⊆T be a subset, which is obtained through independent and random sampling from the training set T with the size of n. For any given weight vector ω, λ≥0 and j∈1,2,…,m+1, we have(9)P∇Jjω;S˜−∇Jjω≥djλ≤exp−2nλ2.Proof.
According to∇ljω,zi,yi=1/1+expω⊤zi−yizi,j and yi∈0,1, we have minaj,−bj≤∇ljω,zi,yi≤max−aj,bj, where aj=min1≤i≤Nzi,j, bj=max1≤i≤Nzi,j, j=1,2,…,m+1. Moreover,(10)max1≤i≤N∇ljω,zi,yi−min1≤i≤N∇ljω,zi,yi≤max−aj,bj−minaj,−bj≤maxbj−aj,−2aj,2bj≤2aj+bj=dj.
Therefore, for anyj∈1,2,…,m+1, it follows from inequality (10) and the Hoeffding inequality [32] that(11)P∇Jjω;S−∇Jjω≥djλ≤exp−2djλ2max1≤i≤N∇ljω;Xi−min1≤i≤N∇ljω;Xi2n≤exp−2nλ2.
For the convenience of the following developments, we make the following remarks. Lettj be the number of iterations until the AS algorithm satisfies the inequality ∇J^jω;S˜nt≥1+1/εαjt, where S˜nt⊆T is a subset obtained through random and independent sampling from T with the size of nt=t×s, j=1,2,…,m+1. Note that tj is a random variable depending on the sample drawn from T. Let tj1 be the smallest integer meeting the following inequalities:(12)αjtj1<ε∇Jjω1+2εandαjtj1−1≥ε∇Jjω1+2ε.Since αjt is a strictly decreasing function with t and Jjω is fixed under any given ω, hence, tj1 is uniquely determined for j=1,2,…,m+1.Lemma 3.
Iftj<tj1, then we have ∇Jjω;S˜tj∇Jjω≥1−εJjω2 with probability >1−δ/2+2m, where j=1,2,…,m+1, 0<δ<1/2.Proof.
We get from the AS algorithm the estimate∇Jjω;S˜tj satisfying ∇Jjω;S˜tj≥1+1/εαjtj at the tj-th step. When tj<tj1, we have αjtj1≥αjtj1−1≥ε∇Jjω/1+2ε and(13)∇Jjω;S˜tj≥1+1εαjtj≥1+1εαjtj1−1≥1+1/εε∇Jjω1+2ε=1−ε1+2ε∇Jjω≥1−ε∇Jjω.
Thus, it always holds that∇Jjω;S˜tj∇Jjω≥1−ε∇Jjω2 as long as ∇Jjω;S˜tj and ∇Jjω have the same sign. On the other hand, if ∇Jjω;S˜tj and ∇Jjω have different signs, then they are quite different as ∇Jjω;S˜tj≥αjtj1+1/ε. Next, we show the difference is large enough, and it is conducive to prove the probability that this situation occurs is small.
In the following, we will give this difference from two cases. If∇Jjω;S˜tj≤0 and ∇Jjω≥0, then ∇Jjω;S˜tj−∇Jjω=−∇Jjω;S˜tj+∇Jjω≥−∇Jjω;S˜tj=∇Jjω;S˜tj≥αtj1+1/ε≥αtj. If ∇Jjω;S˜tj≥0 and ∇Jjω≤0, then ∇Jjω;S˜tj−∇Jjω=∇Jjω;S˜tj−∇Jjω≥∇Jjω;S˜tj=∇Jjω;S˜tj≥αjtj1+1/ε≥αjtj. Therefore, the difference is ∇Jjω;S˜tj−∇Jjω≥αjtj as ∇Jjω;S˜tj∇Jjω≤0. Then, the probability that ∇Jjω;S˜tj and ∇Jjω have different labels is estimated using Lemma 2 as follows:(14)P∇Jjω;S˜tj∇Jjω≤0≤P∇Jjω;S˜tj−∇Jjω≥αjtj=P∇Jjω;S˜tj−∇Jjω≥djln6m+6/δ2tj·s3≤exp−2tj·sln6m+6/δ2tj·s3=δ6m+6exp1tj·s3≤δ6m+6exp1≤δ2m+2.Lemma 4.
Ptj<tj1≥1−δ/2+2m for j=1,2,…,m+1.Proof.
Whentj<tj1, it always holds that ∇Jjω;S˜tj<1+1/εαjtj1 and αjtj1<ε1+2/ε∇Jjω. Then,(15)Ptj<tj1≤P∇Jjω;S˜tj<1+1εαjtj1=P∇Jjω;S˜tj<1+2εεαjtj1−αjtj1≤P∇Jjω;S˜tj<∇Jjω|−αjtj1=P∇Jjω|−∇Jjω;S˜tj>αjtj1.
It follows from the triangle inequality that∇Jjω−∇Jjω;S˜tj≥∇Jjω−∇Jjω;S˜tj>αjtj1. Thus, P∇Jjω−∇Jjω;S˜tj≥αjtj1≥PJjω−∇Jjω;S˜tj>αjtj1. Combining inequality (15) with Lemma 2, we have(16)Ptj<tj1=1−Ptj≥tj1≥1−P∇Jjω;S˜tj−∇Jjω>αjtj1≥1−P∇Jjω−∇Jjω;S˜tj>αjtj1≥1−δ2+2m.Combining with Lemma 3 and 4, we can easily get the following theorem.Theorem 1.
For any1/2>ε>0 and 1/2>δ>0, the estimation ∇Jjω;S˜tj generated by AS algorithm satisfies the following inequality:(17)P∇Jjω;S˜tj∇Jjω≥1−εJjω2≥1−δ1+m,where j=1,2,…,m+1.Theorem 2.
For any1/2>ε>0 and 1/2>δ>0, final estimation ∇J^ω=∇J1ω;S˜t1, …,∇Jm+1ω;S˜tm+1 generated by the AS algorithm satisfies inequality P1−ε∇Jω22≤∇J^ω⊤∇Jω≥1−δ.Proof.
Combining Theorem1 and Lemma 1, then we have(18)P1−ε∇Jω22≤∇J^ω⊤∇Jω=P∑j=1m+1∇Jjω;S˜tj∇Jjω≥1−ε∑j=1m+1∇Jjω2≥P∩j=1m+1∇Jjω;S˜tj∇Jjω≥1−ε∇Jjω2≥∑j=1m+1P∇Jjω;S˜tj∇Jjω≥1−ε∇Jjω2−m≥∑j=1m+11−δm+1−m=1−δ.
In Theorem2, the obtained estimation using the AS algorithm could keep the decreasing value of the current objective function in each iteration calculation of GD. Therefore, it could guarantee that LLR-AS algorithm gets the optimal solution of the convex optimization problem (2), and this conclusion is verified in the experiment. Besides the optimal solution, we also pay attention to the number of sampled examples (NSE). AS algorithm is an iterative sampling algorithm that samples a fixed number of samples each time, then we can estimate NSE using the total number of iterations. We have already proved that AS algorithm terminates finally within t1 steps from Lemmas 3 and 4, where t1=max1≤j≤m+1tj1. It is well known that t1 could be the minimum number of iterations satisfying condition (12); then, we could assume that α1j≈ε∇Jjω/1+2ε, where j=1,2⋯,m+1. Finally, we can estimate NSE as follows:(19)NSE=s·t1≤s·max1≤j≤m+131+2εdj2ε∇Jjω2ln6m+6δ+1.
Next, we will discuss the above formula. The effect of the parameterm and δ on formula (19) is small because they are in the logarithmic function. Namely, the AS algorithm has a low possibility to sample too many examples. So, it is useful for sampling examples from the large-scale data.


## 3.1. Preliminary

The logistic regression classier is generated using the posterior probabilities of two labels ({0, 1}) denoted by a linear function inx∈Rm, where the sum of these two posterior probabilities is one. The form of this model is that(1)logPy=1|xPy=0|x=ω⊤z,where the weight vector ω=ω1,ω2,…,ωm+1⊤∈Rm+1, z=1,x⊤∈Rm+1, ω⊤z is the inner product between ω and z, y is the label of x, and m is the dimensional size.LetTS=x1,y1,x2,y2,…,xN,yN be a training set, each training instance xii=1,2,…,N is represented by a m-dimensional vector xi1,xi2,…,xim⊤, and its label is yi∈0,1. Therefore, the optimal weight vector is obtained by minimizing the following problem:(2)minω∈Rm+1−1N∑i=1Nyiω⊤zi−log1+expω⊤zi,where zi=zi1,zi2,…,zi,m+1⊤=1,xi1,…,xim⊤. Clearly, the objective function Jω=−1/N∑i=1Nyiω⊤zi−log1+expω⊤zi is a strictly convex function. Let ω∗ be the optimal solution of the optimization problem (2), and the classifier of logistic regression is Ifz,0.5:(3)Ifz,0.5=1iffz≥0.5,0otherwise,where fz=1/1+expω∗⊤z. With these notations, the predicted label y^ for a given test instance x0=x01,x02,…,x0m can be derived from the function Ifz0,0.5, where z0=1,x01,x02,…,x0m⊤.The GD algorithm is an iterative optimization algorithm, where it updates the current solution with the solution solved in the previous step and the gradient of the current objective function each time. Owning to its advantage of simple implementation and effectiveness, GD is widely adopted to solve the unconstraint optimization problem (2). Let ωk be the optimal solution at the kth iteration. The gradient estimation ∇Jωk is obtained as follows:(4)∇Jωk=∇J1ωk,∇J2ωk,…,∇Jm+1ωk=1N∑i=1N∇l1ωk,zi,yi,…,1N∑i=1N∇lm+1ωk,zi,yi=1N∑i=1N∇l1ωk,zi,yi,…,∇lm+1ωk,zi,yi=1N∑i=1N∇lωk,zi,yi,where ∇ljωk,zi,yi=−expωk⊤zi/1+expωk⊤zi−yizij, j=1,2,…,m+1. Problem (2) is a convex optimization problem and the obtained solution is the global optimal solution. The computation of the gradient estimation needs all the data, so that it is very time-consuming for large-scale data. To overcome the deficiency, we propose a novel algorithm for fast training logistic regression via adaptive sampling (LLR-AS).


## 3.2. LLR-AS Algorithm

In fact, the LLR-AS algorithm is also a type of minibatch gradient descent algorithm. It randomly samples a subset fromT to compute the gradient estimation ∇J^ω to replace ∇Jω computed with all the data in every iteration. According to the feature of the GD algorithm, its convergence efficiency is closely related to the quality of the gradient estimation ∇Jω. The value of Jω is difficult to be smaller all the time and the final execution time becomes longer if the objective function value cannot keep decreasing using the gradient estimation in each iteration. Because gradient is the fastest direction in which the current objective function value decreases, a similar direction between ∇Jω,∇Jω, and ∇J^ω should be considered to achieve this aim. In the following, a rule is given to obtain an estimate ∇J^ω with high quality. This results in the following inequality:(5)1−ε∇Jω22≤∇J^ω⊤∇Jωwhere X2 is a norm of the vector X and the relaxation parameter 0≤ε<1/2. According to the conclusion [29], the vector −∇J^ω must keep the objective function Jω value decreasing in the iterative calculation of GD if it satisfies inequality (5). Moreover, parameter ε controls the descent speed for the function Jω. The larger the value of ε, the smaller the directional derivative, and the lower the descent speed for the function Jω. Because the optimization problem (2) is convex, the LLR-AS algorithm can get the globally optimal solution. With the rule, our algorithm is outlined in Algorithm 1.Algorithm 1:LLR-AS algorithm.
Input: Dataset TS=x1,y1,x1,y2,…,xN,yN, the stepsize λ>0.Output: The optimal vector ω∗.(1)
Initialize:ω0∈Rm+1 and k=0repeat(2)
Obtain the vector∇J^ωk satisfying inequality (5)(3)
ωk+1=ωk−λ∇J^ωk(4)
k⟶k+1until a convergence test is satisfied;(5)
Returnω∗=ωk


## 3.3. AS Algorithm

In order to get the vector∇J^ω meeting inequality (5), a novel adaptive sampling algorithm is proposed in the following: let T=z1,y1z2,y2,…,zN,yN be the new dataset generated by the set TS, where zi=1,xi1,…,xim⊤, i=1,2,…,N. According to inequality (5), the simplest way of getting the vector ∇J^ω is to sample a subset S from the set T and obtain an estimate from ∇Jω;S=1/S∑zi,yi∈S∇lωk,zi,yi as ∇J^ω, where S is the size of subset S. However, the sample size is difficult to determine for different tasks though there exist some theoretical and empirical results in the literature such as PAC [30] learning theory and learning curves [31]. Moreover, the theoretical results are usually worst-case and learning curves are average-case, so they are not necessarily consistent with each other [19]. To tackle this difficulty, we propose an adaptive sampling algorithm. Our method obtains the aim subset S by continually sampling examples from T until the gradient estimation on the sampled subset S satisfying the stopping rule. It decides the sample size through the information of the sampled examples and solves the difficulty of predetermining the sample size. Therefore, the key issue becomes the problem that how to design the stopping rule for our adaptive sampling to satisfy inequality (5).From a statistical point of view, the estimation of gradient∇Jω satisfying inequality (5) is a m+1-dimensional vector estimating problem. However, the existing sampling procedures mainly focus on a one-dimensional estimate problem, and they cannot be directly applied to the multidimensional problem. Although there exists a close relationship between the components of the gradient ∇Jω, each component of ∇Jω can still be seen as a one-dimensional estimating subproblem. Therefore, the multivariate estimation problem (5) can be solved by solving these one-dimensional estimating subproblems. Inequality (5) is equivalent to ∑j=1m+1∇Jjω;S×∇Jjω−1−ε∇Jjω2≥0 according to the formula of vector inner product. In other words, inequality (5) must hold if each component ∇Jjω;S of ∇Jjω;S simultaneously satisfies its own inequality 1−ε∇Jjω2≤∇Jjω;S×∇Jjω, where j=1,2,…,m+1. So, the problem of seeking for the gradient estimation satisfying inequality (5) can be approximatively divided into m+1 subproblems, where each subproblem is solved by ∇Jjω;S×∇Jjω−1−ε∇Jjω2≥0 at the same time.Each component∇Jjω of ∇Jω can be considered as the exception of the population T, and ∇Jjω,S is the estimation of ∇Jjω on the subset S sampled from the population T. According to the central limit theorem, the difference of value between ∇Jjω,S and ∇Jjω continually becomes small with the enlarging sampled subset S. There exists a critical value of the sample subset S size for satisfying inequality 1−ε∇Jjω2≤∇Jjω;S×∇Jjω with a large probability. Therefore, both the early stopping rule and consecutive sampling are adopted to get the sampled subset S and estimation over it. Given the aim vector ω and the objective function Jω, each element ∇Jjω can be a constant value and its upper bound of the absolute value could be estimated by a function αjt,s=djln6m+1/δ/2t3s3j=1,2,…,m+1, where dj=2min1≤i≤Nzi,j+max1≤i≤Nzi,j and t and s are the size of sampled subset per sampling and total number of sampling, 0<δ<1/2. In the next section, we will show that such a function can achieve this goal, and the stopping condition can be finally transformed into ∇Jjω;S≥1+1/εαjt,s for each component ∇Jjω. However, it needs to recompute ∇Jω;S on the sampled subset S and test whether all the components ∇Jjω;S satisfy their own stopping rule during each round of sampling, and this computational burden could cost much more execution time to achieve this aim. Two improvements are made to solve this problem (Algorithm 2).Algorithm 2:AS algorithm.
Input: Dataset T=z1,y2,z1,y2,…,zN,yN, a weight vector ω, parameters ε∈0,1/2,δ∈0,1/2.Output: The vector ∇J^ω=∇J^1ω,∇J^2ω,/p…,∇J^m+1ω.(1)
Initialize:Λ=1,2,…,m+1, S˜0=∅, t=0(2)
Computedj=2min1≤i≤Nzi,j+max1≤i≤Nzi,j, where j=1,2,…,m+1whileΛ≥1do(3)
t=t+1(4)
Sample a random subsetSt with the size of s from the set T, S˜t=S˜t−1∪Stfor eachj∈Λdo(5)
Compute∇Jjω;S˜tif∇Jjω;S˜t≥1+1/εdjln6m+6/δ/2s3t3then(6)
Λ=Λ−j(7)
∇J^jω=∇Jjω;S˜tendendend(8)
Return∇J^ω=∇J^1ω,∇J^2ω,/p…,∇J^m+1ωLetS˜t−1 be the cumulatively sampled subset obtained by the first t−1 iterations sampling. The set S˜t=S˜t−1∪St, where St is the sampled subset at the t−th iteration. We have(6)∇Jω;S˜t=1S˜t∑zi,yi∈S˜tlω,zi,yi=1S˜t∑zi,yi∈S˜t−1lω,zi,yi+1S˜t∑zi,yi∈St−1lω,zi,yi=S˜t−1S˜t1S˜t−1∑zi,yi∈S˜t−1lω,zi,yi+1S˜t∑zi,yi∈St−1lω,zi,yi=S˜t−1S˜t−1+St∇Jω;S˜t+1S˜t−1+St∑zi,yi∈St−1lω,zi,yi.The computation of ∇Jω;S˜t can be largely reduced using the result on the set S˜t−1 at each iteration. On the other hand, we adopt an asynchronous way to get each component of ∇Jω;S˜. If one or more components satisfy their own stopping rules, then they are directly the corresponding components of the final result without considering in subsequent iteration. Algorithm 2 describes the multivariate adaptive sampling method.


## 3.4. The Effectiveness of AS Algorithm

In this subsection, we study the effectiveness of AS algorithm. To derive our main result, we need the following lemmas and theorems.Lemma 1.
LetA1,A2,…,Am+1 be m+1 random events. For any positive integer m,(7)P∩i=1m+1Ai≥∑i=1m+1PAi−m.Proof.
The lemma is proved by mathematical induction.(1)
Considerm=1 [32]. We have A1∩A2c⊆A2c and PA1∩A2c≤PA2c for m=1, where Ac is the complementary set of the set A. Meanwhile, PA1∩A2+PA2c≥PA1∩A2+PA1∩A2c=PA1. So, PA1∩A2≥PA1−PA2c=PA1−1−PA2=PA1+PA2−1.(2)
Assume that the inequality holds form=k−1; that is, P∩i=1kAi≥∑i=1kPAi−k+1.(3)
Considerm=k. We have(8)P∩j=1k+1Ai=P∩j=1kAi∩Ak+1≥P∩j=1kAi+PAk+1−1≥∑j=1kPAj−k+1+PAk+1−1≥∑j=1k+1PAj−k.
The lemma follows immediately from mathematical induction.Lemma 2.
LetS˜⊆T be a subset, which is obtained through independent and random sampling from the training set T with the size of n. For any given weight vector ω, λ≥0 and j∈1,2,…,m+1, we have(9)P∇Jjω;S˜−∇Jjω≥djλ≤exp−2nλ2.Proof.
According to∇ljω,zi,yi=1/1+expω⊤zi−yizi,j and yi∈0,1, we have minaj,−bj≤∇ljω,zi,yi≤max−aj,bj, where aj=min1≤i≤Nzi,j, bj=max1≤i≤Nzi,j, j=1,2,…,m+1. Moreover,(10)max1≤i≤N∇ljω,zi,yi−min1≤i≤N∇ljω,zi,yi≤max−aj,bj−minaj,−bj≤maxbj−aj,−2aj,2bj≤2aj+bj=dj.
Therefore, for anyj∈1,2,…,m+1, it follows from inequality (10) and the Hoeffding inequality [32] that(11)P∇Jjω;S−∇Jjω≥djλ≤exp−2djλ2max1≤i≤N∇ljω;Xi−min1≤i≤N∇ljω;Xi2n≤exp−2nλ2.
For the convenience of the following developments, we make the following remarks. Lettj be the number of iterations until the AS algorithm satisfies the inequality ∇J^jω;S˜nt≥1+1/εαjt, where S˜nt⊆T is a subset obtained through random and independent sampling from T with the size of nt=t×s, j=1,2,…,m+1. Note that tj is a random variable depending on the sample drawn from T. Let tj1 be the smallest integer meeting the following inequalities:(12)αjtj1<ε∇Jjω1+2εandαjtj1−1≥ε∇Jjω1+2ε.Since αjt is a strictly decreasing function with t and Jjω is fixed under any given ω, hence, tj1 is uniquely determined for j=1,2,…,m+1.Lemma 3.
Iftj<tj1, then we have ∇Jjω;S˜tj∇Jjω≥1−εJjω2 with probability >1−δ/2+2m, where j=1,2,…,m+1, 0<δ<1/2.Proof.
We get from the AS algorithm the estimate∇Jjω;S˜tj satisfying ∇Jjω;S˜tj≥1+1/εαjtj at the tj-th step. When tj<tj1, we have αjtj1≥αjtj1−1≥ε∇Jjω/1+2ε and(13)∇Jjω;S˜tj≥1+1εαjtj≥1+1εαjtj1−1≥1+1/εε∇Jjω1+2ε=1−ε1+2ε∇Jjω≥1−ε∇Jjω.
Thus, it always holds that∇Jjω;S˜tj∇Jjω≥1−ε∇Jjω2 as long as ∇Jjω;S˜tj and ∇Jjω have the same sign. On the other hand, if ∇Jjω;S˜tj and ∇Jjω have different signs, then they are quite different as ∇Jjω;S˜tj≥αjtj1+1/ε. Next, we show the difference is large enough, and it is conducive to prove the probability that this situation occurs is small.
In the following, we will give this difference from two cases. If∇Jjω;S˜tj≤0 and ∇Jjω≥0, then ∇Jjω;S˜tj−∇Jjω=−∇Jjω;S˜tj+∇Jjω≥−∇Jjω;S˜tj=∇Jjω;S˜tj≥αtj1+1/ε≥αtj. If ∇Jjω;S˜tj≥0 and ∇Jjω≤0, then ∇Jjω;S˜tj−∇Jjω=∇Jjω;S˜tj−∇Jjω≥∇Jjω;S˜tj=∇Jjω;S˜tj≥αjtj1+1/ε≥αjtj. Therefore, the difference is ∇Jjω;S˜tj−∇Jjω≥αjtj as ∇Jjω;S˜tj∇Jjω≤0. Then, the probability that ∇Jjω;S˜tj and ∇Jjω have different labels is estimated using Lemma 2 as follows:(14)P∇Jjω;S˜tj∇Jjω≤0≤P∇Jjω;S˜tj−∇Jjω≥αjtj=P∇Jjω;S˜tj−∇Jjω≥djln6m+6/δ2tj·s3≤exp−2tj·sln6m+6/δ2tj·s3=δ6m+6exp1tj·s3≤δ6m+6exp1≤δ2m+2.Lemma 4.
Ptj<tj1≥1−δ/2+2m for j=1,2,…,m+1.Proof.
Whentj<tj1, it always holds that ∇Jjω;S˜tj<1+1/εαjtj1 and αjtj1<ε1+2/ε∇Jjω. Then,(15)Ptj<tj1≤P∇Jjω;S˜tj<1+1εαjtj1=P∇Jjω;S˜tj<1+2εεαjtj1−αjtj1≤P∇Jjω;S˜tj<∇Jjω|−αjtj1=P∇Jjω|−∇Jjω;S˜tj>αjtj1.
It follows from the triangle inequality that∇Jjω−∇Jjω;S˜tj≥∇Jjω−∇Jjω;S˜tj>αjtj1. Thus, P∇Jjω−∇Jjω;S˜tj≥αjtj1≥PJjω−∇Jjω;S˜tj>αjtj1. Combining inequality (15) with Lemma 2, we have(16)Ptj<tj1=1−Ptj≥tj1≥1−P∇Jjω;S˜tj−∇Jjω>αjtj1≥1−P∇Jjω−∇Jjω;S˜tj>αjtj1≥1−δ2+2m.Combining with Lemma 3 and 4, we can easily get the following theorem.Theorem 1.
For any1/2>ε>0 and 1/2>δ>0, the estimation ∇Jjω;S˜tj generated by AS algorithm satisfies the following inequality:(17)P∇Jjω;S˜tj∇Jjω≥1−εJjω2≥1−δ1+m,where j=1,2,…,m+1.Theorem 2.
For any1/2>ε>0 and 1/2>δ>0, final estimation ∇J^ω=∇J1ω;S˜t1, …,∇Jm+1ω;S˜tm+1 generated by the AS algorithm satisfies inequality P1−ε∇Jω22≤∇J^ω⊤∇Jω≥1−δ.Proof.
Combining Theorem1 and Lemma 1, then we have(18)P1−ε∇Jω22≤∇J^ω⊤∇Jω=P∑j=1m+1∇Jjω;S˜tj∇Jjω≥1−ε∑j=1m+1∇Jjω2≥P∩j=1m+1∇Jjω;S˜tj∇Jjω≥1−ε∇Jjω2≥∑j=1m+1P∇Jjω;S˜tj∇Jjω≥1−ε∇Jjω2−m≥∑j=1m+11−δm+1−m=1−δ.
In Theorem2, the obtained estimation using the AS algorithm could keep the decreasing value of the current objective function in each iteration calculation of GD. Therefore, it could guarantee that LLR-AS algorithm gets the optimal solution of the convex optimization problem (2), and this conclusion is verified in the experiment. Besides the optimal solution, we also pay attention to the number of sampled examples (NSE). AS algorithm is an iterative sampling algorithm that samples a fixed number of samples each time, then we can estimate NSE using the total number of iterations. We have already proved that AS algorithm terminates finally within t1 steps from Lemmas 3 and 4, where t1=max1≤j≤m+1tj1. It is well known that t1 could be the minimum number of iterations satisfying condition (12); then, we could assume that α1j≈ε∇Jjω/1+2ε, where j=1,2⋯,m+1. Finally, we can estimate NSE as follows:(19)NSE=s·t1≤s·max1≤j≤m+131+2εdj2ε∇Jjω2ln6m+6δ+1.
Next, we will discuss the above formula. The effect of the parameterm and δ on formula (19) is small because they are in the logarithmic function. Namely, the AS algorithm has a low possibility to sample too many examples. So, it is useful for sampling examples from the large-scale data.


## 4. Experiments

### 4.1. Experimental Setup

Two representative gradient descend methods were selected in this study: a common gradient descend with all the data and a dynamic sample gradient algorithm [25]. These two algorithms are used to solve the logistic regression, and they are named LR-GD and LR-DSG. Seven benchmark datasets are selected for making a fair comparison between our proposal and others [33, 34]; their information is shown in Table 1. All of these selected datasets have larger than 50000 instances.Table 1
Summary of datasets.
DatasetSizeFeaturesClassCifa6000030722Cod-rna5953582Covtype581012542Ijcnn1141691222Mnist350000852Skin-nonskin24505732Susy5000000182Owning to the simplicity and successful application, we select the classification accuracy (Acc) and training time as the performance measure. LLR-AS and LR-DSG are both accelerated algorithms of LR-GD, so we compare their relative speeds (RS) and the difference in their solutions. RS is a ratio of training time between LLR-AS and each of the others. For estimating these three performance measures Acc and RS, we used a 10-fold cross-validation method. To compare the difference between our method and others under a performance measure, we adopt the Wilcoxon signed-rank test (WSRT) [35]. The reason for selecting WSRT is that it does not require a strict data distribution hypothesis and has stable performance. It is empirically considered to be stronger than other tests [36]. The null hypothesis of WSRT denotes that there exists no significant difference between our algorithm and each one of the others under a performance measure, while the alternative is that there exists a significant difference.In the following experiments, we fix an initial sampleS0 of size 1% of the total training set and θ=0.5 for DSG algorithm according to [25], and ε=0.5, δ=0.1 for LLR-AS. The stopping condition for these three iterative algorithms is that the Euclidean distance between the current and the previous solution is smaller than 0.001 and the maximum iteration 5000. The significance level of 0.05 is used. All the experiments are executed in Python 3.8 on the same computer of Intel Xeon E5-2650 CPU and 32 GB of RAM.


### 4.2. Experimental Results and Analysis

In this section, we adopt the classification ability and training efficiency as two important measurements to evaluate the performance of these algorithms and give a detailed comparing analysis and the reason for the experiment result.


#### 4.2.1. Classification Performance Analysis

Under the theoretical hypothesis of logistic regression, its classification performance of logistic regression depends on the solution obtained by gradient descend. LR-GD algorithm is trained by gradient descend with all the training data; then its weight vector is the optimal solution, as well as its classification performance. In the following, we compare the predictive ability of these algorithms from the difference in the obtained solution vector and classification accuracy.(1) The analysis of the difference on solution vector: Pearson correlation coefficient ρ∈0,1 is chosen to evaluate the difference between two solution vectors for its high effectiveness. Its value is inversely proportional to the difference. The larger the value of ρ, the smaller difference between the two vectors. The correlation coefficients between the solution vector obtained by LLP-AS and each one of LR-GD and LR-DSG algorithms on the test data of each dataset are computed, all the statistics results on different datasets are listed in Table 2.Table 2
The correlation coefficientρ between different solutions on seven datasets.
Datasetρ (LLR-AS, LR-GD)ρ (LLR-AS, LR-DSG)StdMeanStdMeanCifa0.0060.9350.0030.986Cod-rna0.0030.9980.0030.991Covtype0.0070.9790.0480.984Ijcnn10.0010.9960.0041.000Mnist0.0000.9950.0030.985Skin-nonskin0.0001.0000.0040.979Susy0.0010.9980.0090.998Average0.9860.989Median0.9960.986Table2 shows that the correlation coefficient ρ of the weight vector between LLR-AS and LR-GD is nearly closed to 1 on almost all datasets except for the Cifa dataset, and this same result can be obtained between LLR-AS and LR-DSG. The mean of correlation coefficients on all the datasets is 0.986 and 0.989, and their medians are 0.996 and 0.986. Moreover, it can get a detailed comparison from the descriptive statistics. The following can be seen: (1) there exists a negligible difference of solution vector between the LLR-AS algorithm and each one of these two algorithms, and then the LLR-AS algorithm can get nearly the same solution vector as the LR-GD algorithm. (2) The standard deviation of the correlation coefficients on each data is tiny; then, this experiment result validates that the proposed algorithm could get a robust solution.Own to the properties of the obtained vector estimation at each iteration, the LLR-AS algorithm performs multiple iterations to continually minimize the objective function value. Meanwhile, the original optimization problem has a unique optimal solution because it is a convex optimization. So, the LLR-AS algorithm is able to guarantee convergence and obtain the optimal solution as the LR-GD algorithm. Furthermore, Theorem2 has also been verified by this experiment result, and the gradient estimation is stable for different datasets.(2) The analysis of the difference in classification. the classification accuracy Acc is adopted to compare the performance between the proposed algorithm and two state-of-the-art approaches, and WSRT is performed to test whether there exists a significant difference among them. Table 3 lists the descriptive statistics of classification accuracy of each algorithm obtained by 10-fold cross-validation, and their results are also plotted in Figure 1.Table 3
Acc of three algorithms on seven datasets.
DatasetLR-GDLR-DSGLLR-ASStdMeanStdMeanStdMeanCifa0.0030.9040.0030.9010.0030.904Cod-rna0.1420.7820.0050.8880.0600.851Covtype0.1030.6650.1170.6820.1450.655Ijcnn10.0040.7940.0040.7940.0110.841Mnist0.0020.7110.0030.7130.0050.722Skin-nonskin0.0020.5830.0020.5730.0050.589Susy0.0020.6800.0020.6810.0050.689Average0.7310.7470.750Median0.7110.7110.722Figure 1
Acc of different algorithms on different datasets.The result on each dataset plotted in Figure1 shows that the LLR-AS achieves better classification accuracy than the LR-GD algorithm and LR-DSG algorithm on the Ijcnn1 dataset, and it has not the worst classification accuracy on the rest of the datasets. To assess the overall classification performance on all the datasets, the mean and median of the result of each algorithm on eight datasets are computed in the last two rows of Table 3. Their mean values of classification accuracy are 0.750, 0.731, and 0.747, and the median values are 0.722, 0.711, and 0.713. Therefore, there exists a negligible difference in classification accuracy among these algorithms. Finally, the obtained p values using WSRT between LLR-AS and each one of the other algorithms are 0.1563 and 0.688, both larger than the given significant level of 0.05. Then, it gets that (1) the LLR-AS algorithm has no significant difference in classification accuracy with the LR-GD algorithm and LR-DSG algorithm on the selected datasets. (2) The LLR-AS algorithm has a stable classification performance because its standard deviation of the classification accuracy of the LLR-AS algorithm on every data is relatively small.The reason for the similar classification result is that the classification performance of logistic regression depends on the solution vector, and the LLR-AS algorithm has no significant different solution vector from the LR-GD algorithm and the LR-DSG algorithm. Moreover, the obtained solution vector of the algorithm is robust according to Theorem2, and the small standard deviation of the correlation coefficient on different datasets also verifies this fact.


#### 4.2.2. Training Efficiency Analysis

Besides the classification performance, training speed is another important measurement to evaluate algorithm training performance. The relative speedRS can evaluate the accelerating extent of the LR-DSG algorithm. Table 4 lists RS on different datasets.Table 4
The relative speed on seven datasets.
DatasetLLR-AS vs LR-GDLLR-AS vs LR-DSGStdMeanStdMeanCifa3.049109.3390.0031.461Cod-rna3.456111.2600.2033.067Covtype9.659270.8980.3483.945Ijcnn14.64164.1320.0041.070Mnist5.01483.3770.0031.082Skin-nonskin1.03620.2820.0041.243Susy4.88791.9420.0091.031Average107.3191.843It finds that the value ofRS between the LLR-AS algorithm and LR-GD algorithm is larger than 20 on all the datasets from Table 4, and its value is larger than 109 on Cifa, Cod-rna, and Covtype dataset. So, the LLR-AS algorithm can largely reduce the training time of the LR-GD algorithm. On the other hand, the value of RS between the LLR-AS algorithm and LR-DSG algorithm is larger than one on these seven datasets, and its average value of RS on all the datasets is 1.843. Therefore, the LLR-AS algorithm indeed needs less training time than the LR-DSG algorithm.There exist three reasons for explaining that the proposed mechanism can achieve a good result on training efficiency. (1) The obtained gradient estimation could be a descent direction of the current objective function at each iteration. Thus, the total number of iterations that is positively correlated with the training efficiency will reduce. (2) The divide-and-conquer approach is adopted to compute each component of the gradient vector at each iteration, and it can be executed in a parallel environment. (3) It is proved that there exists a low possibility to sample too many examples to estimate the gradient; the time of estimating gradient could become short at each iteration. Therefore, the proposed algorithm has a better performance than other representative algorithms owning to the above three merits.


## 4.1. Experimental Setup

Two representative gradient descend methods were selected in this study: a common gradient descend with all the data and a dynamic sample gradient algorithm [25]. These two algorithms are used to solve the logistic regression, and they are named LR-GD and LR-DSG. Seven benchmark datasets are selected for making a fair comparison between our proposal and others [33, 34]; their information is shown in Table 1. All of these selected datasets have larger than 50000 instances.Table 1
Summary of datasets.
DatasetSizeFeaturesClassCifa6000030722Cod-rna5953582Covtype581012542Ijcnn1141691222Mnist350000852Skin-nonskin24505732Susy5000000182Owning to the simplicity and successful application, we select the classification accuracy (Acc) and training time as the performance measure. LLR-AS and LR-DSG are both accelerated algorithms of LR-GD, so we compare their relative speeds (RS) and the difference in their solutions. RS is a ratio of training time between LLR-AS and each of the others. For estimating these three performance measures Acc and RS, we used a 10-fold cross-validation method. To compare the difference between our method and others under a performance measure, we adopt the Wilcoxon signed-rank test (WSRT) [35]. The reason for selecting WSRT is that it does not require a strict data distribution hypothesis and has stable performance. It is empirically considered to be stronger than other tests [36]. The null hypothesis of WSRT denotes that there exists no significant difference between our algorithm and each one of the others under a performance measure, while the alternative is that there exists a significant difference.In the following experiments, we fix an initial sampleS0 of size 1% of the total training set and θ=0.5 for DSG algorithm according to [25], and ε=0.5, δ=0.1 for LLR-AS. The stopping condition for these three iterative algorithms is that the Euclidean distance between the current and the previous solution is smaller than 0.001 and the maximum iteration 5000. The significance level of 0.05 is used. All the experiments are executed in Python 3.8 on the same computer of Intel Xeon E5-2650 CPU and 32 GB of RAM.


## 4.2. Experimental Results and Analysis

In this section, we adopt the classification ability and training efficiency as two important measurements to evaluate the performance of these algorithms and give a detailed comparing analysis and the reason for the experiment result.


### 4.2.1. Classification Performance Analysis

Under the theoretical hypothesis of logistic regression, its classification performance of logistic regression depends on the solution obtained by gradient descend. LR-GD algorithm is trained by gradient descend with all the training data; then its weight vector is the optimal solution, as well as its classification performance. In the following, we compare the predictive ability of these algorithms from the difference in the obtained solution vector and classification accuracy.(1) The analysis of the difference on solution vector: Pearson correlation coefficient ρ∈0,1 is chosen to evaluate the difference between two solution vectors for its high effectiveness. Its value is inversely proportional to the difference. The larger the value of ρ, the smaller difference between the two vectors. The correlation coefficients between the solution vector obtained by LLP-AS and each one of LR-GD and LR-DSG algorithms on the test data of each dataset are computed, all the statistics results on different datasets are listed in Table 2.Table 2
The correlation coefficientρ between different solutions on seven datasets.
Datasetρ (LLR-AS, LR-GD)ρ (LLR-AS, LR-DSG)StdMeanStdMeanCifa0.0060.9350.0030.986Cod-rna0.0030.9980.0030.991Covtype0.0070.9790.0480.984Ijcnn10.0010.9960.0041.000Mnist0.0000.9950.0030.985Skin-nonskin0.0001.0000.0040.979Susy0.0010.9980.0090.998Average0.9860.989Median0.9960.986Table2 shows that the correlation coefficient ρ of the weight vector between LLR-AS and LR-GD is nearly closed to 1 on almost all datasets except for the Cifa dataset, and this same result can be obtained between LLR-AS and LR-DSG. The mean of correlation coefficients on all the datasets is 0.986 and 0.989, and their medians are 0.996 and 0.986. Moreover, it can get a detailed comparison from the descriptive statistics. The following can be seen: (1) there exists a negligible difference of solution vector between the LLR-AS algorithm and each one of these two algorithms, and then the LLR-AS algorithm can get nearly the same solution vector as the LR-GD algorithm. (2) The standard deviation of the correlation coefficients on each data is tiny; then, this experiment result validates that the proposed algorithm could get a robust solution.Own to the properties of the obtained vector estimation at each iteration, the LLR-AS algorithm performs multiple iterations to continually minimize the objective function value. Meanwhile, the original optimization problem has a unique optimal solution because it is a convex optimization. So, the LLR-AS algorithm is able to guarantee convergence and obtain the optimal solution as the LR-GD algorithm. Furthermore, Theorem2 has also been verified by this experiment result, and the gradient estimation is stable for different datasets.(2) The analysis of the difference in classification. the classification accuracy Acc is adopted to compare the performance between the proposed algorithm and two state-of-the-art approaches, and WSRT is performed to test whether there exists a significant difference among them. Table 3 lists the descriptive statistics of classification accuracy of each algorithm obtained by 10-fold cross-validation, and their results are also plotted in Figure 1.Table 3
Acc of three algorithms on seven datasets.
DatasetLR-GDLR-DSGLLR-ASStdMeanStdMeanStdMeanCifa0.0030.9040.0030.9010.0030.904Cod-rna0.1420.7820.0050.8880.0600.851Covtype0.1030.6650.1170.6820.1450.655Ijcnn10.0040.7940.0040.7940.0110.841Mnist0.0020.7110.0030.7130.0050.722Skin-nonskin0.0020.5830.0020.5730.0050.589Susy0.0020.6800.0020.6810.0050.689Average0.7310.7470.750Median0.7110.7110.722Figure 1
Acc of different algorithms on different datasets.The result on each dataset plotted in Figure1 shows that the LLR-AS achieves better classification accuracy than the LR-GD algorithm and LR-DSG algorithm on the Ijcnn1 dataset, and it has not the worst classification accuracy on the rest of the datasets. To assess the overall classification performance on all the datasets, the mean and median of the result of each algorithm on eight datasets are computed in the last two rows of Table 3. Their mean values of classification accuracy are 0.750, 0.731, and 0.747, and the median values are 0.722, 0.711, and 0.713. Therefore, there exists a negligible difference in classification accuracy among these algorithms. Finally, the obtained p values using WSRT between LLR-AS and each one of the other algorithms are 0.1563 and 0.688, both larger than the given significant level of 0.05. Then, it gets that (1) the LLR-AS algorithm has no significant difference in classification accuracy with the LR-GD algorithm and LR-DSG algorithm on the selected datasets. (2) The LLR-AS algorithm has a stable classification performance because its standard deviation of the classification accuracy of the LLR-AS algorithm on every data is relatively small.The reason for the similar classification result is that the classification performance of logistic regression depends on the solution vector, and the LLR-AS algorithm has no significant different solution vector from the LR-GD algorithm and the LR-DSG algorithm. Moreover, the obtained solution vector of the algorithm is robust according to Theorem2, and the small standard deviation of the correlation coefficient on different datasets also verifies this fact.


### 4.2.2. Training Efficiency Analysis

Besides the classification performance, training speed is another important measurement to evaluate algorithm training performance. The relative speedRS can evaluate the accelerating extent of the LR-DSG algorithm. Table 4 lists RS on different datasets.Table 4
The relative speed on seven datasets.
DatasetLLR-AS vs LR-GDLLR-AS vs LR-DSGStdMeanStdMeanCifa3.049109.3390.0031.461Cod-rna3.456111.2600.2033.067Covtype9.659270.8980.3483.945Ijcnn14.64164.1320.0041.070Mnist5.01483.3770.0031.082Skin-nonskin1.03620.2820.0041.243Susy4.88791.9420.0091.031Average107.3191.843It finds that the value ofRS between the LLR-AS algorithm and LR-GD algorithm is larger than 20 on all the datasets from Table 4, and its value is larger than 109 on Cifa, Cod-rna, and Covtype dataset. So, the LLR-AS algorithm can largely reduce the training time of the LR-GD algorithm. On the other hand, the value of RS between the LLR-AS algorithm and LR-DSG algorithm is larger than one on these seven datasets, and its average value of RS on all the datasets is 1.843. Therefore, the LLR-AS algorithm indeed needs less training time than the LR-DSG algorithm.There exist three reasons for explaining that the proposed mechanism can achieve a good result on training efficiency. (1) The obtained gradient estimation could be a descent direction of the current objective function at each iteration. Thus, the total number of iterations that is positively correlated with the training efficiency will reduce. (2) The divide-and-conquer approach is adopted to compute each component of the gradient vector at each iteration, and it can be executed in a parallel environment. (3) It is proved that there exists a low possibility to sample too many examples to estimate the gradient; the time of estimating gradient could become short at each iteration. Therefore, the proposed algorithm has a better performance than other representative algorithms owning to the above three merits.


## 4.2.1. Classification Performance Analysis

Under the theoretical hypothesis of logistic regression, its classification performance of logistic regression depends on the solution obtained by gradient descend. LR-GD algorithm is trained by gradient descend with all the training data; then its weight vector is the optimal solution, as well as its classification performance. In the following, we compare the predictive ability of these algorithms from the difference in the obtained solution vector and classification accuracy.(1) The analysis of the difference on solution vector: Pearson correlation coefficient ρ∈0,1 is chosen to evaluate the difference between two solution vectors for its high effectiveness. Its value is inversely proportional to the difference. The larger the value of ρ, the smaller difference between the two vectors. The correlation coefficients between the solution vector obtained by LLP-AS and each one of LR-GD and LR-DSG algorithms on the test data of each dataset are computed, all the statistics results on different datasets are listed in Table 2.Table 2
The correlation coefficientρ between different solutions on seven datasets.
Datasetρ (LLR-AS, LR-GD)ρ (LLR-AS, LR-DSG)StdMeanStdMeanCifa0.0060.9350.0030.986Cod-rna0.0030.9980.0030.991Covtype0.0070.9790.0480.984Ijcnn10.0010.9960.0041.000Mnist0.0000.9950.0030.985Skin-nonskin0.0001.0000.0040.979Susy0.0010.9980.0090.998Average0.9860.989Median0.9960.986Table2 shows that the correlation coefficient ρ of the weight vector between LLR-AS and LR-GD is nearly closed to 1 on almost all datasets except for the Cifa dataset, and this same result can be obtained between LLR-AS and LR-DSG. The mean of correlation coefficients on all the datasets is 0.986 and 0.989, and their medians are 0.996 and 0.986. Moreover, it can get a detailed comparison from the descriptive statistics. The following can be seen: (1) there exists a negligible difference of solution vector between the LLR-AS algorithm and each one of these two algorithms, and then the LLR-AS algorithm can get nearly the same solution vector as the LR-GD algorithm. (2) The standard deviation of the correlation coefficients on each data is tiny; then, this experiment result validates that the proposed algorithm could get a robust solution.Own to the properties of the obtained vector estimation at each iteration, the LLR-AS algorithm performs multiple iterations to continually minimize the objective function value. Meanwhile, the original optimization problem has a unique optimal solution because it is a convex optimization. So, the LLR-AS algorithm is able to guarantee convergence and obtain the optimal solution as the LR-GD algorithm. Furthermore, Theorem2 has also been verified by this experiment result, and the gradient estimation is stable for different datasets.(2) The analysis of the difference in classification. the classification accuracy Acc is adopted to compare the performance between the proposed algorithm and two state-of-the-art approaches, and WSRT is performed to test whether there exists a significant difference among them. Table 3 lists the descriptive statistics of classification accuracy of each algorithm obtained by 10-fold cross-validation, and their results are also plotted in Figure 1.Table 3
Acc of three algorithms on seven datasets.
DatasetLR-GDLR-DSGLLR-ASStdMeanStdMeanStdMeanCifa0.0030.9040.0030.9010.0030.904Cod-rna0.1420.7820.0050.8880.0600.851Covtype0.1030.6650.1170.6820.1450.655Ijcnn10.0040.7940.0040.7940.0110.841Mnist0.0020.7110.0030.7130.0050.722Skin-nonskin0.0020.5830.0020.5730.0050.589Susy0.0020.6800.0020.6810.0050.689Average0.7310.7470.750Median0.7110.7110.722Figure 1
Acc of different algorithms on different datasets.The result on each dataset plotted in Figure1 shows that the LLR-AS achieves better classification accuracy than the LR-GD algorithm and LR-DSG algorithm on the Ijcnn1 dataset, and it has not the worst classification accuracy on the rest of the datasets. To assess the overall classification performance on all the datasets, the mean and median of the result of each algorithm on eight datasets are computed in the last two rows of Table 3. Their mean values of classification accuracy are 0.750, 0.731, and 0.747, and the median values are 0.722, 0.711, and 0.713. Therefore, there exists a negligible difference in classification accuracy among these algorithms. Finally, the obtained p values using WSRT between LLR-AS and each one of the other algorithms are 0.1563 and 0.688, both larger than the given significant level of 0.05. Then, it gets that (1) the LLR-AS algorithm has no significant difference in classification accuracy with the LR-GD algorithm and LR-DSG algorithm on the selected datasets. (2) The LLR-AS algorithm has a stable classification performance because its standard deviation of the classification accuracy of the LLR-AS algorithm on every data is relatively small.The reason for the similar classification result is that the classification performance of logistic regression depends on the solution vector, and the LLR-AS algorithm has no significant different solution vector from the LR-GD algorithm and the LR-DSG algorithm. Moreover, the obtained solution vector of the algorithm is robust according to Theorem2, and the small standard deviation of the correlation coefficient on different datasets also verifies this fact.


## 4.2.2. Training Efficiency Analysis

Besides the classification performance, training speed is another important measurement to evaluate algorithm training performance. The relative speedRS can evaluate the accelerating extent of the LR-DSG algorithm. Table 4 lists RS on different datasets.Table 4
The relative speed on seven datasets.
DatasetLLR-AS vs LR-GDLLR-AS vs LR-DSGStdMeanStdMeanCifa3.049109.3390.0031.461Cod-rna3.456111.2600.2033.067Covtype9.659270.8980.3483.945Ijcnn14.64164.1320.0041.070Mnist5.01483.3770.0031.082Skin-nonskin1.03620.2820.0041.243Susy4.88791.9420.0091.031Average107.3191.843It finds that the value ofRS between the LLR-AS algorithm and LR-GD algorithm is larger than 20 on all the datasets from Table 4, and its value is larger than 109 on Cifa, Cod-rna, and Covtype dataset. So, the LLR-AS algorithm can largely reduce the training time of the LR-GD algorithm. On the other hand, the value of RS between the LLR-AS algorithm and LR-DSG algorithm is larger than one on these seven datasets, and its average value of RS on all the datasets is 1.843. Therefore, the LLR-AS algorithm indeed needs less training time than the LR-DSG algorithm.There exist three reasons for explaining that the proposed mechanism can achieve a good result on training efficiency. (1) The obtained gradient estimation could be a descent direction of the current objective function at each iteration. Thus, the total number of iterations that is positively correlated with the training efficiency will reduce. (2) The divide-and-conquer approach is adopted to compute each component of the gradient vector at each iteration, and it can be executed in a parallel environment. (3) It is proved that there exists a low possibility to sample too many examples to estimate the gradient; the time of estimating gradient could become short at each iteration. Therefore, the proposed algorithm has a better performance than other representative algorithms owning to the above three merits.


## 5. Summary

A novel algorithm for fast training logistic regression via adaptive sampling has been proposed to effectively handle the massive dataset in this paper. The proposed algorithm solves the difficulty that the sample size needs to be fixed before sampling, and it also offers an idea of dividing the multivariate estimation problem into several easy-to-solve subproblems. Experimental results on real datasets demonstrate that LLR-AS has obtained a similar classification performance with less execution time in comparison with other representative algorithms. Moreover, this proposed algorithm can deal with the multiclassification problem using the one-vs-all scheme, and paper [37] has shown that this scheme is as accurate as any other approach. The proposed algorithm solves the binary classification problem, but it can be used for the multiclassification problem using the one-vs-all scheme. It needs to train several different classifiers, where each classifier is obtained by distinguishing the instances of the same class from the others in the rest classes. When given an unlabeled instance, the final output is the largest result among the results of all the classifiers.Though the proposed algorithm has a good performance for large-scale data, there exist two limitations for dealing with various kinds of real datasets. The gradient estimation needs all the features of the data at each iteration, so that it may take a great challenge of its training efficiency for high dimensional data [38–40]. Furthermore, this algorithm does not consider the label distribution of the data, and then its performance on imbalanced data could decrease. In the future, we will study how to combine sampling and feature selection to scale up machine learning algorithms and design an effective mechanism to deal with the class imbalance problem.



---

*Source: 9991859-2021-05-30.xml*